- Email: [email protected]
A low computation cost method for seizure prediction
Accepted Manuscript Title: A Low Computation Cost Method for Seizure Prediction Author: Yanli Zhang Weidong Zhou Qi Yuan Qi Wu PII: DOI: Reference:
S0920-1211(14)00165-X http://dx.doi.org/doi:10.1016/j.eplepsyres.2014.06.007 EPIRES 5169
To appear in:
Epilepsy Research
Received date: Revised date: Accepted date:
30-8-2013 12-6-2014 17-6-2014
Please cite this article as: Zhang, Y., Zhou, W., Yuan, Q., Wu, Q.,A Low Computation Cost Method for Seizure Prediction, Epilepsy Research (2014), http://dx.doi.org/10.1016/j.eplepsyres.2014.06.007 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
A Low Computation Cost Method for Seizure Prediction Yanli Zhang a,b,c, Weidong Zhou a,c, Qi Yuan a,c, Qi Wu a,c
ip t
School of Information Science and Engineering, Shandong University, Jinan 250100, China b
School of Information and Electronics Engineering, Shandong Institute of Business
Suzhou Institute, Shandong University, Suzhou 215123, China
pt
ed
M
an
c
us
and Technology, Yantai 264005,China
cr
a
Ac ce
The corresponding author: Prof. Weidong Zhou
School of Information Science and Engineering Shandong University
27 Shanda Road, Jinan, P.R.China. Tel.:+86-531-88361551 E-mail Address: [email protected]
Page 1 of 19
HIGHLIGHTS The Higuchi fractal dimension (HFD) changes of EEG recordings can serve as a precursor of epileptic seizures and be used for seizure prediction. The seizure prediction algorithm by combining HFD with Bayesian linear discriminant analysis (BLDA) achieves a high performance. Both HFD and BLDA classifier have a low computational complexity and are suitable for
ip t
real-time seizure prediction.
cr
ABSTRACT
The dynamic changes of electroencephalograph (EEG) signals in the period prior to epileptic
us
seizures play a major role in the seizure prediction. This paper proposes a low computation seizure prediction algorithm that combines a fractal dimension with a machine learning algorithm. The presented seizure prediction algorithm extracts the Higuchi fractal dimension (HFD) of
an
EEG signals as features to classify the patient's preictal or interictal state with Bayesian linear discriminant analysis (BLDA) as a classifier. The outputs of BLDA are smoothed by a Kalman filter for reducing possible sporadic and isolated false alarms and then the final prediction results
M
are produced using a thresholding procedure. The algorithm was evaluated on the intracranial EEG recordings of 21 patients in the Freiburg EEG database. For seizure occurrence period of 30 min and 50 min, our algorithm obtained an average
ed
sensitivity of 86.95% and 89.33%, an average false prediction rate of 0.20/h, and an average prediction time of 24.47 min and 39.39 min, respectively. The results confirm that the changes of HFD can serve as a precursor of ictal activities and be used for distinguishing between interictal
pt
and preictal epochs. Both HFD and BLDA classifier have a low computational complexity. All of
Ac ce
these make the proposed algorithm suitable for real-time seizure prediction. Keywords: EEG, seizure prediction, Higuchi fractal dimension, Bayesian linear discriminant analysis, Kalman filtering
1 Introduction
Epilepsy affects approximately 1% of the word’s population by the sudden and recurrent occurrence of seizures. The epileptic seizures arise from abnormally excessive neuronal spontaneous, synchronized discharges in the cerebral cortex, and there are a significant number of epileptic patients at risk of serious injury or death. Therefore, the ability to effectively predict impending seizures would significantly improve the quality of life of epileptic patients, and provide novel methods for the prevention and treatment of epilepsy. A seizure warning device would augment the therapeutic options or allow behavioral adjustments of those epilepsy patients. Additionally, various intervention systems could be applied to suppress seizures in advance of
Page 2 of 19
their clinical manifestation, e.g., by delivering short-acting antiepileptic drugs or by applying electric stimulation (Schelter et al., 2007). Research on seizure prediction through electroencephalograph (EEG) recordings began in the 1970s. Since then, a large number of studies have been carried out and great efforts have been spent on seizure prediction through EEG monitoring. It has long been observed that the transition from the interictal state (far from seizures) to the ictal state (seizure) is not sudden and may be preceded from minutes to hours by clinical, metabolic or electrical changes, that is there exists a
ip t
preictal state (Lehnertz et al., 2001).
Most current seizure prediction approaches can be summarized into two steps. The first is
cr
extracting linear or nonlinear measurements of scalp or intracranial EEG signals recorded from epileptic patients. The second is classifying them into a preictal or interictal state using statistical
us
analysis or other machine learning techniques (Mirowski et al., 2009). Among these measures, those taken from the theory of chaotic dynamics, including the largest Lyapunov exponent (Iasemidis et al., 1990), correlation dimension (Lehnertz et al., 1998), dynamic similarity index
an
(Le Van Quyen et al., 2001), entropy (Drongelen et al., 2003; Li et al., 2007), and phase synchronization (Mormann et al., 2003; Kuhlmann et al., 2010), have shown higher seizure predictability power than linear measurements such as statistical moments, power spectra
M
features and so on (Mormann et al., 2007). However, to calculate much of these nonlinear measurements, phase space reconstruction is required which increases the algorithm complexity and calculation burden significantly. A nonlinear measurement which can be calculated directly in
ed
the time domain is desired in a real-time seizure prediction system, and Higuchi fractal dimension (HFD) is such an appropriate measurement to quantify signal’s dimensional complexity (Higuchi, 1988), which can be used to study changes in the dynamical behavior of the brain, due to the
pt
following reasons. Though many algorithms are available to compute FD in the time domain, like those proposed by Katz (Katz, 1988) and Petrosian (Petrosian, 1995), HFD provides more
Ac ce
accurate estimation of the fractal dimension than the methods of Katz and, compared with Petrosian’s algorithm, HFD does not depend on a binary sequence and it is less sensitive to noise in many cases (Esteller et al., 2001). HFD has already been used to quantify the complexity of EEG recordings (Khoa et al., 2012; Schneider et al., 2009; Bao et al., 2008). In this paper, HFD is applied to characterize the nonlinear dynamical properties underlying EEG signals for seizure prediction. A good classifier is essential for an excellent seizure prediction method. Bayesian linear discriminant analysis (BLDA) algorithm has shown its good performance for binary classification in the area of brain-computer interface (Hoffmann et al., 2008), due to the fact that it avoids over-fitting phenomenon and time-consuming cross-validation. In this work, this method is employed to classify the HFD features of intracranial EEG recordings in the preictal and interictal state. 2 Material and methods 2.1 EEG dataset and preprocessing
Page 3 of 19
The EEG recordings used in this work came from the Epilepsy Center of the University Hospital of Freiburg, Germany. This database contains intracranial EEG signals of 21 patients recorded during an invasive pre-surgical epilepsy monitoring. These EEG data were obtained using a Neurofile NT digital video EEG system with 128 channels, 256 Hz sampling rate, and a 16-bit analogue-to-digital converter. Six contacts of all implanted grid, strip, and depth electrodes were selected before any analyses by visual inspection of the raw data by a certified epileptologist (Schelter et al., 2006). Three are focal channels (located near the epileptic focus) and three
ip t
extra-focal.
The EEG recordings for each patient contain separate ictal and interictal datasets. The ictal
cr
dataset contains between 2 and 5 epileptic seizures notated by the experienced epileptologists and
at least 50 minutes of preictal data for most seizures. The onset and offset times of seizures had
us
been previously determined by the experienced epileptologists based on identification of epileptic patterns preceding clinical manifestation of seizures in IEEG recordings (Aarabi et al., 2012). Approximate 24h interictal EEG recordings without seizure activity are available in interictal
an
datasets.
Only three focal channels of 21 patients were chosen for seizure prediction in this work. Firstly the moving-window technique was used to split long term EEG recordings into 4s epochs band of EEG epochs to 0.5-30 Hz.
M
without overlap. Then the Chebyshev bandpass digital filters were used to limit the frequency
2.2 Extraction of Higuchi fractal dimension
ed
As the reconstruction of the attractor phase space is not necessary, the calculation of Higuchi fractal dimension is simpler and faster than other classical measures derived from chaos theory. A short introduction to the computation of HFD is included in Appendix A.
pt
For each epoch, the HFD of the EEG segment from each of the three selected channels were calculated and stacked into a three dimensional feature vector.
Ac ce
2.3 Machine learning classification
Once features had been extracted from EEG recordings, classification as preictal or interictal was carried out. The classifier used in this work was Bayesian linear discriminant analysis (Hoffmann et al., 2008). As an extension of Fisher's linear discriminant analysis (FLDA), Bayesian linear discriminant analysis (BLDA) employs regularization to avoid over-fitting to high dimensional and possibly noisy datasets (Hoffmann et al., 2008). A brief description of BLDA is presented in Appendix B. For each patient from the Freiburg database, the available EEG recordings were divided into a training set and a testing set. We set apart the preictal samples of a randomly selected seizure and the same amount of interictal samples chosen randomly, serving as training data. The preictal training data came from the whole 50 min preceding the selected seizure. The remaining samples of the patient were testing data. In summary, we trained the classifier for each patient using the training dataset, and evaluated these classifiers on their testing data.
Page 4 of 19
2.4 Kalman filtering and thresholding To reduce the fluctuations in the outputs of the BLDA and remove the possible sporadic and isolated false alarms, a second-order discrete-time Kalman filter (Chisci et al., 2010), was used for smoothing the BLDA outputs. Kalman filter is an optimal recursive data processing algorithm, which produces the best estimates of the true values of measurements. The output of the BLDA was modeled as the result of a noisy measurement process, i.e., as xk yk vk . Here vk represents
ip t
zero-mean measurement noise with standard deviation σv and Kalman filter was applied to provide an optimal estimation for yk by the measured xk. Introducing the state vector yk then the observation equation and the state equation could be expressed as
(1)
us
yk 1 1 TP yk y 0 1 y wk k k 1
cr
y xk 1 0 k vk y k
y k ,
(2)
white noise with covariance
TP2 2 . 2 wTP
w2
M
2 TP3 w 3 Q 2 TP2 w 2
an
Where Tp was the prediction interval and wk was the process disturbance assumed to be zero-mean
(3)
Kalman filter.
ed
The standard deviations w and v of y k and vk were the only design parameters of the After Kalman filtering, a thresholding method was performed and the outputs of the Kalman
pt
filter were compared with a threshold. If the output for an EEG epoch was greater than the threshold, the epoch was labeled as preictal, otherwise it was marked as an interictal epoch.
Ac ce
In order to set the threshold for each patient, all the preictal training samples of the patient were input into the trained BLDA classifier and the median of the BLDA outputs was represented using the letter m. Then the threshold was given by threshold k m
(4)
The parameter k could be adjusted and thus the threshold was specific for each patient. In the adjustment process of parameter k, there was a trade-off between a high seizure prediction capability and a low false-positive rate. First, the BLDA outputs obtained for the training samples were collected. Then the value of k, minimizing the number of falsely classified training samples, was selected and used to calculate the threshold for each patient. To further reduce the number of false alarm and guarantee a high sensitivity, the alarm for an imminent seizure would be given only when at least two consecutive EEG epochs were marked as preictal. And the seizure was considered to have a successful prediction as long as one alarm had been given for the seizure, while no matter whether the output of the Kalman filter after the alarm
Page 5 of 19
was larger than the threshold or not. 2.5 Performance assessment The performance of the proposed seizure prediction method was assessed on all the testing data of 21 patients in Freiburg database with the following parameters (Winterhalder et al., 2003). Seizure occurrence period (SOP): The period during which a seizure is expected after an alarm. An alarm without a following seizure onset during the SOP will be regarded as a false prediction. Seizure prediction horizon (SPH): The time interval between the alarm and the start of the SOP.
ip t
For a true prediction, after the alarm there is no seizure occurring during the SPH, but a seizure occurs during the SOP. of total seizure events.
us
False prediction rate (FPR): The number of false alarms per hour.
cr
Sensitivity: The measure is defined as the ratio of the number of true predictions to the number
Prediction time: The time interval of the true alarm prior to the onset of the predicted seizure.
an
3 Results 3.1 Overall performance
In this work, the SOP was set to 30 and 50 min, and the SPH was set to 2 min, which is an
M
appropriate value to enable, for instance, drugs administration or warning of the patient (Schelter B et al., 2006). Table 1 lists the results of performance evaluation of this proposed prediction method. The seizure prediction algorithm was tested on the EEG recordings of 21 patients in the
ed
Freiburg EEG database, containing a total of 66 preictal periods and 518.8 h of interictal data. On the whole, an average sensitivity of 86.95% and 89.33% with a same average FPR of 0.20/h were obtained, respectively for SOP = 30 min and 50 min. It can be noticed from Table 1 that there
pt
were 12 patients having 100% sensitivity for SOP = 30 min and 13 patients for SOP = 50 min. In other words, one seizure was missed for 9 or 8 patients with these SOP values. The FPRs of 10
Ac ce
patients were less than 0.15/h for SOP = 30 and 11 patients for SOP =50 min. Among them, four patients had all seizures predicted without false prediction for both SOPs. More than two thirds of all the patients, 15 patients, achieved average prediction time greater than twenty minutes for SOP = 30 min and thirty-five minutes for SOP=50 min. 3.2 Seizure predictability of HFD Investigating the characteristic changes in the preictal state with respect to the interictal state is significant for developing an effective seizure prediction algorithm. In this work, the HFD was used to estimate the changes in the dynamical characteristics underlying the EEG signals, through quantifying the degree of complexity and irregularity of the EEG recordings. As an example, Fig.1 presents two single-channel EEG recordings, which have relatively high and low HFD value respectively. Fig.2 presents the time evolution of the HFD of one-channel EEG recordings from patient 3 and patient 9 and illustrates the proposed seizure prediction method. In Fig.2, the EEG contains two seizures which have been correctly predicted by the proposed algorithm. Compared to the HFD of interictal EEG, the HFD shows a general increase in the preictal period for patient 3
Page 6 of 19
and a general drop for patient 9. To determine the significance of preictal changes of HFD for each patient, the HFD of preictal EEG
was
statistically
compared
with
those
extracted
in
interictal
period
using
one-way analysis of variance. A significance level of P<0.05 was set. Results indicated that 39.1% of the seizures studied in this paper showed a significant drop in HFD prior to the seizure, while 49.4% showed an increase in the preictal state. In the remaining 11.5% of the seizures, no significant changes were observed in the preictal state compared to the interictal state.
ip t
The proposed algorithm might miss the seizures in whose preictal state the changes of HFD
were unapparent or in an opposite direction relative to other seizures of the same patient.
cr
To illustrate this point, Table 2 provides the data of seven patients, including the mean values of interictal and preictal HFD, the minimum and maximum of the outputs of BLDA and Kalman
us
filter, the classification threshold, and the number of epochs classified as preictal by the algorithm. For a seizure event, if none of its preictal epochs were detected correctly, it was not predicted for both SOPs by our proposed algorithm. It can be noted from Table 2 that for six patients (excepting
an
patient 9), the missed seizures by the algorithm have the preictal HFD approximately equal to the interictal HFD or in an opposite change direction relative to other seizures. The missed seizure of patient 9 was maybe because of the heterogeneity of seizures.
M
3.3 Property of Kalman filter
Kalman filter was used to smooth the outputs of the BLDA classifier and to reduce the possible false seizure predictions. In order to assess the effectiveness of the Kalman filter, we
ed
additionally tested other smoothing strategies on the Freiburg long-term EEG recordings, such as moving average filter (Temko et al., 2011). The moving average filter is a low-pass filter used to reduce random noise. The value of a sample is calculated by the average of its neighbors: N 1 x m n 2 N 1 n N
pt y m
(5)
Ac ce
Where x is the input signal, y denotes the output (filtered) signal, and 2N+1 is the order of the filter, that is the number of points used in the moving average. Table 3 reports a comparison between the Kalman filter and the moving average filter. In
comparison to the moving average filter, Kalman filter produced a higher average sensitivity of 89.33% and also an improvement on the average prediction time, meanwhile the average false prediction rate reduced to 0.20 per hour. Besides, Kalman filter did not introduce an undesirable time delay for forecasting purposes like the moving average filter. 4 Discussion HFD was applied to characterize the nonlinear dynamical properties underlying EEG signals for seizure prediction in this work. HFD measures the complexity and irregularity of the EEG signals. The observed changes in HFD are essential for a prediction technique to distinguish between interictal and preictal epochs, regardless of the dimensions drop or increase in the preictal
Page 7 of 19
state. The measure of HFD can serve as a precursor of ictal activity for epileptic patients, and the ability to predict epileptic seizures prior to their occurrences was well illustrated by the results in Table1. 39.1% of the seizures studied in this paper showed a significant drop in HFD prior to the seizure, 49.4% showed an increase and the remaining 11.5% of the seizures had no significant changes observed in the preictal state compared to the interictal state. These results were concordant with the findings of Aarabi et al. (2012) who studied the correlation dimension and reported that 44.9% of the seizures they analyzed showed a significant drop in the values of the
ip t
correlation dimension in the preictal state, 44.9% showed an increase and 10.2% of the seizures were no significant changes observed between the preictal and interictal states. Compared with the
cr
correlation dimension widely used in seizure prediction study, HFD has no need of reconstructing
phase space for EEG signals and can be calculated directly in the time domain. Hence it has a low
us
computational burden which is suitable for real-time prediction.
To lower the amount of calculation, we only chosen three focal channels of 21 patients for seizure prediction in this work and the EEG recordings of three extra-focal channels were not
an
analyzed. Many studies have reported significant preictal changes of univariate nonlinear measures in channels corresponding to the epileptogenic zone (Lehnertz and Elger, 1998; Maiwald of seizure prediction in this work.
M
et al., 2004), and restricting the analysis to three focal channels was without loss in performance Table 4 presents a comparison on the prediction results between the proposed method in this paper and some other methods for SOP = 30 min and 50 min. All the prediction methods in Table
ed
4 were evaluated on the same Freiburg EEG dataset. Winterhalder et al. introduced a performance assessment criterion of seizure prediction methods and investigated changes in the dynamical similarity index (2003) and the synchronization (2006) of the preictal states. Maiwald et al. (2004)
pt
compared the performance of three different patient specific prediction methods including the dynamic similarity index, the accumulated energy, and the effective correlation dimension.
Ac ce
Feldwisch-Drentrup et al. (2010) developed two combinational seizure prediction system by using logical ‘‘AND’’ and ‘‘OR’’ combinations on the mean phase coherence and the dynamic similarity index of EEG recordings. Though these methods had a low FPR, the sensitivity made them insufficient for clinical applications. Compared to these previous works, the presented method improved the sensitivity.
In the recent work of Aarabi et al. (2012), five univariate measures including correlation
dimension, correlation entropy, noise level, Lempel-Ziv complexity, and largest Lyapunov exponent as well as one bivariate measure, nonlinear interdependence, were extracted from intracranial EEG recordings. The spatio-temporal information was then integrated by using patient-specific rules established from the characteristics changes of the preictal states to obtain maximum sensitivity and specificity for seizure prediction. They reported an average sensitivity of 79.9% with an average FPR of 0.17/h for SOP of 30 min and SPH of 10s. For SOP of 50 min, a higher average sensitivity (90.2%) and lower average FPR (0.11/h) were obtained. Compared with
Page 8 of 19
Aarabi’s system, our proposed algorithm provided comparable sensitivity, albeit with a little higher false positive rate. It is worth mentioning that the performance results of Aarabi’s prediction system were calculated only on 11 patients of the EEG dataset, whereas our proposed method was tested on the EEG data of all the 21 patients. In addition, the seizure prediction method of Aarabi needed to extract six features of EEG recordings, which would yield very high computational burden and make the method unrealizable in a real-time prediction system. In contrast, the seizure prediction algorithm presented here, combining the HFD with the BLDA, has
ip t
lower computational complexity because the calculation of HFD does not require phase space reconstruction and the BLDA’s parameter estimation avoids time consuming cross-validation.
cr
To evaluate the performance of our seizure prediction algorithm, we also compared the sensitivities with the ones obtained using the random and the periodical prediction methods as
us
described by Winterhalder (2003). In terms of the sensitivity, our method performed better than the random prediction method (9.52% and 15.35% for SOP = 30 and 50 min, respectively) and the periodical prediction method (10% and 16.67% for SOP = 30 and 50 min, respectively).
an
Although the results are satisfactory compared to that of the other methods in the literature, the proposed seizure prediction algorithm still has a long way from clinical applicability yet, especially, the average FPR of 0.2/h is a little high. The patients in the EEG database used in the
M
work were during monitoring for presurgical evaluation and had a mean seizure frequency of about 3 seizures per day which was higher than that of normal conditions. The proposed seizure
5 Conclusion
ed
prediction algorithm would be improved in the FPR and tested in clinical practice.
In this work, we have investigated the changes of HFD in the preictal state and developed a
pt
patient-specific seizure prediction algorithm with low computational complexity. The results demonstrate that the changes of HFD can serve as a precursor of ictal activities and the proposed
Ac ce
seizure prediction algorithm, by combining HFD features with the patient specific BLDA classifier, achieved high sensitivity and lower false prediction rate when tested on the Freiburg data set. The selection of HFD feature and BLDA classifier makes the prediction have a low computational complexity, and the utility of Kalman filter introduces no further time delay when it is used to smooth the outputs of BLDA. All of these characteristics make the proposed algorithm potentially suitable for real-time seizure prediction. Acknowledgment This work was supported by the Key Program of Natural Science Foundation of Shandong Province (No. ZR2013FZ002), the Program of Science and Technology of Suzhou under Grant ZXY2013030 and the Independent Innovation Foundation of Shandong University under Grant 2012DX008.
Appendix A. Higuchi fractal dimension (HFD)
Page 9 of 19
Given a one dimensional time series x = {x[1], x[2], ..., x[N]}, the algorithm to compute the HFD can be described as follows (Higuchi, 1988 ). The first step is to form k new time series xmk , which is defined by N m xmk x m , x m k , x m 2k ,, x m int k , k
(A.1)
for m=1,2, ..., k. Where m and k are integers, representing the initial time value and the discrete int(r) is the integer part of a real number r. The length of each new time series is then calculated as
cr
int N k m 1 N 1 Lm k x m ik x m i 1 k , N m k i 1 int k k
ip t
time interval between points, respectively. N is the total length of the original time series x, and
us
N 1 is a normalization factor. N m k int k
an
Where
(A.2)
Thus the average length of all time series having the same interval k, is computed as 1 k Lm k . k m 1
M
Lk
The average length L (k) is proportional to k
–D
(A.3)
, where D is the fractal dimension by
ed
Higuchi's method. Therefore, the procedure described above repeats kmax times for k =1,2, ..., kmax, and then the least-square fitting method is used to determine a straight line that best fits the given
pt
data set of {ln(L(k)), ln(1/k)} (k =1,2, ..., kmax,). The straight line determined by the least squares
Ac ce
method is the best-fit curve that has the minimal sum of the deviations squared (least square error) from the given set of data. The slope of the straight line is the Higuchi fractal dimension of time series x. In order to choose an appropriate value of the parameter kmax, HFD values were plotted against a range of kmax. The point at which the FD plateaus is considered a saturation point and that kmax value should be selected (Gómez et,al.,2009). Appendix B. Bayesian linear discriminant analysis (BLDA) The following is a brief description of the key ideas of Bayesian linear discriminant analysis (Hoffmann et al., 2008). As an extension of Fisher's linear discriminant analysis (FLDA), BLDA can be regarded as least squares regression in a Bayesian framework. Assuming that feature vector x and regression target, that is class label t satisfy t T x n ,
(B.1)
Page 10 of 19
Where n denotes the additive white Gaussian noise, and is the regression weight. Then the likelihood function of obey a Gaussian distribution, N
2 2 T p D , exp X t 2 2
(B.2)
Where X is a matrix with rows corresponding to training feature vectors, N dimensional column vector t contains the regression targets of all the training examples, D denotes the pair {X, t},
1
G
2 2 1 T p exp I 2 2 2
(B.3)
cr
Then the prior distribution for the regression weight is specified as
ip t
and represents the inverse variance of the noise.
Where is a very small value selected randomly.
(B.4)
an
0 0 0 0 I 0 0
us
Here G denotes the number of features, and I is a G+1 dimensional, diagonal, square matrix
p D , p
p D , p d
ed
p , , D
M
The posterior distribution of the weight can be obtained based on Bayes rule, (B.5)
Thus the linear discriminant function is given by y m T xˆ
(B.6)
pt
Where xˆ is a new input vector waiting for being classified, m is the mean of the posterior Gaussian distribution, which is computed as
Ac ce
m XX T I Xt 1
(B.7)
The parameters and which yield the optimal separating hyperplane can be selected with a simple iterative algorithm described by MacKay (1992).
References
Aarabi, A., He, B., 2012. A rule-based seizure prediction method for focal neocortical epilepsy. Clin. Neurophysiol. 123, 1111-1122. Bao, F.S., Lie, D.Y.C., Zhang, Y., 2008. A new approach to automated epileptic diagnosis using EEG and probabilistic neural network. International Conference on Tools with Artificial Intelligence. 2, 482-486. Chisci, L., Mavino, A., Perferi, G., Sciandrone, M., Anile, C., Colicchio, G., Fuggetta, F., 2010. Real-time epileptic seizure prediction using AR models and support vector machines. IEEE
Page 11 of 19
Trans. Biomed. Eng. 57(5), 1124-1132. Drongelen, W.v., Nayak, S., Frim, D.M., Kohrman, M.H., Towle, V.L., Lee, H.C., et al., 2003. Seizure anticipation in pediatric epilepsy: use of Kolmogorov entropy. Pediatr Neurol. 29(3), 207-213. Esteller, R., Vachtsevanos, G., Echauz, J., Litt, B., 2001. A comparison of waveform fractal dimension algorithms, IEEE Trans. Circuits-I. 48(2), 177-183. Feldwisch-Drentrup, H., Schelter, B., Jachan, M., Nawrath, J., Timmer, J., Schulze-Bonhage, A.,
ip t
2010. Joining the benefits: combining epileptic seizure prediction methods. Epilepsia 51(8), 1598-1606.
cr
Gómez, C., Mediavilla, A., Hornero, R., Abásolo, D., Fernández, A., 2009. Use of the Higuchi’s
fractal dimension for the analysis of MEG recordings from Alzheimer’s disease patients. Med.
us
Eng. Phys. 31 (3), 306–313.
Higuchi, T., 1988. Approach to an irregular time series on the basis of the fractal theory. Physica D. 31(2), 277-283.
an
Hoffmann, U., Vesin, J.M., Ebrahimi, T., Diserens, K., 2008. An efficient P300-based brain-computer interface for disabled subjects. J. Neurosci. Method 167(1), 115-125. Iasemidis, L.D., Sackellares, J.C., Zaveri, H.P., Williams, W.J., 1990. Phase space topography and
M
the Lyapunov exponent of electrocorticograms in partial seizure. Brain Topogr. 2(3), 187-201. Katz, M., 1988. Fractals and the analysis of waveforms. Comput. Biol. Med. 18, 145-156. Khoa, T.Q.D., Ha, V.Q., Toi, V.V., 2012. Higuchi Fractal Properties of Onset Epilepsy
ed
Electroencephalogram. Comput. Math. Method M. 2012, 461426. Kuhlmann, L., Freestone, D., Lai, A., et al., 2010. Patient-specific bivariate-synchrony-based seizure prediction for short prediction horizons. Epilepsy Res. 91, 214-231.
pt
Le Van Quyen M., Martinerie, J., Navarro, V., Boon, P., D'Havé, M., Adam, C., Renault, B., Varela, F., Baulac, M., 2001. Anticipation of epileptic seizures from standard EEG recordings. Lancet
Ac ce
357, 183-188.
Lehnertz, K., Elger, C.E., 1998. Can epileptic seizures be predicted? Evidence from nonlinear time series analysis of brain electrical activity. Phys. Rev. Lett. 80(22), 5019-5022. Lehnertz, K., Andrzejak, R.G., Arnhold, J., Kreuz, T., Mormann, F., et al., 2001. Nonlinear EEG Analysis in Epilepsy: Its Possible Use for Interictal Focus Localization, Seizure Anticipation, and Prevention. J Clin Neurophysiol. 18 (3), 209-222. Li, X.L., Ouyang, G.X., Richards, Douglas.A., 2007. Predictability analysis of absence seizures with permutation entropy. Epilepsy Res. 77, 70-74. MacKay, D.J.C., 1992. Bayesian interpolation. Neural Comput. 4, 415-447. Maiwald, T., Winterhalder, M., Aschenbrenner-Scheibe, R., Voss, H.U., Schulze-Bonhage, A., Timmer, J., 2004. Comparison of three nonlinear seizure prediction methods by means of the seizure prediction characteristic. Physica D. 194, 357-368. Mirowski, P., Madhavan, D., LeCun, Y., Kuzniecky, R., 2009. Classification of patterns of EEG
Page 12 of 19
synchronization for seizure prediction. Clin. Neurophysiol. 120, 1927-1940. Mormann, F., Andrzejak, R.G., Kreuz, T., Rieke, C., David, P., Elger, C.E., Lehnertz, K., 2003. Automated detection of a preseizure state based on a decrease in synchronization in intracranial electroencephalogram recordings from epilepsy patients. Phys. Rev. E 67, 021912. Mormann, F., Andrzejak, R.G., Elger, C.E., Lehnertz, K., 2007. Seizure prediction: the long and winding road. Brain 130, 314–333. Osterhage, H., Mormann, F., Staniek, M., Lehnertz, K., 2007. Measuring synchronization in the
ip t
epileptic brain: a comparison of different approaches. Int. J. Bifurcat. Chaos. 17(10), 3539-3544.
cr
Petrosian, A., 1995. Kolmogorov complexity of finite sequences and recognition of different preictal EEG patterns. Proc. IEEE Symp. Comp. Med. Syst. 212-217.
us
Schneider, M., Mustaro, P.N., Lima, C.A.M., 2009. Automatic recognition of epileptic seizure in EEG via support vector machine and dimension fractal. Proceedings of International Joint Conference on Neural Networks. 2841-2845.
an
Schelter, B., Winterhalder, M., Maiwald, T., Brandt, A., Schad, A., Timmer, J., Schulze-Bonhage, A., 2006. Do false predictions of seizures depend on the state of vigilance? A report from two seizure-prediction methods and proposed remedies. Epilepsia 47(12), 2058-2070. horizons. Epilepsy Res. 73, 213-217.
M
Schelter, B., Winterhalder, M.,et al., 2007. Seizure prediction: The impact of long prediction Temko, A., Thomas, E., Marnane, W., Lightbody, G., Boylan, G., 2011. EEG-based neonatal Winterhalder,
M.,
ed
seizure detection with Support Vector Machines. Clin. Neurophysiol. 122, 464-473. Maiwald,
T.,Voss,
H.U.,
Aschenbrenner-Scheibe,
R.,
Timmer,
J.,
Schulze-Bonhage, A., 2003. The seizure prediction characteristic: a general framework to
pt
assess and compare seizure prediction methods. Epilepsy Behav. 4, 318-325. Winterhalder, M., Schelter, B., Maiwald, T., Brandt, A., Schad, A., Schulze-Bonhage, A., Timmer,
Ac ce
J., 2006. Spatio-temporal patient–individual assessment of synchronization changes for epileptic seizure prediction. Clin. Neurophysiol. 117, 2399-2413.
Page 13 of 19
ip t cr us
Fig. 1 Two single-channel EEG recordings. The signals in (a) and (b) have relatively high and low
Ac ce
pt
ed
M
an
HFD values respectively.
Page 14 of 19
ip t us
cr Ac ce
pt
ed
M
an
(a)
(b)
Fig. 2 The HFD of one-channel EEG recordings from patient 3 (a) and patient 9 (b), and the illustration of the seizure prediction method. In the top left panels of (a) and (b), the EEG recordings respectively contain a seizure event with a preictal period of approximately 50 minutes. The two vertical lines indicate the seizure onset and offset respectively, and the seizure onset is labeled as time zero. The top right panels of (a) and (b) respectively illustrate 1h interictal EEG. The middle panels of (a) and (b) present the HFD time-traces of the EEG recordings given in the top panels. The bottom panels illustrate the proposed seizure prediction method. The blue curves and the green profiles respectively display the outputs of BLDA and Kalman filter. The classification thresholds are represented by red horizontal lines. For each patient, the threshold is a constant, that is, the red lines in both the left and right graphs of the bottom panels are at the same level. When the outputs of Kalman filter have at least two consecutive EEG epochs greater than the threshold, an alarm in blue column bar form for the impending seizure will be given.
Page 15 of 19
Table 1 Results for each patient using the proposed prediction method SOP=30min, SPH=2 min
Number of
Number
seizure
of true
events
alarms
1
3
2
SOP=50min, SPH=2 min Average
Number
prediction
of true
time(min)
alarms
0.25
19.78
100.00
0.33
3
75.00
4
4
5
4
6
Average
Sensitivity
FPR
(%)
(/h)
3
100.00
0.33
36.33
29.13
2
100.00
0.21
47.40
0.04
29.84
3
75.00
0.13
36.49
100.00
0.00
28.57
4
100.00
4
100.00
0.46
24.72
4
100.00
2
1
50.00
0.21
29.93
2
100.00
7
2
2
100.00
0.00
16.00
2
8
1
1
100.00
0.04
29.73
1
9
4
3
75.00
0.27
16.69
3
cr
Patient
10
4
4
100.00
0.24
26.07
11
3
2
67.00
0.12
16.10
12
3
2
67.00
0.04
29.30
13
1
1
100.00
0.00
14
3
3
100.00
15
3
2
16
4
17 18
FPR
(%)
(/h)
3
100.00
2
2
3
4
4
prediction time(min)
ip t
Sensitivity
0.00
47.12
0.33
49.72
0.42
33.63
0.00
26.27
100.00
0.36
49.60
75.00
0.27
34.54
4
100.00
0.28
41.47
2
67.00
0.12
49.93
3
100.00
0.12
19.71
12.73
1
100.00
0.00
49.87
0.40
28.49
3
100.00
0.68
39.18
67.00
0.13
24.00
2
67.00
0.08
45.00
4
100.00
0.33
27.97
4
100.00
0.38
49.64
4
4
100.00
0.00
28.48
4
100.00
0.00
43.20
4
3
75.00
0.12
27.98
3
75.00
0.00
43.38
19
3
3
100.00
0.38
27.40
2
67.00
0.08
17.07
20
4
3
75.00
0.58
24.51
3
75.00
0.42
34.16
21
4
3
75.00
0.36
16.44
3
75.00
0.08
33.51
86.95
0.20
24.47
89.33
0.20
39.39
an
M
ed
pt
Ac ce
Average
us
100.00
Page 16 of 19
Table 2 The data of seven patients having one seizure missed prediction for both SOPs
15
18
20
1.0663
1.0425
1.0750
1.0292
1.0467
minimum
maximum
1
1.0774
-1.9164
2.1352
-0.2972
0.7435
2
2
1.0764
-2.3967
1.7881
-1.0182
0.4676
0
3
1.1314
-3.6602
4.2210
-2.6464
1.9088
165
4
1.1822
-4.9540
4.3037
-0.9915
2.0962
276
1
1.0791
-3.7769
3.9563
-1.0265
1.7531
0
2
1.0843
-3.9209
3.7343
-2.4761
1.8906
2
3
1.0661
-1.6516
3.7225
1.1045
2.6769
449
4
1.0766
-4.3212
4.0349
-0.5685
1.8849
9
1
1.0777
-1.6283
1.3617
-1.0683
0.6245
0
2
1.0564
-1.3398
2.1196
-0.3303
1.1204
2
3
1.0639
-1.8746
1.9509
-0.6656
1.1517
8
1
1.0430
-1.3812
2.5734
-0.4781
1.5258
0
2
1.0563
-0.9429
4.2652
-0.3772
2.5456
28
3
1.0537
-0.6365
4.8526
0.7022
3.2622
78
1
1.1962
-1.2933
8.5827
-1.2989
5.2133
263
2
1.1796
-0.8237
4.0800
-0.7299
2.3789
48
3
1.0700
-1.0790
-0.2776
-0.8911
-0.3715
0
4
1.2097
-0.8523
3.5874
-0.6217
1.8732
72
1
1.0344
-1.4375
4.5179
-0.3503
1.7115
22
2
1.0270
-1.9806
2.9048
-1.0831
0.5910
0
1.8153
0.9182
2.4942
1.6382
1.3734
1.7369
Ac ce
ip t
maximum
Outputs of BLDA
cr
1.1136
minimum
3
1.0323
-1.9291
3.5050
-1.1559
1.4498
51
4
1.0378
-1.4811
8.0002
0.0076
2.0560
146
1
1.0470
-1.2382
4.5815
-0.5447
1.8642
2
2
1.0500
-1.1504
2.7341
-0.4145
0.3564
0
3
1.0602
-0.5575
4.6726
0.3230
1.7787
4
4
1.0577
-1.5914
2.7027
-0.3947
2.0336
18
pt
21
0.7375
number of epochs classified as preictal
mean of preictal HFD
us
11
1.0639
Outputs of Kalman filter
Seizure number
an
9
threshold
M
3
mean of interictal HFD
ed
Patient
Page 17 of 19
Table 3 Comparison between Kalman filter and moving average filter SOP=50min, SPH=2min Postprocessing method
Number of seizure events
Number of true alarms
Average Sensitivity (%)
Average
Average prediction
FPR (/h)
time(min)
66
57
87.38
0.26
30.33
Kalman filter
66
58
89.33
0.20
39.39
Ac ce
pt
ed
M
an
us
cr
ip t
Moving average filter
Page 18 of 19
Table 4 Comparison between the proposed method and other methods SOP=30min Year
feature
Sensitivity
FPR
Sensitivity
FPR
(%)
(/h)
(%)
(/h)
2003
dynamical similarity index
21-42
0.15
32-67
0.15
Maiwald
2004
dynamical similarity index
21-42
0.15
52
0.15
correlation dimension
13-30
0.15
46
0.15
accumulated energy
18-31
0.15
34
0.15
mean phase coherence, lag
60
0.15
70
0.15
35.2
0.15
45
43.2
0.15
79.9
Winterhalder
2006
synchronization index 2010
mean phase coherence or
Drentrup
dynamical similarity index mean phase coherence and
cr
Feldwisch-
dynamical similarity index Aarabi
2012
five
univariate
measures
and one bivariate measure Present
Higuchi fractal dimension
86.95
0.17
0.20
53
0.15 0.15
90.2
0.11
89.33
0.20
Ac ce
pt
ed
M
an
This work
ip t
Winterhalder
us
Authors
SOP=50min
Page 19 of 19
S0920-1211(14)00165-X http://dx.doi.org/doi:10.1016/j.eplepsyres.2014.06.007 EPIRES 5169
To appear in:
Epilepsy Research
Received date: Revised date: Accepted date:
30-8-2013 12-6-2014 17-6-2014
Please cite this article as: Zhang, Y., Zhou, W., Yuan, Q., Wu, Q.,A Low Computation Cost Method for Seizure Prediction, Epilepsy Research (2014), http://dx.doi.org/10.1016/j.eplepsyres.2014.06.007 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
A Low Computation Cost Method for Seizure Prediction Yanli Zhang a,b,c, Weidong Zhou a,c, Qi Yuan a,c, Qi Wu a,c
ip t
School of Information Science and Engineering, Shandong University, Jinan 250100, China b
School of Information and Electronics Engineering, Shandong Institute of Business
Suzhou Institute, Shandong University, Suzhou 215123, China
pt
ed
M
an
c
us
and Technology, Yantai 264005,China
cr
a
Ac ce
The corresponding author: Prof. Weidong Zhou
School of Information Science and Engineering Shandong University
27 Shanda Road, Jinan, P.R.China. Tel.:+86-531-88361551 E-mail Address: [email protected]
Page 1 of 19
HIGHLIGHTS The Higuchi fractal dimension (HFD) changes of EEG recordings can serve as a precursor of epileptic seizures and be used for seizure prediction. The seizure prediction algorithm by combining HFD with Bayesian linear discriminant analysis (BLDA) achieves a high performance. Both HFD and BLDA classifier have a low computational complexity and are suitable for
ip t
real-time seizure prediction.
cr
ABSTRACT
The dynamic changes of electroencephalograph (EEG) signals in the period prior to epileptic
us
seizures play a major role in the seizure prediction. This paper proposes a low computation seizure prediction algorithm that combines a fractal dimension with a machine learning algorithm. The presented seizure prediction algorithm extracts the Higuchi fractal dimension (HFD) of
an
EEG signals as features to classify the patient's preictal or interictal state with Bayesian linear discriminant analysis (BLDA) as a classifier. The outputs of BLDA are smoothed by a Kalman filter for reducing possible sporadic and isolated false alarms and then the final prediction results
M
are produced using a thresholding procedure. The algorithm was evaluated on the intracranial EEG recordings of 21 patients in the Freiburg EEG database. For seizure occurrence period of 30 min and 50 min, our algorithm obtained an average
ed
sensitivity of 86.95% and 89.33%, an average false prediction rate of 0.20/h, and an average prediction time of 24.47 min and 39.39 min, respectively. The results confirm that the changes of HFD can serve as a precursor of ictal activities and be used for distinguishing between interictal
pt
and preictal epochs. Both HFD and BLDA classifier have a low computational complexity. All of
Ac ce
these make the proposed algorithm suitable for real-time seizure prediction. Keywords: EEG, seizure prediction, Higuchi fractal dimension, Bayesian linear discriminant analysis, Kalman filtering
1 Introduction
Epilepsy affects approximately 1% of the word’s population by the sudden and recurrent occurrence of seizures. The epileptic seizures arise from abnormally excessive neuronal spontaneous, synchronized discharges in the cerebral cortex, and there are a significant number of epileptic patients at risk of serious injury or death. Therefore, the ability to effectively predict impending seizures would significantly improve the quality of life of epileptic patients, and provide novel methods for the prevention and treatment of epilepsy. A seizure warning device would augment the therapeutic options or allow behavioral adjustments of those epilepsy patients. Additionally, various intervention systems could be applied to suppress seizures in advance of
Page 2 of 19
their clinical manifestation, e.g., by delivering short-acting antiepileptic drugs or by applying electric stimulation (Schelter et al., 2007). Research on seizure prediction through electroencephalograph (EEG) recordings began in the 1970s. Since then, a large number of studies have been carried out and great efforts have been spent on seizure prediction through EEG monitoring. It has long been observed that the transition from the interictal state (far from seizures) to the ictal state (seizure) is not sudden and may be preceded from minutes to hours by clinical, metabolic or electrical changes, that is there exists a
ip t
preictal state (Lehnertz et al., 2001).
Most current seizure prediction approaches can be summarized into two steps. The first is
cr
extracting linear or nonlinear measurements of scalp or intracranial EEG signals recorded from epileptic patients. The second is classifying them into a preictal or interictal state using statistical
us
analysis or other machine learning techniques (Mirowski et al., 2009). Among these measures, those taken from the theory of chaotic dynamics, including the largest Lyapunov exponent (Iasemidis et al., 1990), correlation dimension (Lehnertz et al., 1998), dynamic similarity index
an
(Le Van Quyen et al., 2001), entropy (Drongelen et al., 2003; Li et al., 2007), and phase synchronization (Mormann et al., 2003; Kuhlmann et al., 2010), have shown higher seizure predictability power than linear measurements such as statistical moments, power spectra
M
features and so on (Mormann et al., 2007). However, to calculate much of these nonlinear measurements, phase space reconstruction is required which increases the algorithm complexity and calculation burden significantly. A nonlinear measurement which can be calculated directly in
ed
the time domain is desired in a real-time seizure prediction system, and Higuchi fractal dimension (HFD) is such an appropriate measurement to quantify signal’s dimensional complexity (Higuchi, 1988), which can be used to study changes in the dynamical behavior of the brain, due to the
pt
following reasons. Though many algorithms are available to compute FD in the time domain, like those proposed by Katz (Katz, 1988) and Petrosian (Petrosian, 1995), HFD provides more
Ac ce
accurate estimation of the fractal dimension than the methods of Katz and, compared with Petrosian’s algorithm, HFD does not depend on a binary sequence and it is less sensitive to noise in many cases (Esteller et al., 2001). HFD has already been used to quantify the complexity of EEG recordings (Khoa et al., 2012; Schneider et al., 2009; Bao et al., 2008). In this paper, HFD is applied to characterize the nonlinear dynamical properties underlying EEG signals for seizure prediction. A good classifier is essential for an excellent seizure prediction method. Bayesian linear discriminant analysis (BLDA) algorithm has shown its good performance for binary classification in the area of brain-computer interface (Hoffmann et al., 2008), due to the fact that it avoids over-fitting phenomenon and time-consuming cross-validation. In this work, this method is employed to classify the HFD features of intracranial EEG recordings in the preictal and interictal state. 2 Material and methods 2.1 EEG dataset and preprocessing
Page 3 of 19
The EEG recordings used in this work came from the Epilepsy Center of the University Hospital of Freiburg, Germany. This database contains intracranial EEG signals of 21 patients recorded during an invasive pre-surgical epilepsy monitoring. These EEG data were obtained using a Neurofile NT digital video EEG system with 128 channels, 256 Hz sampling rate, and a 16-bit analogue-to-digital converter. Six contacts of all implanted grid, strip, and depth electrodes were selected before any analyses by visual inspection of the raw data by a certified epileptologist (Schelter et al., 2006). Three are focal channels (located near the epileptic focus) and three
ip t
extra-focal.
The EEG recordings for each patient contain separate ictal and interictal datasets. The ictal
cr
dataset contains between 2 and 5 epileptic seizures notated by the experienced epileptologists and
at least 50 minutes of preictal data for most seizures. The onset and offset times of seizures had
us
been previously determined by the experienced epileptologists based on identification of epileptic patterns preceding clinical manifestation of seizures in IEEG recordings (Aarabi et al., 2012). Approximate 24h interictal EEG recordings without seizure activity are available in interictal
an
datasets.
Only three focal channels of 21 patients were chosen for seizure prediction in this work. Firstly the moving-window technique was used to split long term EEG recordings into 4s epochs band of EEG epochs to 0.5-30 Hz.
M
without overlap. Then the Chebyshev bandpass digital filters were used to limit the frequency
2.2 Extraction of Higuchi fractal dimension
ed
As the reconstruction of the attractor phase space is not necessary, the calculation of Higuchi fractal dimension is simpler and faster than other classical measures derived from chaos theory. A short introduction to the computation of HFD is included in Appendix A.
pt
For each epoch, the HFD of the EEG segment from each of the three selected channels were calculated and stacked into a three dimensional feature vector.
Ac ce
2.3 Machine learning classification
Once features had been extracted from EEG recordings, classification as preictal or interictal was carried out. The classifier used in this work was Bayesian linear discriminant analysis (Hoffmann et al., 2008). As an extension of Fisher's linear discriminant analysis (FLDA), Bayesian linear discriminant analysis (BLDA) employs regularization to avoid over-fitting to high dimensional and possibly noisy datasets (Hoffmann et al., 2008). A brief description of BLDA is presented in Appendix B. For each patient from the Freiburg database, the available EEG recordings were divided into a training set and a testing set. We set apart the preictal samples of a randomly selected seizure and the same amount of interictal samples chosen randomly, serving as training data. The preictal training data came from the whole 50 min preceding the selected seizure. The remaining samples of the patient were testing data. In summary, we trained the classifier for each patient using the training dataset, and evaluated these classifiers on their testing data.
Page 4 of 19
2.4 Kalman filtering and thresholding To reduce the fluctuations in the outputs of the BLDA and remove the possible sporadic and isolated false alarms, a second-order discrete-time Kalman filter (Chisci et al., 2010), was used for smoothing the BLDA outputs. Kalman filter is an optimal recursive data processing algorithm, which produces the best estimates of the true values of measurements. The output of the BLDA was modeled as the result of a noisy measurement process, i.e., as xk yk vk . Here vk represents
ip t
zero-mean measurement noise with standard deviation σv and Kalman filter was applied to provide an optimal estimation for yk by the measured xk. Introducing the state vector yk then the observation equation and the state equation could be expressed as
(1)
us
yk 1 1 TP yk y 0 1 y wk k k 1
cr
y xk 1 0 k vk y k
y k ,
(2)
white noise with covariance
TP2 2 . 2 wTP
w2
M
2 TP3 w 3 Q 2 TP2 w 2
an
Where Tp was the prediction interval and wk was the process disturbance assumed to be zero-mean
(3)
Kalman filter.
ed
The standard deviations w and v of y k and vk were the only design parameters of the After Kalman filtering, a thresholding method was performed and the outputs of the Kalman
pt
filter were compared with a threshold. If the output for an EEG epoch was greater than the threshold, the epoch was labeled as preictal, otherwise it was marked as an interictal epoch.
Ac ce
In order to set the threshold for each patient, all the preictal training samples of the patient were input into the trained BLDA classifier and the median of the BLDA outputs was represented using the letter m. Then the threshold was given by threshold k m
(4)
The parameter k could be adjusted and thus the threshold was specific for each patient. In the adjustment process of parameter k, there was a trade-off between a high seizure prediction capability and a low false-positive rate. First, the BLDA outputs obtained for the training samples were collected. Then the value of k, minimizing the number of falsely classified training samples, was selected and used to calculate the threshold for each patient. To further reduce the number of false alarm and guarantee a high sensitivity, the alarm for an imminent seizure would be given only when at least two consecutive EEG epochs were marked as preictal. And the seizure was considered to have a successful prediction as long as one alarm had been given for the seizure, while no matter whether the output of the Kalman filter after the alarm
Page 5 of 19
was larger than the threshold or not. 2.5 Performance assessment The performance of the proposed seizure prediction method was assessed on all the testing data of 21 patients in Freiburg database with the following parameters (Winterhalder et al., 2003). Seizure occurrence period (SOP): The period during which a seizure is expected after an alarm. An alarm without a following seizure onset during the SOP will be regarded as a false prediction. Seizure prediction horizon (SPH): The time interval between the alarm and the start of the SOP.
ip t
For a true prediction, after the alarm there is no seizure occurring during the SPH, but a seizure occurs during the SOP. of total seizure events.
us
False prediction rate (FPR): The number of false alarms per hour.
cr
Sensitivity: The measure is defined as the ratio of the number of true predictions to the number
Prediction time: The time interval of the true alarm prior to the onset of the predicted seizure.
an
3 Results 3.1 Overall performance
In this work, the SOP was set to 30 and 50 min, and the SPH was set to 2 min, which is an
M
appropriate value to enable, for instance, drugs administration or warning of the patient (Schelter B et al., 2006). Table 1 lists the results of performance evaluation of this proposed prediction method. The seizure prediction algorithm was tested on the EEG recordings of 21 patients in the
ed
Freiburg EEG database, containing a total of 66 preictal periods and 518.8 h of interictal data. On the whole, an average sensitivity of 86.95% and 89.33% with a same average FPR of 0.20/h were obtained, respectively for SOP = 30 min and 50 min. It can be noticed from Table 1 that there
pt
were 12 patients having 100% sensitivity for SOP = 30 min and 13 patients for SOP = 50 min. In other words, one seizure was missed for 9 or 8 patients with these SOP values. The FPRs of 10
Ac ce
patients were less than 0.15/h for SOP = 30 and 11 patients for SOP =50 min. Among them, four patients had all seizures predicted without false prediction for both SOPs. More than two thirds of all the patients, 15 patients, achieved average prediction time greater than twenty minutes for SOP = 30 min and thirty-five minutes for SOP=50 min. 3.2 Seizure predictability of HFD Investigating the characteristic changes in the preictal state with respect to the interictal state is significant for developing an effective seizure prediction algorithm. In this work, the HFD was used to estimate the changes in the dynamical characteristics underlying the EEG signals, through quantifying the degree of complexity and irregularity of the EEG recordings. As an example, Fig.1 presents two single-channel EEG recordings, which have relatively high and low HFD value respectively. Fig.2 presents the time evolution of the HFD of one-channel EEG recordings from patient 3 and patient 9 and illustrates the proposed seizure prediction method. In Fig.2, the EEG contains two seizures which have been correctly predicted by the proposed algorithm. Compared to the HFD of interictal EEG, the HFD shows a general increase in the preictal period for patient 3
Page 6 of 19
and a general drop for patient 9. To determine the significance of preictal changes of HFD for each patient, the HFD of preictal EEG
was
statistically
compared
with
those
extracted
in
interictal
period
using
one-way analysis of variance. A significance level of P<0.05 was set. Results indicated that 39.1% of the seizures studied in this paper showed a significant drop in HFD prior to the seizure, while 49.4% showed an increase in the preictal state. In the remaining 11.5% of the seizures, no significant changes were observed in the preictal state compared to the interictal state.
ip t
The proposed algorithm might miss the seizures in whose preictal state the changes of HFD
were unapparent or in an opposite direction relative to other seizures of the same patient.
cr
To illustrate this point, Table 2 provides the data of seven patients, including the mean values of interictal and preictal HFD, the minimum and maximum of the outputs of BLDA and Kalman
us
filter, the classification threshold, and the number of epochs classified as preictal by the algorithm. For a seizure event, if none of its preictal epochs were detected correctly, it was not predicted for both SOPs by our proposed algorithm. It can be noted from Table 2 that for six patients (excepting
an
patient 9), the missed seizures by the algorithm have the preictal HFD approximately equal to the interictal HFD or in an opposite change direction relative to other seizures. The missed seizure of patient 9 was maybe because of the heterogeneity of seizures.
M
3.3 Property of Kalman filter
Kalman filter was used to smooth the outputs of the BLDA classifier and to reduce the possible false seizure predictions. In order to assess the effectiveness of the Kalman filter, we
ed
additionally tested other smoothing strategies on the Freiburg long-term EEG recordings, such as moving average filter (Temko et al., 2011). The moving average filter is a low-pass filter used to reduce random noise. The value of a sample is calculated by the average of its neighbors: N 1 x m n 2 N 1 n N
pt y m
(5)
Ac ce
Where x is the input signal, y denotes the output (filtered) signal, and 2N+1 is the order of the filter, that is the number of points used in the moving average. Table 3 reports a comparison between the Kalman filter and the moving average filter. In
comparison to the moving average filter, Kalman filter produced a higher average sensitivity of 89.33% and also an improvement on the average prediction time, meanwhile the average false prediction rate reduced to 0.20 per hour. Besides, Kalman filter did not introduce an undesirable time delay for forecasting purposes like the moving average filter. 4 Discussion HFD was applied to characterize the nonlinear dynamical properties underlying EEG signals for seizure prediction in this work. HFD measures the complexity and irregularity of the EEG signals. The observed changes in HFD are essential for a prediction technique to distinguish between interictal and preictal epochs, regardless of the dimensions drop or increase in the preictal
Page 7 of 19
state. The measure of HFD can serve as a precursor of ictal activity for epileptic patients, and the ability to predict epileptic seizures prior to their occurrences was well illustrated by the results in Table1. 39.1% of the seizures studied in this paper showed a significant drop in HFD prior to the seizure, 49.4% showed an increase and the remaining 11.5% of the seizures had no significant changes observed in the preictal state compared to the interictal state. These results were concordant with the findings of Aarabi et al. (2012) who studied the correlation dimension and reported that 44.9% of the seizures they analyzed showed a significant drop in the values of the
ip t
correlation dimension in the preictal state, 44.9% showed an increase and 10.2% of the seizures were no significant changes observed between the preictal and interictal states. Compared with the
cr
correlation dimension widely used in seizure prediction study, HFD has no need of reconstructing
phase space for EEG signals and can be calculated directly in the time domain. Hence it has a low
us
computational burden which is suitable for real-time prediction.
To lower the amount of calculation, we only chosen three focal channels of 21 patients for seizure prediction in this work and the EEG recordings of three extra-focal channels were not
an
analyzed. Many studies have reported significant preictal changes of univariate nonlinear measures in channels corresponding to the epileptogenic zone (Lehnertz and Elger, 1998; Maiwald of seizure prediction in this work.
M
et al., 2004), and restricting the analysis to three focal channels was without loss in performance Table 4 presents a comparison on the prediction results between the proposed method in this paper and some other methods for SOP = 30 min and 50 min. All the prediction methods in Table
ed
4 were evaluated on the same Freiburg EEG dataset. Winterhalder et al. introduced a performance assessment criterion of seizure prediction methods and investigated changes in the dynamical similarity index (2003) and the synchronization (2006) of the preictal states. Maiwald et al. (2004)
pt
compared the performance of three different patient specific prediction methods including the dynamic similarity index, the accumulated energy, and the effective correlation dimension.
Ac ce
Feldwisch-Drentrup et al. (2010) developed two combinational seizure prediction system by using logical ‘‘AND’’ and ‘‘OR’’ combinations on the mean phase coherence and the dynamic similarity index of EEG recordings. Though these methods had a low FPR, the sensitivity made them insufficient for clinical applications. Compared to these previous works, the presented method improved the sensitivity.
In the recent work of Aarabi et al. (2012), five univariate measures including correlation
dimension, correlation entropy, noise level, Lempel-Ziv complexity, and largest Lyapunov exponent as well as one bivariate measure, nonlinear interdependence, were extracted from intracranial EEG recordings. The spatio-temporal information was then integrated by using patient-specific rules established from the characteristics changes of the preictal states to obtain maximum sensitivity and specificity for seizure prediction. They reported an average sensitivity of 79.9% with an average FPR of 0.17/h for SOP of 30 min and SPH of 10s. For SOP of 50 min, a higher average sensitivity (90.2%) and lower average FPR (0.11/h) were obtained. Compared with
Page 8 of 19
Aarabi’s system, our proposed algorithm provided comparable sensitivity, albeit with a little higher false positive rate. It is worth mentioning that the performance results of Aarabi’s prediction system were calculated only on 11 patients of the EEG dataset, whereas our proposed method was tested on the EEG data of all the 21 patients. In addition, the seizure prediction method of Aarabi needed to extract six features of EEG recordings, which would yield very high computational burden and make the method unrealizable in a real-time prediction system. In contrast, the seizure prediction algorithm presented here, combining the HFD with the BLDA, has
ip t
lower computational complexity because the calculation of HFD does not require phase space reconstruction and the BLDA’s parameter estimation avoids time consuming cross-validation.
cr
To evaluate the performance of our seizure prediction algorithm, we also compared the sensitivities with the ones obtained using the random and the periodical prediction methods as
us
described by Winterhalder (2003). In terms of the sensitivity, our method performed better than the random prediction method (9.52% and 15.35% for SOP = 30 and 50 min, respectively) and the periodical prediction method (10% and 16.67% for SOP = 30 and 50 min, respectively).
an
Although the results are satisfactory compared to that of the other methods in the literature, the proposed seizure prediction algorithm still has a long way from clinical applicability yet, especially, the average FPR of 0.2/h is a little high. The patients in the EEG database used in the
M
work were during monitoring for presurgical evaluation and had a mean seizure frequency of about 3 seizures per day which was higher than that of normal conditions. The proposed seizure
5 Conclusion
ed
prediction algorithm would be improved in the FPR and tested in clinical practice.
In this work, we have investigated the changes of HFD in the preictal state and developed a
pt
patient-specific seizure prediction algorithm with low computational complexity. The results demonstrate that the changes of HFD can serve as a precursor of ictal activities and the proposed
Ac ce
seizure prediction algorithm, by combining HFD features with the patient specific BLDA classifier, achieved high sensitivity and lower false prediction rate when tested on the Freiburg data set. The selection of HFD feature and BLDA classifier makes the prediction have a low computational complexity, and the utility of Kalman filter introduces no further time delay when it is used to smooth the outputs of BLDA. All of these characteristics make the proposed algorithm potentially suitable for real-time seizure prediction. Acknowledgment This work was supported by the Key Program of Natural Science Foundation of Shandong Province (No. ZR2013FZ002), the Program of Science and Technology of Suzhou under Grant ZXY2013030 and the Independent Innovation Foundation of Shandong University under Grant 2012DX008.
Appendix A. Higuchi fractal dimension (HFD)
Page 9 of 19
Given a one dimensional time series x = {x[1], x[2], ..., x[N]}, the algorithm to compute the HFD can be described as follows (Higuchi, 1988 ). The first step is to form k new time series xmk , which is defined by N m xmk x m , x m k , x m 2k ,, x m int k , k
(A.1)
for m=1,2, ..., k. Where m and k are integers, representing the initial time value and the discrete int(r) is the integer part of a real number r. The length of each new time series is then calculated as
cr
int N k m 1 N 1 Lm k x m ik x m i 1 k , N m k i 1 int k k
ip t
time interval between points, respectively. N is the total length of the original time series x, and
us
N 1 is a normalization factor. N m k int k
an
Where
(A.2)
Thus the average length of all time series having the same interval k, is computed as 1 k Lm k . k m 1
M
Lk
The average length L (k) is proportional to k
–D
(A.3)
, where D is the fractal dimension by
ed
Higuchi's method. Therefore, the procedure described above repeats kmax times for k =1,2, ..., kmax, and then the least-square fitting method is used to determine a straight line that best fits the given
pt
data set of {ln(L(k)), ln(1/k)} (k =1,2, ..., kmax,). The straight line determined by the least squares
Ac ce
method is the best-fit curve that has the minimal sum of the deviations squared (least square error) from the given set of data. The slope of the straight line is the Higuchi fractal dimension of time series x. In order to choose an appropriate value of the parameter kmax, HFD values were plotted against a range of kmax. The point at which the FD plateaus is considered a saturation point and that kmax value should be selected (Gómez et,al.,2009). Appendix B. Bayesian linear discriminant analysis (BLDA) The following is a brief description of the key ideas of Bayesian linear discriminant analysis (Hoffmann et al., 2008). As an extension of Fisher's linear discriminant analysis (FLDA), BLDA can be regarded as least squares regression in a Bayesian framework. Assuming that feature vector x and regression target, that is class label t satisfy t T x n ,
(B.1)
Page 10 of 19
Where n denotes the additive white Gaussian noise, and is the regression weight. Then the likelihood function of obey a Gaussian distribution, N
2 2 T p D , exp X t 2 2
(B.2)
Where X is a matrix with rows corresponding to training feature vectors, N dimensional column vector t contains the regression targets of all the training examples, D denotes the pair {X, t},
1
G
2 2 1 T p exp I 2 2 2
(B.3)
cr
Then the prior distribution for the regression weight is specified as
ip t
and represents the inverse variance of the noise.
Where is a very small value selected randomly.
(B.4)
an
0 0 0 0 I 0 0
us
Here G denotes the number of features, and I is a G+1 dimensional, diagonal, square matrix
p D , p
p D , p d
ed
p , , D
M
The posterior distribution of the weight can be obtained based on Bayes rule, (B.5)
Thus the linear discriminant function is given by y m T xˆ
(B.6)
pt
Where xˆ is a new input vector waiting for being classified, m is the mean of the posterior Gaussian distribution, which is computed as
Ac ce
m XX T I Xt 1
(B.7)
The parameters and which yield the optimal separating hyperplane can be selected with a simple iterative algorithm described by MacKay (1992).
References
Aarabi, A., He, B., 2012. A rule-based seizure prediction method for focal neocortical epilepsy. Clin. Neurophysiol. 123, 1111-1122. Bao, F.S., Lie, D.Y.C., Zhang, Y., 2008. A new approach to automated epileptic diagnosis using EEG and probabilistic neural network. International Conference on Tools with Artificial Intelligence. 2, 482-486. Chisci, L., Mavino, A., Perferi, G., Sciandrone, M., Anile, C., Colicchio, G., Fuggetta, F., 2010. Real-time epileptic seizure prediction using AR models and support vector machines. IEEE
Page 11 of 19
Trans. Biomed. Eng. 57(5), 1124-1132. Drongelen, W.v., Nayak, S., Frim, D.M., Kohrman, M.H., Towle, V.L., Lee, H.C., et al., 2003. Seizure anticipation in pediatric epilepsy: use of Kolmogorov entropy. Pediatr Neurol. 29(3), 207-213. Esteller, R., Vachtsevanos, G., Echauz, J., Litt, B., 2001. A comparison of waveform fractal dimension algorithms, IEEE Trans. Circuits-I. 48(2), 177-183. Feldwisch-Drentrup, H., Schelter, B., Jachan, M., Nawrath, J., Timmer, J., Schulze-Bonhage, A.,
ip t
2010. Joining the benefits: combining epileptic seizure prediction methods. Epilepsia 51(8), 1598-1606.
cr
Gómez, C., Mediavilla, A., Hornero, R., Abásolo, D., Fernández, A., 2009. Use of the Higuchi’s
fractal dimension for the analysis of MEG recordings from Alzheimer’s disease patients. Med.
us
Eng. Phys. 31 (3), 306–313.
Higuchi, T., 1988. Approach to an irregular time series on the basis of the fractal theory. Physica D. 31(2), 277-283.
an
Hoffmann, U., Vesin, J.M., Ebrahimi, T., Diserens, K., 2008. An efficient P300-based brain-computer interface for disabled subjects. J. Neurosci. Method 167(1), 115-125. Iasemidis, L.D., Sackellares, J.C., Zaveri, H.P., Williams, W.J., 1990. Phase space topography and
M
the Lyapunov exponent of electrocorticograms in partial seizure. Brain Topogr. 2(3), 187-201. Katz, M., 1988. Fractals and the analysis of waveforms. Comput. Biol. Med. 18, 145-156. Khoa, T.Q.D., Ha, V.Q., Toi, V.V., 2012. Higuchi Fractal Properties of Onset Epilepsy
ed
Electroencephalogram. Comput. Math. Method M. 2012, 461426. Kuhlmann, L., Freestone, D., Lai, A., et al., 2010. Patient-specific bivariate-synchrony-based seizure prediction for short prediction horizons. Epilepsy Res. 91, 214-231.
pt
Le Van Quyen M., Martinerie, J., Navarro, V., Boon, P., D'Havé, M., Adam, C., Renault, B., Varela, F., Baulac, M., 2001. Anticipation of epileptic seizures from standard EEG recordings. Lancet
Ac ce
357, 183-188.
Lehnertz, K., Elger, C.E., 1998. Can epileptic seizures be predicted? Evidence from nonlinear time series analysis of brain electrical activity. Phys. Rev. Lett. 80(22), 5019-5022. Lehnertz, K., Andrzejak, R.G., Arnhold, J., Kreuz, T., Mormann, F., et al., 2001. Nonlinear EEG Analysis in Epilepsy: Its Possible Use for Interictal Focus Localization, Seizure Anticipation, and Prevention. J Clin Neurophysiol. 18 (3), 209-222. Li, X.L., Ouyang, G.X., Richards, Douglas.A., 2007. Predictability analysis of absence seizures with permutation entropy. Epilepsy Res. 77, 70-74. MacKay, D.J.C., 1992. Bayesian interpolation. Neural Comput. 4, 415-447. Maiwald, T., Winterhalder, M., Aschenbrenner-Scheibe, R., Voss, H.U., Schulze-Bonhage, A., Timmer, J., 2004. Comparison of three nonlinear seizure prediction methods by means of the seizure prediction characteristic. Physica D. 194, 357-368. Mirowski, P., Madhavan, D., LeCun, Y., Kuzniecky, R., 2009. Classification of patterns of EEG
Page 12 of 19
synchronization for seizure prediction. Clin. Neurophysiol. 120, 1927-1940. Mormann, F., Andrzejak, R.G., Kreuz, T., Rieke, C., David, P., Elger, C.E., Lehnertz, K., 2003. Automated detection of a preseizure state based on a decrease in synchronization in intracranial electroencephalogram recordings from epilepsy patients. Phys. Rev. E 67, 021912. Mormann, F., Andrzejak, R.G., Elger, C.E., Lehnertz, K., 2007. Seizure prediction: the long and winding road. Brain 130, 314–333. Osterhage, H., Mormann, F., Staniek, M., Lehnertz, K., 2007. Measuring synchronization in the
ip t
epileptic brain: a comparison of different approaches. Int. J. Bifurcat. Chaos. 17(10), 3539-3544.
cr
Petrosian, A., 1995. Kolmogorov complexity of finite sequences and recognition of different preictal EEG patterns. Proc. IEEE Symp. Comp. Med. Syst. 212-217.
us
Schneider, M., Mustaro, P.N., Lima, C.A.M., 2009. Automatic recognition of epileptic seizure in EEG via support vector machine and dimension fractal. Proceedings of International Joint Conference on Neural Networks. 2841-2845.
an
Schelter, B., Winterhalder, M., Maiwald, T., Brandt, A., Schad, A., Timmer, J., Schulze-Bonhage, A., 2006. Do false predictions of seizures depend on the state of vigilance? A report from two seizure-prediction methods and proposed remedies. Epilepsia 47(12), 2058-2070. horizons. Epilepsy Res. 73, 213-217.
M
Schelter, B., Winterhalder, M.,et al., 2007. Seizure prediction: The impact of long prediction Temko, A., Thomas, E., Marnane, W., Lightbody, G., Boylan, G., 2011. EEG-based neonatal Winterhalder,
M.,
ed
seizure detection with Support Vector Machines. Clin. Neurophysiol. 122, 464-473. Maiwald,
T.,Voss,
H.U.,
Aschenbrenner-Scheibe,
R.,
Timmer,
J.,
Schulze-Bonhage, A., 2003. The seizure prediction characteristic: a general framework to
pt
assess and compare seizure prediction methods. Epilepsy Behav. 4, 318-325. Winterhalder, M., Schelter, B., Maiwald, T., Brandt, A., Schad, A., Schulze-Bonhage, A., Timmer,
Ac ce
J., 2006. Spatio-temporal patient–individual assessment of synchronization changes for epileptic seizure prediction. Clin. Neurophysiol. 117, 2399-2413.
Page 13 of 19
ip t cr us
Fig. 1 Two single-channel EEG recordings. The signals in (a) and (b) have relatively high and low
Ac ce
pt
ed
M
an
HFD values respectively.
Page 14 of 19
ip t us
cr Ac ce
pt
ed
M
an
(a)
(b)
Fig. 2 The HFD of one-channel EEG recordings from patient 3 (a) and patient 9 (b), and the illustration of the seizure prediction method. In the top left panels of (a) and (b), the EEG recordings respectively contain a seizure event with a preictal period of approximately 50 minutes. The two vertical lines indicate the seizure onset and offset respectively, and the seizure onset is labeled as time zero. The top right panels of (a) and (b) respectively illustrate 1h interictal EEG. The middle panels of (a) and (b) present the HFD time-traces of the EEG recordings given in the top panels. The bottom panels illustrate the proposed seizure prediction method. The blue curves and the green profiles respectively display the outputs of BLDA and Kalman filter. The classification thresholds are represented by red horizontal lines. For each patient, the threshold is a constant, that is, the red lines in both the left and right graphs of the bottom panels are at the same level. When the outputs of Kalman filter have at least two consecutive EEG epochs greater than the threshold, an alarm in blue column bar form for the impending seizure will be given.
Page 15 of 19
Table 1 Results for each patient using the proposed prediction method SOP=30min, SPH=2 min
Number of
Number
seizure
of true
events
alarms
1
3
2
SOP=50min, SPH=2 min Average
Number
prediction
of true
time(min)
alarms
0.25
19.78
100.00
0.33
3
75.00
4
4
5
4
6
Average
Sensitivity
FPR
(%)
(/h)
3
100.00
0.33
36.33
29.13
2
100.00
0.21
47.40
0.04
29.84
3
75.00
0.13
36.49
100.00
0.00
28.57
4
100.00
4
100.00
0.46
24.72
4
100.00
2
1
50.00
0.21
29.93
2
100.00
7
2
2
100.00
0.00
16.00
2
8
1
1
100.00
0.04
29.73
1
9
4
3
75.00
0.27
16.69
3
cr
Patient
10
4
4
100.00
0.24
26.07
11
3
2
67.00
0.12
16.10
12
3
2
67.00
0.04
29.30
13
1
1
100.00
0.00
14
3
3
100.00
15
3
2
16
4
17 18
FPR
(%)
(/h)
3
100.00
2
2
3
4
4
prediction time(min)
ip t
Sensitivity
0.00
47.12
0.33
49.72
0.42
33.63
0.00
26.27
100.00
0.36
49.60
75.00
0.27
34.54
4
100.00
0.28
41.47
2
67.00
0.12
49.93
3
100.00
0.12
19.71
12.73
1
100.00
0.00
49.87
0.40
28.49
3
100.00
0.68
39.18
67.00
0.13
24.00
2
67.00
0.08
45.00
4
100.00
0.33
27.97
4
100.00
0.38
49.64
4
4
100.00
0.00
28.48
4
100.00
0.00
43.20
4
3
75.00
0.12
27.98
3
75.00
0.00
43.38
19
3
3
100.00
0.38
27.40
2
67.00
0.08
17.07
20
4
3
75.00
0.58
24.51
3
75.00
0.42
34.16
21
4
3
75.00
0.36
16.44
3
75.00
0.08
33.51
86.95
0.20
24.47
89.33
0.20
39.39
an
M
ed
pt
Ac ce
Average
us
100.00
Page 16 of 19
Table 2 The data of seven patients having one seizure missed prediction for both SOPs
15
18
20
1.0663
1.0425
1.0750
1.0292
1.0467
minimum
maximum
1
1.0774
-1.9164
2.1352
-0.2972
0.7435
2
2
1.0764
-2.3967
1.7881
-1.0182
0.4676
0
3
1.1314
-3.6602
4.2210
-2.6464
1.9088
165
4
1.1822
-4.9540
4.3037
-0.9915
2.0962
276
1
1.0791
-3.7769
3.9563
-1.0265
1.7531
0
2
1.0843
-3.9209
3.7343
-2.4761
1.8906
2
3
1.0661
-1.6516
3.7225
1.1045
2.6769
449
4
1.0766
-4.3212
4.0349
-0.5685
1.8849
9
1
1.0777
-1.6283
1.3617
-1.0683
0.6245
0
2
1.0564
-1.3398
2.1196
-0.3303
1.1204
2
3
1.0639
-1.8746
1.9509
-0.6656
1.1517
8
1
1.0430
-1.3812
2.5734
-0.4781
1.5258
0
2
1.0563
-0.9429
4.2652
-0.3772
2.5456
28
3
1.0537
-0.6365
4.8526
0.7022
3.2622
78
1
1.1962
-1.2933
8.5827
-1.2989
5.2133
263
2
1.1796
-0.8237
4.0800
-0.7299
2.3789
48
3
1.0700
-1.0790
-0.2776
-0.8911
-0.3715
0
4
1.2097
-0.8523
3.5874
-0.6217
1.8732
72
1
1.0344
-1.4375
4.5179
-0.3503
1.7115
22
2
1.0270
-1.9806
2.9048
-1.0831
0.5910
0
1.8153
0.9182
2.4942
1.6382
1.3734
1.7369
Ac ce
ip t
maximum
Outputs of BLDA
cr
1.1136
minimum
3
1.0323
-1.9291
3.5050
-1.1559
1.4498
51
4
1.0378
-1.4811
8.0002
0.0076
2.0560
146
1
1.0470
-1.2382
4.5815
-0.5447
1.8642
2
2
1.0500
-1.1504
2.7341
-0.4145
0.3564
0
3
1.0602
-0.5575
4.6726
0.3230
1.7787
4
4
1.0577
-1.5914
2.7027
-0.3947
2.0336
18
pt
21
0.7375
number of epochs classified as preictal
mean of preictal HFD
us
11
1.0639
Outputs of Kalman filter
Seizure number
an
9
threshold
M
3
mean of interictal HFD
ed
Patient
Page 17 of 19
Table 3 Comparison between Kalman filter and moving average filter SOP=50min, SPH=2min Postprocessing method
Number of seizure events
Number of true alarms
Average Sensitivity (%)
Average
Average prediction
FPR (/h)
time(min)
66
57
87.38
0.26
30.33
Kalman filter
66
58
89.33
0.20
39.39
Ac ce
pt
ed
M
an
us
cr
ip t
Moving average filter
Page 18 of 19
Table 4 Comparison between the proposed method and other methods SOP=30min Year
feature
Sensitivity
FPR
Sensitivity
FPR
(%)
(/h)
(%)
(/h)
2003
dynamical similarity index
21-42
0.15
32-67
0.15
Maiwald
2004
dynamical similarity index
21-42
0.15
52
0.15
correlation dimension
13-30
0.15
46
0.15
accumulated energy
18-31
0.15
34
0.15
mean phase coherence, lag
60
0.15
70
0.15
35.2
0.15
45
43.2
0.15
79.9
Winterhalder
2006
synchronization index 2010
mean phase coherence or
Drentrup
dynamical similarity index mean phase coherence and
cr
Feldwisch-
dynamical similarity index Aarabi
2012
five
univariate
measures
and one bivariate measure Present
Higuchi fractal dimension
86.95
0.17
0.20
53
0.15 0.15
90.2
0.11
89.33
0.20
Ac ce
pt
ed
M
an
This work
ip t
Winterhalder
us
Authors
SOP=50min
Page 19 of 19