# ARMA Modelling of Sleep Spindles

João Costa, Manuel Duarte Ortigueira, Arnaldo Batista

► To cite this version:

João Costa, Manuel Duarte Ortigueira, Arnaldo Batista. ARMA Modelling of Sleep Spindles. 2nd Doctoral Conference on Computing, Electrical and Industrial Systems (DoCEIS), Feb 2011, Costa de Caparica, Portugal. pp.341-348, 10.1007/978-3-642-19170-1_37. hal-01566572

HAL Id: hal-01566572

https://inria.hal.science/hal-01566572

Submitted on 21 Jul 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# ARMA Modelling of Sleep Spindles

João Caldas da Costa¹, Manuel Duarte Ortigueira², Arnaldo Batista²

¹ Department of Systems and Informatics, EST, IPS, Setubal, Portugal

² UNINOVA and Department of Electrical Engineering, University Nova, Lisbon, Portugal

**Abstract.** Differences in EEG sleep spindles constitute a promising indicator of neurodegenerative disorders. In this paper an ARMA modelling to sleep spindles is proposed and tested. The primary objective is to distinguish, via poles and zeros location, between regular, elderly and dementia subjects. In order to achieve this goal, a model validation has been done.

**Keywords:** sleep spindles, ARMA, EEG

## 1 Introduction

Sleep spindles are particular EEG patterns which occur during the sleep cycle with center frequency in the band 11.5 to 15 Hz. They are used as one of the features to classify the sleep stages [1]. Sleep spindles are promising objective indicators in neurodegenerative disorders [2]. In order to interpret them, their structure needs to be clarified or a suitable model needs to be found. In this work, an autoregressive moving average (ARMA) model for sleep spindles is used to detect meaningful differences when applied to spindles from different types of people. More clearly, we wish to distinguish normal, elderly and dementia subjects based on the location of poles and zeros.

In [3] automated spindle detection by using autoregressive (AR) modelling for feature extraction is proposed. It is concluded that AR model parameters provide a good representation of the EEG data. It is expected that even better results can be obtained from the use of ARMA models.

In order to validate the ARMA models, a system is created and its response to white noise is obtained. Then, a model is estimated using the previous response and the estimated model is compared with the original one. It is also studied the best order of the model in order to represent a sleep spindle.

## 2 Contribution to sustainability

Sustainability is to promote the best for people and environment, both now and in the future [8], and contributions to early diagnosis of diseases can lead to a better tomorrow. This paper comes with this perspective. The objective is an early detection of changes in brain to prevent or, at least, mitigate the influence of certain diseases.
# 3 Sleep Spindles

It is commonly referred in literature that sleep spindles are the most interesting hallmark of stage 2 sleep electroencephalograms (EEG) [1]. A sleep spindle is a burst of brain activity visible on an EEG and it consists of 11-15 Hz waves with duration between 0.5s and 2s in healthy adults, they are bilateral and synchronous in their appearance, with amplitude up to 30 µV (Fig.1).

<figure><img src="image_2.png" /><figcaption>Fig. 1. EEG signal showing a sleep spindle.</figcaption></figure>

The spindle is characterized by progressively increasing, then gradually decreasing amplitude, which gives the waveform its characteristic name [4]. It is now reliable that sleep spindles are originated in the thalamus and can be recorded as potential changes at the cortical surface [5].

Sleep spindles are affected by brain pathology, as well as by normal and pathological aging (e.g., dementia) [1]. With normal aging, sleep spindles are less numerous and less well formed. In dementia, the sleep EEG patterns suggest accelerated aging [6].

Sleep EEG measures seem promising as objective indicators in neurodegenerative disorders, including dementia, where sleep changes appear to be an exaggeration of changes that come normally with aging.

# 4 ARMA models and “Itakura-Saito” distance

## 4.1 ARMA model

In signal processing, autoregressive moving average (ARMA) models are typically applied to correlated time series data. Given a time series, we can consider it as the output of an ARMA system driven by white noise. The ARMA model is a tool for understanding and, whenever necessary, predicting future values in time series. The model consists of two parts, an autoregressive (AR) part and a moving average (MA) part. The model is usually referred to as ARMA(p,q) where p is the order of the autoregressive part and q is the order of the moving average part.
Compared with the pure MA or AR models, ARMA models more suitable for describing the characteristics of a given process with minimum number of parameters using both poles and zeros, rather than just poles or zeros [7].

