- Email: [email protected]
Analysis of EEG signals during epileptic and alcoholic states using AR modeling techniques
Disponible en ligne sur www.sciencedirect.com
ITBM-RBM 29 (2008) 44–52
Original article
Analysis of EEG signals during epileptic and alcoholic states using AR modeling techniques Analyse des signaux EEG pendant les états épileptiques et alcooliques utilisant les techniques AR O. Faust b , R.U. Acharya a,∗ , A.R. Allen b , C.M. Lin a a
Electrical and Computer Engineering Division, Ngee Ann Polytechnic, 535, Clementi Road, 599489 Singapore, Singapore b School of Engineering, University of Aberdeen, Aberdeen AB24 3UE Scotland, UK Received 19 June 2007; accepted 6 November 2007 Available online 26 December 2007
Abstract Electroencephalogram (EEG) analysis remains problematic due to limited understanding of the signal origin, which leads to the difficulty of designing evaluation methods. In spite of these shortcomings, the EEG is a valuable tool in the evaluation of some neurological disorders as well as in the evaluation of overall cerebral activity. In most studies, which use quantitative EEG analysis, the properties of measured EEG are computed by applying power spectral density (PSD) estimation for selected representative EEG samples. The sample for which the PSD is calculated is assumed to be stationary. This work deals with a comparative study of the PSD obtained from normal, epileptic and alcoholic EEG signals. The power density spectra were calculated using fast Fourier transform (FFT) by Welch’s method, auto regressive (AR) method by Yule–Walker and Burg’s method. The results are tabulated for these different classes of EEG signals. © 2007 Elsevier Masson SAS. All rights reserved. Résumé L’analyse des électroencéphalogrammes (EEG) reste problématique à cause de la compréhension limitée de l’origine du signal. Cela rend difficile la conception de méthodes d’évaluation. L’EEG est cependant un outil précieux pour l’analyse de désordres neurologiques et d’activité cérébrale totale. Dans la plupart des études qui utilisent l’EEG quantitative, les propriétés des signaux sont évaluées par la densité spectrale (PSD) d’échantillons représentatifs et supposés stationnaires. Cet article propose une étude comparative de la PSD des signaux EEG obtenus sur des sujets normaux, épileptiques et alcooliques. Les PSD ont été calculées en utilisant les transformées de Fourier rapides (FFT) obtenues soit par la méthode Welch, soit par la méthode régressive (AR) de Yule-Walker, soit par la méthode de Burg. Les résultats obtenus sont présentés selon les différentes classes de signaux EEG. © 2007 Elsevier Masson SAS. All rights reserved. Keywords: Electroencephalogram; AR model; FFT; Alcoholic; Epileptic Mots clés : EEG ; Désordres neurologiques ; Activité cérébrale
1. Introduction The electroencephalogram (EEG) is a reliable reflection of the many physiological factors modulating the brain. Despite the
∗
Corresponding author. E-mail address: [email protected] (R.U. Acharya).
1297-9562/$ – see front matter © 2007 Elsevier Masson SAS. All rights reserved. doi:10.1016/j.rbmret.2007.11.003
many applications of EEG in clinical neurophysiology [1], its visual interpretation is very subjective and does not lend itself to statistical analysis. As a result, a number of research groups have proposed methods to quantify the information content of the EEG. Among these are the Fourier transform (FFT), the wavelet transform, chaos theory, entropy and subband wavelet entropy [2–4]. In addition, EEG signal modeling is important to achieve a better understanding of the physical mechanisms generating
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
these signals and to identify the causes of EEG signal changes [5]. Pattern recognition and classification of EEG abnormalities can be achieved through the analysis of the estimated model parameters. Modeling can also be used for predicting the future neurological outcome and for data compression. Simulations based on EEG signal model can be used to demonstrate the effectiveness of a certain quantitative analysis method or EEG feature extraction. Since the early days of automatic EEG processing, representations based on a FFT have been most commonly applied. This approach is based on earlier observations that the EEG spectrum contains some characteristic waveforms that fall primarily within four frequency bands: delta (< 4 Hz), theta (4–8 Hz), alpha (8–14 Hz) and beta (14–30 Hz). Such methods have proved beneficial for various EEG characterizations, but fast FFT, suffer from large noise sensitivity. Furthermore, spectral leakage occurs if windowing is applied. Parametric power spectrum estimation methods such as auto regressive (AR) reduce the spectral leakage problems and gives better frequency resolution. Also, AR method has an advantage over FFT, that it needs shorter duration data records than FFT [6,7]. In addition, it is faster than continued wavelet transform techniques, especially in real time applications (Kiymik et al., 2001 [31]) Over the past two decades, much research has been done with the use of conventional temporal and frequency analyses measures in the detection of epileptic form activity in EEGs and reasonably good results have been obtained [6,7]. The results of different types of epilepsy from five patients were reported [8]. Their fast FFT power spectra revealed the presence of fast oscillations in two frequency bands (40–50 Hz and 80–120 Hz). This activity was found to appear at seizure onset and often on electrocorticograph (ECoG) contacts supposed to be close to the seizure focus. In the same time, it was also reported that high-frequency rhythmic activity in SEEG signals can be observed in frontal lobe epilepsy [9]. AR spectra showed that rapid activities may be observed in many sites over a wide region of the frontal lobe. Here again, these sites are supposed to anatomically belong to the epileptogenic zone. The data obtained from the analysis of the ictal EEG events compared with clinical and interictal EEG features, indicate that an asymmetric EEG pattern (mainly consisting of a rhythmic burst of fast activity) consistently preceded both symmetric and asymmetric spasms, thus suggesting a localized cortical origin of the ictal discharge giving rise to the spasms [10]. Mackenzie et al. reported results on the regional and spectral distribution of EEG rhythms in picrotoxin-induced seizures in the rat [11]. In the former study, AR analysis revealed the presence of a short rapid discharge of average frequency around 20 Hz. In the latter, fast activity was found during seizures mainly in neocortex and in several subcortical areas. Al-Nashash et al. have used a method based on neural network for modeling EEG signals [12]. Furthermore, they have implemented EEG modeling algorithms with adaptive Markov process amplitude [13]. Kannathal et al. have studied the normal, epileptic and alcoholic EEG signals using nonlinear parameters like correlation dimension, largest Lyapunov exponent, Hurst exponent, entropy etc [14]. They showed that, during epilepsy
45
and alcoholic condition, EEG signals become less random. Also, they have analyzed the EEG signals using nonlinear parameters at different mental states similar to reflexological stimulation at the feet, influence of rock music, classical music, etc [15]. They showed that the variability of the EEG signal is higher under normal condition compared to the music and reflexological stimulation. In this paper, we have analyzed normal, epileptic and alcoholic EEG signals using FFT and AR modeling techniques. In the result section we use the ANOVA test to show that the power spectral density (PSDs) of these signals are distinct. Furthermore, the ANOVA test shows that AR modeling is superior to FFT. This statement is supported by results obtained from receiver operating characteristic (ROC) tests. 2. Materials All the EEG signals used for this study were recorded with 128 channel system with 12 bit A/D resolution. The EEG data were obtained from two different sources, one source provided the data for epileptic seizure analysis and the other one provided the data for alcoholic subjects EEG analysis. EEG recordings of control and alcoholic subjects were obtained from the University of California, Irvine Knowledge Discovery in Databases (UCI KDD) archive. We have chosen 30 data sets each from control (normal) and alcoholic group for our analysis. The electrode positions were located at standard sites [16]. Zhang et al. [17] describe in detail the data collection process. Set A consisted of segments taken from surface EEG recordings that were carried out on five healthy volunteers using a standardized electrode placement scheme. Volunteers were relaxed in an awake state with eyes open. These normal and alcoholic EEG signals were recorded with a sampling frequency of 256 Hz. For epileptic data analysis, data were obtained from the EEG database available from the Bonn University with a sampling rate of 173.61 Hz. The EEG signals were taken from the correct epileptonic zone [18]. Two sets, each containing 30 single channel EEG segments of 23.6-second duration, were composed for the study. These segments were selected and cut out from continuous multichannel EEG recordings after visual inspection for artifacts, e.g., due to muscle activity or eye movements. Set E contains EEG data from five patients diagnosed with epilepsy. In this set the EEG signals were recorded during seizure activity. 3. Spectral analysis The approaches for spectrum estimation may be generally categorized into one of the two classes. The first includes the classical or nonparametric methods that begin by estimating the autocorrelation sequence rx (k) from a given set of data. The power spectrum is then estimated by Fourier transforming the estimated autocorrelation sequence. The second class includes the nonclassical or parametric approaches, which are based on using a model for the process in order to estimate the power spectrum. These techniques are briefly explained below.
46
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
3.1. Nonparametric method (FFT) by Welch’s method
in the Yule–Walker equations:
The Welch’s method [19] is one of the most popular methods to estimate the power spectrum at any given time sequence. The sequences are allowed to overlap and a data window is applied to each sequence. This will produce set of modified periodograms that are to be averaged. The data sequences xi (n) can be represented as xi (n) = x(n + iD)
n = 0, 1, 2, . . . , M − 1 and
i = 0, 1, 2, . . . , L − 1
(1)
where iD is the starting point for the ith sequence. Finally, we can form L data segments each of length 2M. The resulting modified periodogram is: 2 M−1 1 (i) (f ) = P˜ xx xi (n)w(n)e−j2πfn MU
(2)
⎡ ⎢ ⎢ ⎢ ⎢ ⎣
rxx (0)
rxx (−1)
rxx (1) .. .
rxx (0) .. .
rxx (p − 1) rxx (p − 2) ⎡ ⎤ rxx (1) ⎢ r (2) ⎥ ⎢ xx ⎥ ⎥ =⎢ ⎢ .. ⎥ ⎣ . ⎦
M−1 1 2 U= w (n) M
(3)
n=0
and w(n) is the window function. The Welch power spectrum is average of these modified periodograms and is given by: W Pxx (f ) =
L−1
1 ˜ (i) Pxx (f ) L
(4)
i=0
3.2. Parametric method The nonparametric method (e.g. FFT method) is not optimal because this method suffers from spectral leakage effects due to windowing. Spectral leakage might mask weak signals present in the data. A parametric (model based) power spectrum estimation method avoids the problem of spectral leakage and provides better frequency resolution than nonparametric or classical methods. The parametric methods assume the signal to be a stationary random process. This process is modeled as the output of a filter with white noise input. These methods estimate the PSD by first estimating the parameters (coefficients) of the linear system (filter) that hypothetically generates the signal [30]. In AR model method, there are poles present in the z-domain and all the zeros are at the origin. In the present study, we have estimated the AR models using Yule–Walker and Burg’s method. A brief description of these methods are given below. 3.2.1. Yule–Walker method The Yule–Walker AR method of spectral estimation computes the AR parameters by forming a biased estimate of the signal’s autocorrelation function and solving the least squares minimization of the forward prediction error [20]. This results
(5)
rxx (p) where rxx is a biased form of the autocorrelation function which ensures that autocorrelation matrix above is positive definite. The biased form of the autocorrelation estimate is calculated as follows:
n=0
where U is the normalization factor for the power in the window function and is selected as:
⎤ ⎡ ⎤ · · · rxx (−p + 1) aˆ p (1) ⎢ ⎥ · · · rxx (−p + 2) ⎥ ⎥ ⎢ aˆ P (2) ⎥ ⎥×⎢ ⎥ .. ⎥ ⎢ .. ⎥ .. . ⎦ ⎣ . ⎦ . aˆ p (p) ··· rxx (0)
rxx (m) =
N−m−1 1 ∗ x (n)x(n + m) N
m ≥ 0.
(6)
N=0
The AR coefficients (ˆap ) can be obtained by solving the above set of p + 1 linear equations (for instance, by using the fast Levinson–Durbin algorithm). The corresponding PSD estimate is: 2 σˆ wp YW (f ) = Pxx
−j2πfk 2 1 + p k=1 ap (k)e
(7)
where σˆ wp is the estimated minimum mean square error for the pth-order predictor calculated as follows: 2 σˆ wp
ˆ pf = rxx (0) =E
p
1 − |ˆak (k)|2
(8)
k=1
3.2.2. Burg’s method The Burg’s method for AR spectral estimation is based on minimizing the forward and backward prediction errors, while satisfying the Levinson–Durbin recursion [21,22]. In contrast to other AR estimation methods, the Burg’s method avoids calculating the autocorrelation function and instead estimates the reflection coefficients directly. The primary advantages of the Burg’s method are resolving closely spaced sinusoids in signals with low noise levels and estimating short data records, in which case the AR, PSD estimates are very close to the true values. In addition, the Burg’s method ensures a stable AR model and is computationally efficient. The accuracy of the Burg’s method is lower for high-order models, long data records and high signalto-noise ratios (which can cause line splitting or the generation of extraneous peaks in the spectrum estimate). The spectral density estimate computed by the Burg’s method is also susceptible to frequency shifts (relative to the true frequency) resulting from the initial phase of noisy sinusoidal signals. This effect is magnified when analyzing short data sequences. Burg’s method differs BU (f ), is from the Yule–Walker method in the way the PSD, Pxx
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
obtained as shown in the following equation:
BU (f ) Pxx
៝p E
= 2 p ap (k)e−j2πfk 1 +
(9)
3.2.3. Akaike information criterion For both Burg and Yule–Walker, the model order was chosen as the one that minimizes the Akaike information criterion (AIC) figure of merit [23,24]: AIC(p) = N · In(λˆ 2 ) + 2p
k=1
where the ap (k) are again the estimates of the AR parameters obtained from the Levinson–Durbin recursion, and ៝ p are reflection coefficients in an equivalent lattice strucE ture, which are chosen to obtain the total least square error.
47
(10)
where N is the number of data samples and λˆ 2 is the estimated white noise variance. To reduce computational costs, we assumed as optimal the value of p that fulfilled the AIC criterion in the first two epochs. One of the most important aspects of the AR method is the selection of the order p. Much work has been done by various researchers on this problem and many experimental results have
Fig. 1. PSD estimation using Welch’s method.
48
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
been given in literature [25]. In this work the order of the AR model is taken as: p = 20. 3.2.4. Peak amplitude detection The maxima of the individual PSDs, shown in Figs. 1–3, were detected on a sample by sample basis. A particular sample was detected as (local) maxima, if its value was larger than the value of both previous and next sample. The detection function builds up a set of ordered pairs. Each pair was composed out of the value and the position of the (local) maxima. In a second step the pair with the largest value, that is the maximum, is removed from the set. The ordered pair with the largest value in the remaining
set is the first local maxima. After the first local maximum is removed from the set, the pair with largest value is the second maxima. In Figs. 1–3 the global, first and second maxima are indicated with circle, cross and triangle, respectively. 4. Results The PSD was calculated for each of the signals using Welch’s, Yule–Walker and Burg’s method. The first plot in Fig. 1 shows the PSDs obtained form the first 10 signals in the control data. The second plot (in Fig. 1) shows the PSDs from the first 10 signals of the alcohol data set and the last plot shows the PSDs
Fig. 2. PSD estimation using Yule–Walker method.
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
49
Fig. 3. PSD estimation using Burg’s method.
from the first 10 signals in epileptic data set. In this figure, all PSDs where estimated using Welch’s method. For each of the PSD functions, a circle marks, a cross and a triangle marks the first, second and third prominent frequencies, respectively. The first, second, and third prominent frequencies represent the first, second, and third maxima of the PSD function, respectively. Figs. 2 and 3 show the PSDs using Yule–Walker and Burg’s method (AR model). It can be observed from these figures, that the prominent peaks are at similar frequencies for the same PSD estimation method, but the peaks are at distinct frequencies in different methods. In the next section we quantify this statement with the ANOVA test.
4.1. ANOVA test Tables 1–3 show the first three peak powers, their corresponding frequencies and their corresponding ratios for the three types of signals using Welch’s, Yule–Walker and Burg’s method, respectively. The first parameter, P1, represents mean and variance of the peak power. The second parameter, F1, symbolizes mean and variance of the frequency where P1 is located. The third parameter, P1/F1, corresponds to the peak power divided by its corresponding frequency. Similarly, P2 represents mean and variance of second peak power, F2 is the mean and variance of the corresponding frequency.
50
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
Table 1 ANOVA results for Welch’s method Parameter
Control
Alcohol
P1 F1 P1/F1 P2 F2 P2/F2 P3 F3 P3/F3
1.717e + 01 2.800e + 00 1.575e + 01 4.441e + 00 8.300e + 00 6.403e-01 3.867e + 00 1.210e + 01 1.039e-01
± ± ± ± ± ± ± ± ±
4.107e + 02 7.476e + 00 4.477e + 02 3.666e + 01 1.263e + 01 6.009e-01 1.003e + 02 1.423e + 01 1.566e + 00
Epileptic
1.348e + 01 3.300e + 00 1.207e + 01 2.534e + 00 8.167e + 00 3.934e-01 1.930e + 00 1.237e + 01 1.938e-01
± ± ± ± ± ± ± ± ±
3.388e + 02 1.008e + 01 3.331e + 02 1.003e + 01 1.587e + 01 1.992e-01 1.182e + 01 1.707e + 01 1.213e-01
P-value
9.912e + 03 3.533e + 00 2.413e + 03 1.128e + 04 9.067e + 00 1.402e + 03 5.825e + 03 1.403e + 01 4.880e + 02
± ± ± ± ± ± ± ± ±
1.792e + 08 3.430e + 00 9.191e + 06 .859e + 08 8.133e + 00 3.367e + 06 4.890e + 07 1.121e + 01 4.583e + 05
9.068e-07 5.500e-01 1.648e-07 4.979e-08 5.623e-01 4.096e-07 4.163e-08 1.039e-01 1.664e-06
Table 2 ANOVA results for Yule–Walker method P/F
Control
P1 F1 P1/F1 P2 F2 P2/F2 P3 F3 P3/F3
2.788e + 01 1.566e + 00 5.524e + 01 4.416e + 00 1.738e + 01 4.455e-02 5.599e-01 3.492e + 01 2.170e-02
Alcohol ± ± ± ± ± ± ± ± ±
5.732e + 02 1.820e + 01 2.346e + 03 6.742e + 01 7.265e + 01 1.017e + 00 4.065e-01 1.841e + 02 6.466e-04
2.780e + 01 3.050e + 00 5.500e + 01 2.120e + 00 2.058e + 01 2.047e-01 3.999e-01 4.168e + 01 1.604e-02
Epileptic ± ± ± ± ± ± ± ± ±
The ratio (P2/F2) represents the ratio of P2 divided by F2. The parameters P3, F3, and P3/F3 indicate third peak power, corresponding frequency and the third peak power divided by its corresponding frequency. In all the tables, the parameters are listed in the first column. The second column contains the parameter values for the control group. The second and third columns hold the parameter values for alcoholic and epileptic groups, respectively. The last column holds the ‘p-value’ which was evaluated using ANOVA test. Our results show very low ‘p-value’, indicating that the results are clinically significant. This test helps to determine whether or not data groups have distinguishable characteristics, even if they show the same mean value. Table 1 shows the PSD results of the three classes using FFT: Welch’s method. In this method, the P1, P2, P3 and F1, F2, F3 are close to each in all the three classes. Tables 2 and 3 show the PSD results of the three classes using AR method by Yule–Walker
1.296e + 03 2.947e + 01 5.251e + 03 1.594e + 01 1.456e + 02 1.804e-01 6.969e-01 3.599e + 02 1.454e-03
2.231e + 04 6.383e + 00 4.040e + 03 5.358e + 03 2.853e + 01 2.772e + 02 2.047e + 02 7.823e + 01 8.257e + 00
P-value ± ± ± ± ± ± ± ± ±
5.504e + 08 8.408e + 00 1.977e + 07 4.910e + 07 7.016e + 02 1.116e + 05 4.580e + 05 9.604e + 02 7.404e + 02
7.251e-10 1.492e-04 4.759e-09 4.046e-07 4.455e-02 4.758e-08 7.074e-02 2.900e-11 6.951e-02
and Burg’s method, respectively. In this method, the P1, P2, P3 and their corresponding frequencies are relatively distinct from each other in all the three classes. However, Burg’s method is more effective as compared to the Yule–Walker method because the ‘p-values’ shown in Table 3 are very low. This can also be evaluated using ROC curve. 4.2. Receiver operating characteristic test Receiver operating characteristic (ROC) test was performed using the parameters F1, F2, and F3. ROC is a one-class classification method, which classifies a target class within a dataset. The dataset consists of target (desired) values or vectors and outliers. To show that the frequency locations of the prominent frequencies are distinct, the alcoholic and epileptic data were treated as target and the control data was treated as outlier. Figs. 4–6 show the results of the ROC test based on PSDs
Table 3 ANOVA results for Burg’s method P/F P1 F1 P1/F1 P2 F2 P2/F2 P3 F3 P3/F3
Control 2.8141 + 01 1.567e + 00 5.576e + 01 4.477e + 00 1.740e + 01 4.855e-01 5.601e-01 3.497e + 01 2.167e-02
Alcohol ± ± ± ± ± ± ± ± ±
5.732e + 02 1.820e + 01 2.415e + 03 7.063e + 01 7.313e + 01 1.068e + 00 4.046e-01 1.848e + 02 6.441e-04
2.494e + 01 4.200e + 00 4.914e + 01 2.124e + 00 2.232e + 01 2.048e-01 4.046e-01 4.227e + 01 1.623e-02
Epileptic ± ± ± ± ± ± ± ± ±
1.309e + 03 4.641e + 01 5.307e + 03 1.736e + 01 1.951e + 02 1.941e-01 7.086e-01 3.857e + 02 1.476e-03
2.264e + 04 6.350e + 00 4.157e + 03 5.497e + 03 2.860e + 01 2.898e + 02 2.617e + 02 7.647e + 01 1.066e + 01
P-value ± ± ± ± ± ± ± ± ±
5.755e + 08 8.761e + 00 2.027e + 07 4.972e + 07 7.553e + 02 1.154e + 05 5.639 + 05 1.070e + 03 9.224e + 02
9.310e-10 1.466e-03 2.782e-09 2.476e-07 6.812e-02 2.137e-08 3.058e-02 7.889e-10 2.914e-02
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
Fig. 4. ROC results using FFT (Welch’s) method.
obtained with the Welch’s, Yule–Walker and Burg’s method, respectively. The area under the ROC curve indicates whether the target and the outliers are distinct. It can be observed that the ROC test for Welch’s method, shown in Fig. 4, yields low values for the area under the curves. In general, this small area indicates that the target and the outliers are similar. So, it is difficult to distinguish the control data from either alcoholic or epileptic data. It is because, in Welch’s method (Fig. 1), first three spikes are mistaken to be the prominent frequencies. The parametric methods perform much better in the ROC test, as shown in Figs. 5 and 6. The area under the ROC curve for the epileptic group is the same for Burg’s and Yule–Walker method. But for the alcohol dataset, Burg’s method performs slightly better than the Yule–Walker method.
Fig. 5. ROC results using Yule–Walker method.
51
Fig. 6. ROC results using Burg’s method.
5. Discussion In this discussion section, we review some of the work which preceded the current study. The discussion provides further insight in the subject, which is necessary to understand relevance and distinctness of our work. In other studies, dynamical properties like correlation dimension, largest Lyapunov exponent, Hurst exponent, entropy, were evaluated for alcoholic, normal, and epileptic EEG signals [14]. The results show that the alcoholic EEG is more complex than the epileptic seizure signal. And they also showed that the normal EEG was more complex compared to the other two types, indicating higher neuronal process in the brain during the normal active state. In the normal stage, the cortex is more active and more neurons are available for processing [26]. And in epileptic stage, the cortex becomes inactive and EEG becomes less random. The neurons in the cerebral hemispheres during the seizure misfire and create abnormal electrical activity. People with epilepsy have seizures that are a bit like an electrical brainstorm. Hence, the number of neurons available for processing the information reduces during the seizures. Authors have showed that the variability of epileptic activity was less as compared to that of nonepileptic activity [27,28]. It was supported by the reduced dimensionality of epileptic seizures as compared to nonepileptic EEG signals. This concept finds support in the observations that neuronal hypersynchrony underlies seizures: a phenomenon during which the number of independent variables required to describe the system was smaller [29]. The algorithm that we used is sensitive towards noise because only three samples are considered for the peak amplitude detection. Therefore, this algorithm works best for signals that contain mostly low frequencies.
52
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
The PSD estimation using Welch’s method yields signals with high frequency content (presence of spikes). Each of these spikes constitute for local maxima. One way to reduce the number of local maxima is to smooth the signal. However, this smoothing may alter the information content in the signal. In this particular case, Welch’s method is good in discriminating closely spaced frequency peaks. The smoothing function might combine two distinct peaks to one. Parametric PSD estimation methods, like Burg’s and Yule–Walker, yield a smooth spectrum. Therefore, the simple maxima detection algorithm works well with these signals. Thus, the nine features extracted from these PSD plots are distinct for control, epileptic, alcoholic groups. These nine features can constitute a feature vector which is fed into a classifier. Such a classifier maps individual feature vectors to one of three different classes. 6. Conclusion The EEG signal can be used to diagnose different mental conditions like active, epilepsy and alcoholic. The changes in EEG signals might be quite prominent, as in the case of an epileptic patient or more hidden (complex), as in the case of alcoholic subject. In the time domain, only a trained eye can detect these different states. This work shows that, these changes in the states are also visible in the spectral domain. We used three different ways (Welch’s, Yule–Walker and Burg’s) to represent the power distribution in the frequency domain and compared their performances in Tables 1–3. The nine parameters evaluated in the results show that, the parametric methods are superior compared to the nonparametric methods. We have distinguished between alcoholic, epileptic and control EEG signals using ROC method. Based on the ROC test we come to the conclusion that Burg’s method yields the most distinguishable parameters. This study can be seen as the groundwork for more sophisticated EEG signal classification. This signal classification is another step towards an automated system that is able to diagnose different mental conditions based on EEG signals. Such a system would significantly improve clinical workflows, because it frees up trained personal from routine jobs. References [1] James C., Lowe D. Using dynamical embedding to isolate seizure components in the ictal EEG, 1st International Conference on Advances in Medical Signal and Information Processing 2000, 476, 2000, 158–165. [2] Geocadin R, Ghodadra R, Kimura T, Lie H, Sherman D, Hanley D, et al. A novel quantitative EEG injury measure of global cerebral ischemia. Clin Neurophysiol 2000;111:1779–87. [3] Bezerianos A., Tong S., Malhorta A., Thakor N. Information measures of brain dynamics, Nonlinear signal and image processing conference, Baltimore, MD, 2001. [4] Rosso O, Blanco S, Yordanova J, Kolev V, Figliola A, Schurmann M, et al. Wavelet entropy: a new tool for analysis of short duration brain electrical signals. J Neurosci Methods 2001;105:65–75.
[5] Bai O, Nakamura M, Ikeda A, Shibasaki H. Nonlinear Markov process amplitude EEG model for nonlinear coupling interaction of spontaneous EEG. IEEE Trans Biomed Eng 2000;47:1141–6. [6] Guler I, Kiymik MK, Akin M, Alkan A. AR spectral analysis of EEG signals by using maximum likelihood estimation. Med Biol 2001;31:441–50. [7] Akin M, Kiymik MK. Application of periodogram and AR spectral analysis to EEG signals. J Med Syst 2000;24(4):247–56. [8] Fisher RS, Webber WR, Lesser RP, Arroyo S, Uematsu S. High-frequency EEG activity at the start of seizures. J Clin Neurophysiol 1992;9(3):441–8. [9] Allen PJ, Fish DR, Smith SJ. Very high-frequency rhythmic activity during SEEG suppression in frontal lobe epilepsy. Electroencephalogr Clin Neurophysiol 1992;82(2):155–62. [10] Panzica F, Franceschetti S, Binelli S, Canafoglia L, Granata T, Avanzini G. Spectral properties of EEG fast activity ictal discharges associated with infantile spasms. Clin Neurophysiol 1999;110(4):593–603. [11] Mackenzie L, Medvedev A, Hiscock JJ, Pope KJ, Willoughby JO. Picrotoxin-induced generalised convulsive seizure in rat: changes in regional distribution and frequency of the power of electroencephalogram rhythms. Clin Neurophysiol 2002;113(4):586–96. [12] Al-Nashash HA, Zalzala AMS, Thakor NV. A neural networks approach to EEG signals modeling. Med Biol 2003;3:1–6. [13] Al-Nashash Ha, Al-Assaf Y, Paul J, Thakor N. EEG signal modeling using adaptive Markov process amplitude. IEEE Trans Biomed Eng 2004;51(5):744–51. [14] Kannathal N, Acharya RU, Lim Choo M, Quek W, Mohamed H, Sadasivan PK. Characterization of EEG: a comparative study. Comput Methods Programs Biomed 2005;80(1):17–23. [15] Kannathal N, Rajendra Acharya U, Fadilah A, Thelma T, Sadasivan PK. Nonlinear analysis of EEG signals at different mental states. Biomedical 2004 16 March;3(7). [16] American Electroencephalographic Association 1990. Standard electrode position nomenclature,
ITBM-RBM 29 (2008) 44–52
Original article
Analysis of EEG signals during epileptic and alcoholic states using AR modeling techniques Analyse des signaux EEG pendant les états épileptiques et alcooliques utilisant les techniques AR O. Faust b , R.U. Acharya a,∗ , A.R. Allen b , C.M. Lin a a
Electrical and Computer Engineering Division, Ngee Ann Polytechnic, 535, Clementi Road, 599489 Singapore, Singapore b School of Engineering, University of Aberdeen, Aberdeen AB24 3UE Scotland, UK Received 19 June 2007; accepted 6 November 2007 Available online 26 December 2007
Abstract Electroencephalogram (EEG) analysis remains problematic due to limited understanding of the signal origin, which leads to the difficulty of designing evaluation methods. In spite of these shortcomings, the EEG is a valuable tool in the evaluation of some neurological disorders as well as in the evaluation of overall cerebral activity. In most studies, which use quantitative EEG analysis, the properties of measured EEG are computed by applying power spectral density (PSD) estimation for selected representative EEG samples. The sample for which the PSD is calculated is assumed to be stationary. This work deals with a comparative study of the PSD obtained from normal, epileptic and alcoholic EEG signals. The power density spectra were calculated using fast Fourier transform (FFT) by Welch’s method, auto regressive (AR) method by Yule–Walker and Burg’s method. The results are tabulated for these different classes of EEG signals. © 2007 Elsevier Masson SAS. All rights reserved. Résumé L’analyse des électroencéphalogrammes (EEG) reste problématique à cause de la compréhension limitée de l’origine du signal. Cela rend difficile la conception de méthodes d’évaluation. L’EEG est cependant un outil précieux pour l’analyse de désordres neurologiques et d’activité cérébrale totale. Dans la plupart des études qui utilisent l’EEG quantitative, les propriétés des signaux sont évaluées par la densité spectrale (PSD) d’échantillons représentatifs et supposés stationnaires. Cet article propose une étude comparative de la PSD des signaux EEG obtenus sur des sujets normaux, épileptiques et alcooliques. Les PSD ont été calculées en utilisant les transformées de Fourier rapides (FFT) obtenues soit par la méthode Welch, soit par la méthode régressive (AR) de Yule-Walker, soit par la méthode de Burg. Les résultats obtenus sont présentés selon les différentes classes de signaux EEG. © 2007 Elsevier Masson SAS. All rights reserved. Keywords: Electroencephalogram; AR model; FFT; Alcoholic; Epileptic Mots clés : EEG ; Désordres neurologiques ; Activité cérébrale
1. Introduction The electroencephalogram (EEG) is a reliable reflection of the many physiological factors modulating the brain. Despite the
∗
Corresponding author. E-mail address: [email protected] (R.U. Acharya).
1297-9562/$ – see front matter © 2007 Elsevier Masson SAS. All rights reserved. doi:10.1016/j.rbmret.2007.11.003
many applications of EEG in clinical neurophysiology [1], its visual interpretation is very subjective and does not lend itself to statistical analysis. As a result, a number of research groups have proposed methods to quantify the information content of the EEG. Among these are the Fourier transform (FFT), the wavelet transform, chaos theory, entropy and subband wavelet entropy [2–4]. In addition, EEG signal modeling is important to achieve a better understanding of the physical mechanisms generating
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
these signals and to identify the causes of EEG signal changes [5]. Pattern recognition and classification of EEG abnormalities can be achieved through the analysis of the estimated model parameters. Modeling can also be used for predicting the future neurological outcome and for data compression. Simulations based on EEG signal model can be used to demonstrate the effectiveness of a certain quantitative analysis method or EEG feature extraction. Since the early days of automatic EEG processing, representations based on a FFT have been most commonly applied. This approach is based on earlier observations that the EEG spectrum contains some characteristic waveforms that fall primarily within four frequency bands: delta (< 4 Hz), theta (4–8 Hz), alpha (8–14 Hz) and beta (14–30 Hz). Such methods have proved beneficial for various EEG characterizations, but fast FFT, suffer from large noise sensitivity. Furthermore, spectral leakage occurs if windowing is applied. Parametric power spectrum estimation methods such as auto regressive (AR) reduce the spectral leakage problems and gives better frequency resolution. Also, AR method has an advantage over FFT, that it needs shorter duration data records than FFT [6,7]. In addition, it is faster than continued wavelet transform techniques, especially in real time applications (Kiymik et al., 2001 [31]) Over the past two decades, much research has been done with the use of conventional temporal and frequency analyses measures in the detection of epileptic form activity in EEGs and reasonably good results have been obtained [6,7]. The results of different types of epilepsy from five patients were reported [8]. Their fast FFT power spectra revealed the presence of fast oscillations in two frequency bands (40–50 Hz and 80–120 Hz). This activity was found to appear at seizure onset and often on electrocorticograph (ECoG) contacts supposed to be close to the seizure focus. In the same time, it was also reported that high-frequency rhythmic activity in SEEG signals can be observed in frontal lobe epilepsy [9]. AR spectra showed that rapid activities may be observed in many sites over a wide region of the frontal lobe. Here again, these sites are supposed to anatomically belong to the epileptogenic zone. The data obtained from the analysis of the ictal EEG events compared with clinical and interictal EEG features, indicate that an asymmetric EEG pattern (mainly consisting of a rhythmic burst of fast activity) consistently preceded both symmetric and asymmetric spasms, thus suggesting a localized cortical origin of the ictal discharge giving rise to the spasms [10]. Mackenzie et al. reported results on the regional and spectral distribution of EEG rhythms in picrotoxin-induced seizures in the rat [11]. In the former study, AR analysis revealed the presence of a short rapid discharge of average frequency around 20 Hz. In the latter, fast activity was found during seizures mainly in neocortex and in several subcortical areas. Al-Nashash et al. have used a method based on neural network for modeling EEG signals [12]. Furthermore, they have implemented EEG modeling algorithms with adaptive Markov process amplitude [13]. Kannathal et al. have studied the normal, epileptic and alcoholic EEG signals using nonlinear parameters like correlation dimension, largest Lyapunov exponent, Hurst exponent, entropy etc [14]. They showed that, during epilepsy
45
and alcoholic condition, EEG signals become less random. Also, they have analyzed the EEG signals using nonlinear parameters at different mental states similar to reflexological stimulation at the feet, influence of rock music, classical music, etc [15]. They showed that the variability of the EEG signal is higher under normal condition compared to the music and reflexological stimulation. In this paper, we have analyzed normal, epileptic and alcoholic EEG signals using FFT and AR modeling techniques. In the result section we use the ANOVA test to show that the power spectral density (PSDs) of these signals are distinct. Furthermore, the ANOVA test shows that AR modeling is superior to FFT. This statement is supported by results obtained from receiver operating characteristic (ROC) tests. 2. Materials All the EEG signals used for this study were recorded with 128 channel system with 12 bit A/D resolution. The EEG data were obtained from two different sources, one source provided the data for epileptic seizure analysis and the other one provided the data for alcoholic subjects EEG analysis. EEG recordings of control and alcoholic subjects were obtained from the University of California, Irvine Knowledge Discovery in Databases (UCI KDD) archive. We have chosen 30 data sets each from control (normal) and alcoholic group for our analysis. The electrode positions were located at standard sites [16]. Zhang et al. [17] describe in detail the data collection process. Set A consisted of segments taken from surface EEG recordings that were carried out on five healthy volunteers using a standardized electrode placement scheme. Volunteers were relaxed in an awake state with eyes open. These normal and alcoholic EEG signals were recorded with a sampling frequency of 256 Hz. For epileptic data analysis, data were obtained from the EEG database available from the Bonn University with a sampling rate of 173.61 Hz. The EEG signals were taken from the correct epileptonic zone [18]. Two sets, each containing 30 single channel EEG segments of 23.6-second duration, were composed for the study. These segments were selected and cut out from continuous multichannel EEG recordings after visual inspection for artifacts, e.g., due to muscle activity or eye movements. Set E contains EEG data from five patients diagnosed with epilepsy. In this set the EEG signals were recorded during seizure activity. 3. Spectral analysis The approaches for spectrum estimation may be generally categorized into one of the two classes. The first includes the classical or nonparametric methods that begin by estimating the autocorrelation sequence rx (k) from a given set of data. The power spectrum is then estimated by Fourier transforming the estimated autocorrelation sequence. The second class includes the nonclassical or parametric approaches, which are based on using a model for the process in order to estimate the power spectrum. These techniques are briefly explained below.
46
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
3.1. Nonparametric method (FFT) by Welch’s method
in the Yule–Walker equations:
The Welch’s method [19] is one of the most popular methods to estimate the power spectrum at any given time sequence. The sequences are allowed to overlap and a data window is applied to each sequence. This will produce set of modified periodograms that are to be averaged. The data sequences xi (n) can be represented as xi (n) = x(n + iD)
n = 0, 1, 2, . . . , M − 1 and
i = 0, 1, 2, . . . , L − 1
(1)
where iD is the starting point for the ith sequence. Finally, we can form L data segments each of length 2M. The resulting modified periodogram is: 2 M−1 1 (i) (f ) = P˜ xx xi (n)w(n)e−j2πfn MU
(2)
⎡ ⎢ ⎢ ⎢ ⎢ ⎣
rxx (0)
rxx (−1)
rxx (1) .. .
rxx (0) .. .
rxx (p − 1) rxx (p − 2) ⎡ ⎤ rxx (1) ⎢ r (2) ⎥ ⎢ xx ⎥ ⎥ =⎢ ⎢ .. ⎥ ⎣ . ⎦
M−1 1 2 U= w (n) M
(3)
n=0
and w(n) is the window function. The Welch power spectrum is average of these modified periodograms and is given by: W Pxx (f ) =
L−1
1 ˜ (i) Pxx (f ) L
(4)
i=0
3.2. Parametric method The nonparametric method (e.g. FFT method) is not optimal because this method suffers from spectral leakage effects due to windowing. Spectral leakage might mask weak signals present in the data. A parametric (model based) power spectrum estimation method avoids the problem of spectral leakage and provides better frequency resolution than nonparametric or classical methods. The parametric methods assume the signal to be a stationary random process. This process is modeled as the output of a filter with white noise input. These methods estimate the PSD by first estimating the parameters (coefficients) of the linear system (filter) that hypothetically generates the signal [30]. In AR model method, there are poles present in the z-domain and all the zeros are at the origin. In the present study, we have estimated the AR models using Yule–Walker and Burg’s method. A brief description of these methods are given below. 3.2.1. Yule–Walker method The Yule–Walker AR method of spectral estimation computes the AR parameters by forming a biased estimate of the signal’s autocorrelation function and solving the least squares minimization of the forward prediction error [20]. This results
(5)
rxx (p) where rxx is a biased form of the autocorrelation function which ensures that autocorrelation matrix above is positive definite. The biased form of the autocorrelation estimate is calculated as follows:
n=0
where U is the normalization factor for the power in the window function and is selected as:
⎤ ⎡ ⎤ · · · rxx (−p + 1) aˆ p (1) ⎢ ⎥ · · · rxx (−p + 2) ⎥ ⎥ ⎢ aˆ P (2) ⎥ ⎥×⎢ ⎥ .. ⎥ ⎢ .. ⎥ .. . ⎦ ⎣ . ⎦ . aˆ p (p) ··· rxx (0)
rxx (m) =
N−m−1 1 ∗ x (n)x(n + m) N
m ≥ 0.
(6)
N=0
The AR coefficients (ˆap ) can be obtained by solving the above set of p + 1 linear equations (for instance, by using the fast Levinson–Durbin algorithm). The corresponding PSD estimate is: 2 σˆ wp YW (f ) = Pxx
−j2πfk 2 1 + p k=1 ap (k)e
(7)
where σˆ wp is the estimated minimum mean square error for the pth-order predictor calculated as follows: 2 σˆ wp
ˆ pf = rxx (0) =E
p
1 − |ˆak (k)|2
(8)
k=1
3.2.2. Burg’s method The Burg’s method for AR spectral estimation is based on minimizing the forward and backward prediction errors, while satisfying the Levinson–Durbin recursion [21,22]. In contrast to other AR estimation methods, the Burg’s method avoids calculating the autocorrelation function and instead estimates the reflection coefficients directly. The primary advantages of the Burg’s method are resolving closely spaced sinusoids in signals with low noise levels and estimating short data records, in which case the AR, PSD estimates are very close to the true values. In addition, the Burg’s method ensures a stable AR model and is computationally efficient. The accuracy of the Burg’s method is lower for high-order models, long data records and high signalto-noise ratios (which can cause line splitting or the generation of extraneous peaks in the spectrum estimate). The spectral density estimate computed by the Burg’s method is also susceptible to frequency shifts (relative to the true frequency) resulting from the initial phase of noisy sinusoidal signals. This effect is magnified when analyzing short data sequences. Burg’s method differs BU (f ), is from the Yule–Walker method in the way the PSD, Pxx
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
obtained as shown in the following equation:
BU (f ) Pxx
៝p E
= 2 p ap (k)e−j2πfk 1 +
(9)
3.2.3. Akaike information criterion For both Burg and Yule–Walker, the model order was chosen as the one that minimizes the Akaike information criterion (AIC) figure of merit [23,24]: AIC(p) = N · In(λˆ 2 ) + 2p
k=1
where the ap (k) are again the estimates of the AR parameters obtained from the Levinson–Durbin recursion, and ៝ p are reflection coefficients in an equivalent lattice strucE ture, which are chosen to obtain the total least square error.
47
(10)
where N is the number of data samples and λˆ 2 is the estimated white noise variance. To reduce computational costs, we assumed as optimal the value of p that fulfilled the AIC criterion in the first two epochs. One of the most important aspects of the AR method is the selection of the order p. Much work has been done by various researchers on this problem and many experimental results have
Fig. 1. PSD estimation using Welch’s method.
48
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
been given in literature [25]. In this work the order of the AR model is taken as: p = 20. 3.2.4. Peak amplitude detection The maxima of the individual PSDs, shown in Figs. 1–3, were detected on a sample by sample basis. A particular sample was detected as (local) maxima, if its value was larger than the value of both previous and next sample. The detection function builds up a set of ordered pairs. Each pair was composed out of the value and the position of the (local) maxima. In a second step the pair with the largest value, that is the maximum, is removed from the set. The ordered pair with the largest value in the remaining
set is the first local maxima. After the first local maximum is removed from the set, the pair with largest value is the second maxima. In Figs. 1–3 the global, first and second maxima are indicated with circle, cross and triangle, respectively. 4. Results The PSD was calculated for each of the signals using Welch’s, Yule–Walker and Burg’s method. The first plot in Fig. 1 shows the PSDs obtained form the first 10 signals in the control data. The second plot (in Fig. 1) shows the PSDs from the first 10 signals of the alcohol data set and the last plot shows the PSDs
Fig. 2. PSD estimation using Yule–Walker method.
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
49
Fig. 3. PSD estimation using Burg’s method.
from the first 10 signals in epileptic data set. In this figure, all PSDs where estimated using Welch’s method. For each of the PSD functions, a circle marks, a cross and a triangle marks the first, second and third prominent frequencies, respectively. The first, second, and third prominent frequencies represent the first, second, and third maxima of the PSD function, respectively. Figs. 2 and 3 show the PSDs using Yule–Walker and Burg’s method (AR model). It can be observed from these figures, that the prominent peaks are at similar frequencies for the same PSD estimation method, but the peaks are at distinct frequencies in different methods. In the next section we quantify this statement with the ANOVA test.
4.1. ANOVA test Tables 1–3 show the first three peak powers, their corresponding frequencies and their corresponding ratios for the three types of signals using Welch’s, Yule–Walker and Burg’s method, respectively. The first parameter, P1, represents mean and variance of the peak power. The second parameter, F1, symbolizes mean and variance of the frequency where P1 is located. The third parameter, P1/F1, corresponds to the peak power divided by its corresponding frequency. Similarly, P2 represents mean and variance of second peak power, F2 is the mean and variance of the corresponding frequency.
50
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
Table 1 ANOVA results for Welch’s method Parameter
Control
Alcohol
P1 F1 P1/F1 P2 F2 P2/F2 P3 F3 P3/F3
1.717e + 01 2.800e + 00 1.575e + 01 4.441e + 00 8.300e + 00 6.403e-01 3.867e + 00 1.210e + 01 1.039e-01
± ± ± ± ± ± ± ± ±
4.107e + 02 7.476e + 00 4.477e + 02 3.666e + 01 1.263e + 01 6.009e-01 1.003e + 02 1.423e + 01 1.566e + 00
Epileptic
1.348e + 01 3.300e + 00 1.207e + 01 2.534e + 00 8.167e + 00 3.934e-01 1.930e + 00 1.237e + 01 1.938e-01
± ± ± ± ± ± ± ± ±
3.388e + 02 1.008e + 01 3.331e + 02 1.003e + 01 1.587e + 01 1.992e-01 1.182e + 01 1.707e + 01 1.213e-01
P-value
9.912e + 03 3.533e + 00 2.413e + 03 1.128e + 04 9.067e + 00 1.402e + 03 5.825e + 03 1.403e + 01 4.880e + 02
± ± ± ± ± ± ± ± ±
1.792e + 08 3.430e + 00 9.191e + 06 .859e + 08 8.133e + 00 3.367e + 06 4.890e + 07 1.121e + 01 4.583e + 05
9.068e-07 5.500e-01 1.648e-07 4.979e-08 5.623e-01 4.096e-07 4.163e-08 1.039e-01 1.664e-06
Table 2 ANOVA results for Yule–Walker method P/F
Control
P1 F1 P1/F1 P2 F2 P2/F2 P3 F3 P3/F3
2.788e + 01 1.566e + 00 5.524e + 01 4.416e + 00 1.738e + 01 4.455e-02 5.599e-01 3.492e + 01 2.170e-02
Alcohol ± ± ± ± ± ± ± ± ±
5.732e + 02 1.820e + 01 2.346e + 03 6.742e + 01 7.265e + 01 1.017e + 00 4.065e-01 1.841e + 02 6.466e-04
2.780e + 01 3.050e + 00 5.500e + 01 2.120e + 00 2.058e + 01 2.047e-01 3.999e-01 4.168e + 01 1.604e-02
Epileptic ± ± ± ± ± ± ± ± ±
The ratio (P2/F2) represents the ratio of P2 divided by F2. The parameters P3, F3, and P3/F3 indicate third peak power, corresponding frequency and the third peak power divided by its corresponding frequency. In all the tables, the parameters are listed in the first column. The second column contains the parameter values for the control group. The second and third columns hold the parameter values for alcoholic and epileptic groups, respectively. The last column holds the ‘p-value’ which was evaluated using ANOVA test. Our results show very low ‘p-value’, indicating that the results are clinically significant. This test helps to determine whether or not data groups have distinguishable characteristics, even if they show the same mean value. Table 1 shows the PSD results of the three classes using FFT: Welch’s method. In this method, the P1, P2, P3 and F1, F2, F3 are close to each in all the three classes. Tables 2 and 3 show the PSD results of the three classes using AR method by Yule–Walker
1.296e + 03 2.947e + 01 5.251e + 03 1.594e + 01 1.456e + 02 1.804e-01 6.969e-01 3.599e + 02 1.454e-03
2.231e + 04 6.383e + 00 4.040e + 03 5.358e + 03 2.853e + 01 2.772e + 02 2.047e + 02 7.823e + 01 8.257e + 00
P-value ± ± ± ± ± ± ± ± ±
5.504e + 08 8.408e + 00 1.977e + 07 4.910e + 07 7.016e + 02 1.116e + 05 4.580e + 05 9.604e + 02 7.404e + 02
7.251e-10 1.492e-04 4.759e-09 4.046e-07 4.455e-02 4.758e-08 7.074e-02 2.900e-11 6.951e-02
and Burg’s method, respectively. In this method, the P1, P2, P3 and their corresponding frequencies are relatively distinct from each other in all the three classes. However, Burg’s method is more effective as compared to the Yule–Walker method because the ‘p-values’ shown in Table 3 are very low. This can also be evaluated using ROC curve. 4.2. Receiver operating characteristic test Receiver operating characteristic (ROC) test was performed using the parameters F1, F2, and F3. ROC is a one-class classification method, which classifies a target class within a dataset. The dataset consists of target (desired) values or vectors and outliers. To show that the frequency locations of the prominent frequencies are distinct, the alcoholic and epileptic data were treated as target and the control data was treated as outlier. Figs. 4–6 show the results of the ROC test based on PSDs
Table 3 ANOVA results for Burg’s method P/F P1 F1 P1/F1 P2 F2 P2/F2 P3 F3 P3/F3
Control 2.8141 + 01 1.567e + 00 5.576e + 01 4.477e + 00 1.740e + 01 4.855e-01 5.601e-01 3.497e + 01 2.167e-02
Alcohol ± ± ± ± ± ± ± ± ±
5.732e + 02 1.820e + 01 2.415e + 03 7.063e + 01 7.313e + 01 1.068e + 00 4.046e-01 1.848e + 02 6.441e-04
2.494e + 01 4.200e + 00 4.914e + 01 2.124e + 00 2.232e + 01 2.048e-01 4.046e-01 4.227e + 01 1.623e-02
Epileptic ± ± ± ± ± ± ± ± ±
1.309e + 03 4.641e + 01 5.307e + 03 1.736e + 01 1.951e + 02 1.941e-01 7.086e-01 3.857e + 02 1.476e-03
2.264e + 04 6.350e + 00 4.157e + 03 5.497e + 03 2.860e + 01 2.898e + 02 2.617e + 02 7.647e + 01 1.066e + 01
P-value ± ± ± ± ± ± ± ± ±
5.755e + 08 8.761e + 00 2.027e + 07 4.972e + 07 7.553e + 02 1.154e + 05 5.639 + 05 1.070e + 03 9.224e + 02
9.310e-10 1.466e-03 2.782e-09 2.476e-07 6.812e-02 2.137e-08 3.058e-02 7.889e-10 2.914e-02
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
Fig. 4. ROC results using FFT (Welch’s) method.
obtained with the Welch’s, Yule–Walker and Burg’s method, respectively. The area under the ROC curve indicates whether the target and the outliers are distinct. It can be observed that the ROC test for Welch’s method, shown in Fig. 4, yields low values for the area under the curves. In general, this small area indicates that the target and the outliers are similar. So, it is difficult to distinguish the control data from either alcoholic or epileptic data. It is because, in Welch’s method (Fig. 1), first three spikes are mistaken to be the prominent frequencies. The parametric methods perform much better in the ROC test, as shown in Figs. 5 and 6. The area under the ROC curve for the epileptic group is the same for Burg’s and Yule–Walker method. But for the alcohol dataset, Burg’s method performs slightly better than the Yule–Walker method.
Fig. 5. ROC results using Yule–Walker method.
51
Fig. 6. ROC results using Burg’s method.
5. Discussion In this discussion section, we review some of the work which preceded the current study. The discussion provides further insight in the subject, which is necessary to understand relevance and distinctness of our work. In other studies, dynamical properties like correlation dimension, largest Lyapunov exponent, Hurst exponent, entropy, were evaluated for alcoholic, normal, and epileptic EEG signals [14]. The results show that the alcoholic EEG is more complex than the epileptic seizure signal. And they also showed that the normal EEG was more complex compared to the other two types, indicating higher neuronal process in the brain during the normal active state. In the normal stage, the cortex is more active and more neurons are available for processing [26]. And in epileptic stage, the cortex becomes inactive and EEG becomes less random. The neurons in the cerebral hemispheres during the seizure misfire and create abnormal electrical activity. People with epilepsy have seizures that are a bit like an electrical brainstorm. Hence, the number of neurons available for processing the information reduces during the seizures. Authors have showed that the variability of epileptic activity was less as compared to that of nonepileptic activity [27,28]. It was supported by the reduced dimensionality of epileptic seizures as compared to nonepileptic EEG signals. This concept finds support in the observations that neuronal hypersynchrony underlies seizures: a phenomenon during which the number of independent variables required to describe the system was smaller [29]. The algorithm that we used is sensitive towards noise because only three samples are considered for the peak amplitude detection. Therefore, this algorithm works best for signals that contain mostly low frequencies.
52
O. Faust et al. / ITBM-RBM 29 (2008) 44–52
The PSD estimation using Welch’s method yields signals with high frequency content (presence of spikes). Each of these spikes constitute for local maxima. One way to reduce the number of local maxima is to smooth the signal. However, this smoothing may alter the information content in the signal. In this particular case, Welch’s method is good in discriminating closely spaced frequency peaks. The smoothing function might combine two distinct peaks to one. Parametric PSD estimation methods, like Burg’s and Yule–Walker, yield a smooth spectrum. Therefore, the simple maxima detection algorithm works well with these signals. Thus, the nine features extracted from these PSD plots are distinct for control, epileptic, alcoholic groups. These nine features can constitute a feature vector which is fed into a classifier. Such a classifier maps individual feature vectors to one of three different classes. 6. Conclusion The EEG signal can be used to diagnose different mental conditions like active, epilepsy and alcoholic. The changes in EEG signals might be quite prominent, as in the case of an epileptic patient or more hidden (complex), as in the case of alcoholic subject. In the time domain, only a trained eye can detect these different states. This work shows that, these changes in the states are also visible in the spectral domain. We used three different ways (Welch’s, Yule–Walker and Burg’s) to represent the power distribution in the frequency domain and compared their performances in Tables 1–3. The nine parameters evaluated in the results show that, the parametric methods are superior compared to the nonparametric methods. We have distinguished between alcoholic, epileptic and control EEG signals using ROC method. Based on the ROC test we come to the conclusion that Burg’s method yields the most distinguishable parameters. This study can be seen as the groundwork for more sophisticated EEG signal classification. This signal classification is another step towards an automated system that is able to diagnose different mental conditions based on EEG signals. Such a system would significantly improve clinical workflows, because it frees up trained personal from routine jobs. References [1] James C., Lowe D. Using dynamical embedding to isolate seizure components in the ictal EEG, 1st International Conference on Advances in Medical Signal and Information Processing 2000, 476, 2000, 158–165. [2] Geocadin R, Ghodadra R, Kimura T, Lie H, Sherman D, Hanley D, et al. A novel quantitative EEG injury measure of global cerebral ischemia. Clin Neurophysiol 2000;111:1779–87. [3] Bezerianos A., Tong S., Malhorta A., Thakor N. Information measures of brain dynamics, Nonlinear signal and image processing conference, Baltimore, MD, 2001. [4] Rosso O, Blanco S, Yordanova J, Kolev V, Figliola A, Schurmann M, et al. Wavelet entropy: a new tool for analysis of short duration brain electrical signals. J Neurosci Methods 2001;105:65–75.
[5] Bai O, Nakamura M, Ikeda A, Shibasaki H. Nonlinear Markov process amplitude EEG model for nonlinear coupling interaction of spontaneous EEG. IEEE Trans Biomed Eng 2000;47:1141–6. [6] Guler I, Kiymik MK, Akin M, Alkan A. AR spectral analysis of EEG signals by using maximum likelihood estimation. Med Biol 2001;31:441–50. [7] Akin M, Kiymik MK. Application of periodogram and AR spectral analysis to EEG signals. J Med Syst 2000;24(4):247–56. [8] Fisher RS, Webber WR, Lesser RP, Arroyo S, Uematsu S. High-frequency EEG activity at the start of seizures. J Clin Neurophysiol 1992;9(3):441–8. [9] Allen PJ, Fish DR, Smith SJ. Very high-frequency rhythmic activity during SEEG suppression in frontal lobe epilepsy. Electroencephalogr Clin Neurophysiol 1992;82(2):155–62. [10] Panzica F, Franceschetti S, Binelli S, Canafoglia L, Granata T, Avanzini G. Spectral properties of EEG fast activity ictal discharges associated with infantile spasms. Clin Neurophysiol 1999;110(4):593–603. [11] Mackenzie L, Medvedev A, Hiscock JJ, Pope KJ, Willoughby JO. Picrotoxin-induced generalised convulsive seizure in rat: changes in regional distribution and frequency of the power of electroencephalogram rhythms. Clin Neurophysiol 2002;113(4):586–96. [12] Al-Nashash HA, Zalzala AMS, Thakor NV. A neural networks approach to EEG signals modeling. Med Biol 2003;3:1–6. [13] Al-Nashash Ha, Al-Assaf Y, Paul J, Thakor N. EEG signal modeling using adaptive Markov process amplitude. IEEE Trans Biomed Eng 2004;51(5):744–51. [14] Kannathal N, Acharya RU, Lim Choo M, Quek W, Mohamed H, Sadasivan PK. Characterization of EEG: a comparative study. Comput Methods Programs Biomed 2005;80(1):17–23. [15] Kannathal N, Rajendra Acharya U, Fadilah A, Thelma T, Sadasivan PK. Nonlinear analysis of EEG signals at different mental states. Biomedical 2004 16 March;3(7). [16] American Electroencephalographic Association 1990. Standard electrode position nomenclature,