As referred, a stationary ARMA process of order (p,q) is considered as the output of a linear time-invariant(LTI) digital filter driven by white noise. The transfer function of the system is given by:

$$H(z)=\frac{\sum_{m=0}^{q} b_m z^{-m}}{\sum_{k=0}^{p} a_k z^{-k}} \quad (1)$$

with $a_0=1$. The process corresponding to this model satisfies the difference equation:

$$x(n) = - \sum_{k=1}^{p} a_k x(n-k) + \sum_{m=0}^{q} b_m w(n-m) \quad (2)$$

where $w(n)$ is the input sequence, a zero-mean white noise and $x(n)$ is the output sequence. The main task in the modeling can be formulated as:

Given a segment of a time series, $x(n)$, n=0,1,2 ..., L-1, estimate the p+q+1 ARMA parameters.

## 4.2 ARMA model validation

In order to validate the accuracy of the ARMA models some tests have been made.

Two models were tested, with 5 and 3 poles respectively and both with one zero. The models were excited with white Gaussian noise. An ARMA model was then estimated based only on the output of the model. The correct model orders are assumed to be known. The procedure was repeated several times (100) and means were calculated.

The original models used had the following transfer functions:

$$H1(z)=\frac{1-0.2 z^{-1}}{1+0.1 z^{-1}+0.02 z^{-2}+0.154 z^{-3}+0.1597 z^{-4}+0.5584 z^{-5}} \quad (3)$$

$$H2(z)=\frac{1+0.9 z^{-1}}{1-z^{-1}+0.66 z^{-2}-0.4 z^{-3}} \quad (4)$$

It can be seen, from (Tables 1 and 2) and from pole-zero map (Fig. 2 and 4) that the estimators produced very accurate results. In (Fig. 3 and 5) the clusters for poles and zeros positions are shown (these are the locations of all the poles and zeros in all the experiments).
<table><caption>**Table 1.** H1(z) - ARMA(5,1) coefficients</caption><thead><tr><th></th><th>a<sub>0</sub></th><th>a<sub>1</sub></th><th>b<sub>0</sub></th><th>b<sub>1</sub></th><th>b<sub>2</sub></th><th>b<sub>3</sub></th><th>b<sub>4</sub></th><th>b<sub>5</sub></th></tr></thead><tbody><tr><td>Original</td><td>1.0000</td><td>-0.2000</td><td>1.0000</td><td>0.1000</td><td>0.0200</td><td>0.1540</td><td>0.1597</td><td>0.5584</td></tr><tr><td>Estimated(mean)</td><td>1.0000</td><td>-0.1973</td><td>1.0000</td><td>0.1027</td><td>0.0211</td><td>0.1514</td><td>0.1633</td><td>0.5602</td></tr><tr><td>Error</td><td>0.0000</td><td>0.0027</td><td>0.0000</td><td>0.0027</td><td>0.0011</td><td>0.0026</td><td>0.0036</td><td>0.0017</td></tr><tr><td>Mean of errors</td><td>0.0000</td><td>0.0481</td><td>0.0000</td><td>0.0402</td><td>0.0233</td><td>0.0198</td><td>0.0215</td><td>0.0221</td></tr><tr><td>Quadratic error</td><td>0.0000</td><td>0.0036</td><td>0.0000</td><td>0.0023</td><td>0.0009</td><td>0.0005</td><td>0.0007</td><td>0.0008</td></tr></tbody></table>

<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Zeros and poles from original and estimated (mean) ARMA(5,1) systems.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Clusters of zeros and poles from estimated ARMA(5,1) system.</figcaption></figure>

<table><caption>**Table 2.** H2(z) - ARMA(5,1) coefficients</caption><thead><tr><th></th><th>a<sub>0</sub></th><th>a<sub>1</sub></th><th>b<sub>0</sub></th><th>b<sub>1</sub></th><th>b<sub>2</sub></th><th>b<sub>3</sub></th></tr></thead><tbody><tr><td>Original</td><td>1.0000</td><td>0.9000</td><td>1.0000</td><td>-1.0000</td><td>0.6600</td><td>-0.4000</td></tr><tr><td>Estimated(mean)</td><td>1.0000</td><td>0.8736</td><td>1.0000</td><td>-1.0043</td><td>0.6737</td><td>-0.4045</td></tr><tr><td>Error</td><td>0.0000</td><td>0.0264</td><td>0.0000</td><td>0.0043</td><td>0.0137</td><td>0.0045</td></tr><tr><td>Mean of Errors</td><td>0.0000</td><td>0.0314</td><td>0.0000</td><td>0.0267</td><td>0.0345</td><td>0.0216</td></tr><tr><td>Quadratic error</td><td>0.0000</td><td>0.0017</td><td>0.0000</td><td>0.0011</td><td>0.0019</td><td>0.0007</td></tr></tbody></table>

<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Zeros and poles from original and estimated (mean) ARMA(3,1) systems.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Clusters of zeros and poles from estimated ARMA(3,1) system.</figcaption></figure>

In (Fig. 6) the Spectra from H2(z) and its corresponding ARMA(3,1) model is shown. It can be seen that both spectra are almost identical.

Tests have also been carried to determine the order of the model to be used in sleep spindle modelling. In (Figs. 4, 5, 7 and 8) the poles and zeros maps of 4 systems with different orders (numerator and denominator) are shown.
<figure><img src="image_7.png" /><figcaption>Fig. 6. Spectra from H2(z) and it's corresponding ARMA(3,1) model.</figcaption></figure>

In (Fig. 6) the Spectra from H2(z) and it's corresponding ARMA(3,1) model is shown. It can be seen that both spectra are almost identical.

Tests have also been carried to determine the order of the model to be used in sleep spindle modelling. In (Figs. 4, 5, 7 and 8) the poles and zeros maps of 4 systems with different orders (numerator and denominator) are shown.

For systems with orders larger than 5 poles and 1 zero, the new poles or zeros tend to “accommodate” themselves to the system with minor differences in the overall model. For example, when one more zero is added, only the position of the other zero suffers notorious change to accommodate the new pole, with small variations in poles positions (Figs. 8 and 10). On the other hand, when we increase simultaneously the pole and zero orders, extra pole/zero pairs appear in very close positions or in reverse positions, revealing the presence of allpass subsystems.

<figure><img src="image_8.png" /><figcaption>Fig. 7. Poles and zeros map of a spindle ARMA model with 5 poles, 0 zeros.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 8. Poles and zeros map of a spindle ARMA model with 5 poles, 1 zeros.</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>Fig. 9. Poles and zeros map of a spindle ARMA model with 6 poles, 1 zeros.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 10. Poles and zeros map of a spindle ARMA model with 5 poles, 2 zeros.</figcaption></figure>

## 4.3 “Itakura-Saito” distances

The “Itakura-Saito” distance is a measure of the perceptual difference between original spectrum P(W) and an approximation,

$\hat{P}(w)$, of that spectrum. It can be used to compare the coefficients of the AR polynomials. $I_i$ is defined as:

$$ D_{IS}(P(w), \hat{P}^i(w)) = \frac{1}{2\pi} \int_{-\pi}^{\pi} \left[ \frac{P(w)}{\hat{P}^i(w)} - \log \frac{P(w)}{\hat{P}^i(w)} \right] dw \quad (5) $$

# 5 Experimental results

Spindles from night sleep of 5 subjects were used. Three sets of spindles from healthy subjects (S1, S2 and S3), a set of spindles from an elderly healthy subject (ELD) and a set of spindles from a dementia patient (DEM). The data used is from a real EEG with 512 Hz sampling rate. It has been pre-processed with a band-pass filter with cutoff frequencies of 5Hz and 22Hz.

For each person, the same procedure has been applied, consisting of:

* Visual identification of the sleep spindles;

* Estimation of an ARMA model with 5 poles and 1 zero, thus, obtaining A and B polynomials; the mean of A and B polynomials obtained from each set of spindles was computed;

* Zeros and Poles map of all systems were obtained;

* For computing the “Itakura-Saito” distances, the real poles and zeros were removed

It is possible to distinguish, either by the analysis of poles map (Fig. 11) or by “Itakura-Saito” distances healthy subjects from dementia/elderly subjects. It is particularly notorious the “zero” position, which in the elderly/dementia subjects is very close to -1. On the other side, normal subjects “zero” is found to be located on the right hand side complex plane (Fig. 8).
<figure><img src="image_12.png" /><figcaption>**Fig. 11.** Zeros and poles from elderly and normal subjects.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Fig. 12.** Poles from elderly, dementia and normal subjects.</figcaption></figure>

However, it is not possible to distinguish between elderly and dementia subjects as the poles and zeros in both cases are very close to each others, as it can be seen from (Fig. 12).

The same result can be obtained by “Itakura-Saito” distance, in (Table 3) the distances between various subjects is showed. It is clearly seen that bigger distances (distance>0.1) are measured between elderly/dementia and normal subjects.

<table><caption>**Table 3.** “Itakura-Saito” distances between complex poles of subjects coefficients</caption><thead><tr><th></th><th>ELD</th><th>S1</th><th>S2</th><th>S3</th></tr></thead><tbody><tr><th>DEM</th><td>0.0109</td><td>0.2235</td><td>0.1420</td><td>0.1044</td></tr><tr><th>S3</th><td>0.1037</td><td>0.0365</td><td>0.0134</td><td></td></tr><tr><th>S2</th><td>0.1580</td><td>0.0777</td><td></td><td></td></tr><tr><th>S1</th><td>0.1992</td><td></td><td></td><td></td></tr></tbody></table>

It can be seen (Fig. 11 and 12) that the pole position give some biomarker for the presence of dementia. Complex poles from normal subject lie in specific areas, different from complex poles from elderly and dementia patients. However, pole location from elderly and dementia patients lie in similar regions inside the unit circle.

<figure><img src="image_14.png" /><figcaption>**Fig. 13.** Spectra from the ARMA model and from the original signal.</figcaption></figure>

In (Fig. 13) signal spectra corresponding ARMA model and to the periodogram estimate are shown. As it can be seen they give us similar information and from
them we can conclude that the center frequency is 12.5 Hz, similar to a sinusoid in
the 11-14 Hz band.

6 Conclusions

ARMA models can make a good representation of sleep spindles. It is showed that it
is possible to distinguish between regular subjects and elderly or dementia subjects.
However, it is not easy, using this method to distinguish between elderly and
dementia subjects. According to [9,10] there is a increased loss of spindles in
dementia patients when compared to elderly healthy subjects. From the experiments
we performed it seems not to be possible to distinguish different abnormalities in the
brain, probably because the effect on each individual spindle is similar. This requires
further research.

This type of spindle modeling opens a door into the perspective of using it in
the automatic spindle detection.

References

1. De Gennaro, L., Ferrara, M. Sleep spindles: an overview. Sleep Med Rev 2003;7:423-40.

2. Ktonas, P.Y., Golemati, S., Xanthopoulos, P., Sakkalis, V., Ortigueira, M.D, et al. Time-frequency analysis methods to quantify the time-varying microstructure of sleep EEG spindles: Possibility for dementia biomarkers? J. of Neuroscience Methods, Vol 185-1, 2009, 133-142.

3. Gorur, D. Halici, U., Aydin, H., Ongun, G., Ozgen, F., Leblebicioglu, K., Sleep Spindles Detection Using Autoregressive Modelling, Proc. of ICANN/ICONIP 2003.

4. Rechtschaffen, A, Kales, A. A manual of standardised terminology, techniques and scoring system for sleep stages of human subjects. Washington, DC: Public Health Service, U.S. Government Printing Office; 1968.

5. Steriade, M., Jones, E.G., and Llinas. Thalamic Oscillations and Signaling. Neuroscience Institute Publications. New York: John Wiley & Sons, 1990.

6. Petit, D., Gagnon, J.F., Fantini, M.L., Ferini-Strambi, L., Montplaisir, J. Sleep and quantitative EEG in neurodegenerative disorders. J Psychosom Res 2004;56:487-96.

7. Kizilkaya, A., Kayran, A. H. ARMA model parameter estimation based on the equivalent MA approach. Digital Signal Processing, Vol 16, Issue 6, 2006)

8. Wikipedia, http://pt.wikipedia.org/wiki/Sustentabilidade

9. Spinosa, M.J., Garzon, E. Sleep spindles: validated concepts and breakthroughs. J. epilepsy clin. neurophysiol. vol.13 no.4 Porto Alegre,2007

10. Reynolds III, C. F., Kupfer, D. J., Taska, L. S., et al. EEG sleep in elderly depressed, demented, and healthy subjects. Biological Psychiatry, vol. 20, no. 4, pp. 431-442, 1985