This item is only available as the following downloads:
1 EEG ANALYSIS OF BRAIN DYNAMICAL BEHAVIOR WITH APPLICATIONS IN EPILEPSY By JUI HONG CHIEN A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGRE E OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2011
2 2011 Jui Hong Chien
3 To my Mom, Dad and Oli
4 ACKNOWLEDGMENTS I would first like to thank Dr. Panos Pardalos for his guidance, wisdom, and his persistent drive towards excellence which helped me and numerous students to repeatedly achieve great accomplishment throughout the years. Most of all, I thank Dr. Deng Shan Shiau for sharing his brilliant ideas, time, energy, and research insights as well as the effort he exerted to assist me through a difficult portion of my graduate schooling. The benefits of his impact on my graduate experience are beyond measure. I am very grateful for the opportunity to have work ed with Dr. J. Chris Sackellares and for the large positive impact and clinica l guidance he has had on my research I would also like to thank Dr. Thomas Mareci Dr. Hans v an Oostrom, and Dr. Steven Roper for serving on my graduate committee and providing valuable feedback. I thank Optima Neuroscience Inc. Th ese studies would not ha ve be en conducted if Optima had not provided a valuable epileptic EEG database and computational machines I am also grateful for the support of Dr. Ryan Kern, Dr. Kevin Kelly as well as Dr. Jonathan Halford for their assistance with computation work and E EG recordings. I am grateful for my colleagues in the Center for Applied Optimization, Petros Xanthoupoulos, Jicong Zhang, Dmytro Korenkevych, and Siqian Shen. I also appreciate Dr. Mingzhou Ding for his guidance and support the beginning of my research e specially when I took his class and attended his group meetings. Special thanks go out Dr. Il ( Memming ) Park for his encouragement and urging about doing research. I also thank Dr. Sungho Oh and Dr. Rajasimhan Rajagovindan as well as Reixin Jiang for the ir companion ship and advice when I confronted challenges and felt alone I am also deeply indebted to Zhang Chong for providing me with her valuable time and support during her stay at UF Her genuine personality
5 strengthens my belief that beaut y consists of being smart, simple and sincere. G ood times with these friends, not limited to research ha ve become an unforgettable and integral part of my experience at UF. My achievement up to this moment would not have been possible without the support from all of my family member s First I want to thank to my mom for raising me and being supportive of all my decisions. My dad who was a tacit urn but brilliant father, is the main reason why I want to be capable of many things and pursue my Ph.D ab r oad I have onl y begun to realize how much he influenced me despite the limited time that we had together His brave and responsible attitude toward his business, entire family, and death influence d me deeply I am also grateful to m y uncles for help ing the family busine ss so that I can focus on my studies after the death of my dad I am also grateful to my aunties for caring for me so well since I was young so that I have become positive and confident enough to face challenges I am indebted to those who contributed to m y family business, the Chien Furniture. Finally my heartfelt thank s go to my intelligent and beautiful fiance, Oli who has no Ph.D. degree yet but taught me how to perceive my emotions and express my feelings to live healthily Her wisdom guided me to fu nction well and her insistent love induced me to become brave, mature and responsible in a relationship The quality of my li f e has been enriched because of her patience understanding and love. She sees those positive influences my dad passed o n to me and cherishes th o se elements in my personality Her willingness to cook and do the laundry let me fully dedicate myself to research. Besides patience, her willingness to be away from me when I need time and space to f ocus was also much appreciated.
6 TABLE OF CONTENTS page ACKNOWLEDGMENTS ................................ ................................ ................................ .. 4 LIST OF TABLES ................................ ................................ ................................ ............ 9 LIST OF FIGURES ................................ ................................ ................................ ........ 10 LIST OF ABBREVIATIONS ................................ ................................ ........................... 12 ABSTRACT ................................ ................................ ................................ ................... 13 CHAPTER 1 INTRODUCTION ................................ ................................ ................................ .... 15 Electroencephalogram ................................ ................................ ............................ 15 Epilepsy ................................ ................................ ................................ .................. 16 Objective and Contributions of the Dissertation ................................ ...................... 18 Organization of Chapters ................................ ................................ ........................ 20 2 THE PROPERTIES AND ANALYSIS OF NON STATIONARY TIME SERIES ....... 21 Foreword ................................ ................................ ................................ ................. 21 Properties of Time Series ................................ ................................ ....................... 23 Stationarity ................................ ................................ ................................ ....... 23 Local Stationarity ................................ ................................ .............................. 23 Parametric Methods for Analyzing Non Stationary Time Series ............................. 25 Non Parametric Methods for Analyzing Non Sta tionary Time Series ..................... 29 Wavelet Transforms ................................ ................................ ......................... 29 Local Stationary Wavelets ................................ ................................ ................ 30 Forecasting Non Stationary Time Series ................................ ................................ 33 Kalman Filter ................................ ................................ ................................ .... 34 Local Stationary Wavelet (LSW) ................................ ................................ ....... 37 3 SIGNAL REGULARITY BASED AUTOMATED SEIZURE PREDICTION ALGORITHM ON SCALP EEG ................................ ................................ ............... 39 Background and Significance ................................ ................................ ................. 39 Methods ................................ ................................ ................................ .................. 46 Seizure Warning Mechanism ................................ ................................ ............ 46 Overview ................................ ................................ ................................ .... 46 Pattern Match Regularity Statistic (PMRS) ................................ ................ 46 PMRS convergence (T index) ................................ ................................ .... 47 Selection of channel groups for seizure warning monitor ing ...................... 49 Detection of PMRS convergence ................................ ............................... 50 EEG Data Characteristics ................................ ................................ ................. 51
7 Subjects and EEG recording specifications ................................ ............... 51 Data selection ................................ ................................ ............................ 52 Statistical Evaluation ................................ ................................ ........................ 53 Estimation of performance statistics ................................ .......................... 53 Performance statistical validation comparison with a random seizure warning scheme ................................ ................................ ...................... 53 Cross validation on performance comparison ................................ ............ 55 Combining p values ................................ ................................ ................... 55 Results ................................ ................................ ................................ .................... 56 Training Results ................................ ................................ ............................... 56 Test Results ................................ ................................ ................................ ..... 57 Discussion ................................ ................................ ................................ .............. 58 Dataset Requirement for Performance Evaluation ................................ ........... 58 Patient Specific Optimization ................................ ................................ ............ 59 Vigilance States as a Possible Confounding factor ................................ .......... 60 Significance Test ................................ ................................ .............................. 61 4 PSYCOGENIC NON EPILEPTIC SEIZURE AND COMPLEX PARTIAL SEIZURE PATIENTS CLASSIFICATION ................................ ............................... 69 Preface ................................ ................................ ................................ ................... 69 Small World Network in EEG Data ................................ ................................ ... 70 Graph Theoretical Analysis on EEG Signals ................................ .................... 70 Preprocess ................................ ................................ ................................ 70 Revealing network structures in EEG data ................................ ................. 71 Results of Studies Using Small World Network Analysis to EEG Data ............ 76 Inter Hemispheric Power Asymmetry ................................ ............................... 79 Method s ................................ ................................ ................................ .................. 80 EEG Data Characteristics ................................ ................................ ................. 80 Subjects and EEG recording specifications ................................ ............... 80 Awake and relaxed state EEG data selection ................................ ............ 81 Functional Network Graph ................................ ................................ ................ 81 Inter Hemispheric Power Asymmetry ................................ ............................... 84 Results ................................ ................................ ................................ .................... 85 Network Measures ................................ ................................ ........................... 85 Inter Hemispheric Power Asymmetry ................................ ............................... 86 Discussion ................................ ................................ ................................ .............. 87 5 CONNECTIVITY TRANSITION FROM INTERICATL TO ICTAL STATES ON NEOCORTICAL EPILEPSY ................................ ................................ .................... 98 Preamble ................................ ................................ ................................ ................ 98 Material ................................ ................................ ................................ ................. 102 Methods ................................ ................................ ................................ ................ 104 Resul ts and Hypothesis Testing ................................ ................................ ........... 106 Discussion ................................ ................................ ................................ ............ 107
8 6 CONCLUSION ................................ ................................ ................................ ...... 116 LIST OF REFERENCES ................................ ................................ ............................. 122 BIOGRAPHICAL SKETCH ................................ ................................ .......................... 136
9 LIST OF TABLES Table page 3 1 The sensitivitie s, false positive rates and p values achieved by the proposed algorithm on training and test dataset.. ................................ ............................... 62 4 1 Small world network index, patients during awake and relaxed state. ................................ ........................... 90 4 2 Mann Whitney U test results of network measures of CPS and PNES patient groups during awake and relaxed state. ................................ ............................. 90 4 3 Mann Whitney U test results of relative power asymmetry of CPS and PNES patient groups during awake and relaxed state. ................................ ................. 91 4 4 sample T test results of relative power asymmetry of CPS and PNES patient groups during awake and relaxed state. ................................ ...... 91 4 5 P values of two way ANOVA (patient groups as fac tor) in each frequency band. ................................ ................................ ................................ .................. 91 4 6 P values of two way ANOVA (patient groups as factor) for each anatomically symmetric pair. ................................ ................................ ................................ ... 91
10 LIST OF FIGURES Figure page 3 1 Dynamic features of three EEG electrode signals ................................ .............. 63 3 2 The group T index plot with warning algorithm parameter indic ations ................ 64 3 3 Flow chart of the seizure warning mechanism ................................ .................... 65 3 4 Recording duration of all 71 subjects in the dataset ................................ ........... 66 3 5 One of training results from 20 trials. ................................ ................................ .. 67 3 6 The perfo rmance comparisons of the 20 th trail ................................ ................... 68 4 1 E ye closed awake and relaxed state EEG signals of a patient having CPS. ...... 92 4 2 E ye closed awake and relaxed state EEG signals of a patient having PNES ... 92 4 3 A weighted adjacency matrix ................................ ................................ .............. 93 4 4 An adjacency matrix after apply a threshold on the weighted adjacency matrix ................................ ................................ ................................ .................. 94 4 5 of individual anatomically symmetric pairs in the delta frequency band. .. 95 4 6 of individual anatomically symmetric pairs i n the theta frequency band. ................................ ................................ ................................ .................. 95 4 7 of individual anatomically symmetric pairs in the alpha frequency band. ................................ ................................ ................................ .................. 96 4 8 of individual anatomically symmetric pairs in the beta frequency band. ................................ ................................ ................................ .................. 96 4 9 of individual anatomically symmetric pairs in the gamma frequency band. ................................ ................................ ................................ .................. 97 5 1 Electrode placement of the neocortical epilepsy patient ................................ ... 110 5 2 of the segment with a seizure onset at the 120 minute ........................... 111 5 3 of the segment without any seizure activity ................................ ............. 112 5 4 Averaged degree over 300 calculation windows of nodes in the network ......... 113 5 5 The average of over 30 minutes before and after a seizure onset .......... 114
11 5 6 The average of over 30 minutes before an d after an imaginary seizure onset that a lso located at 120 minute point ................................ ...................... 115
12 LIST OF ABBREVIATIONS AED Anti epileptic drug ANOVA Analysis of variance AR Autoregressive ARMA Autoregressive moving average ARIMA Autoregressiv e integrated moving average CC Cross correlation CPS Complex partial seizure EEG Electroencephalogram EMU E pilepsy monitoring unit EWS Evolutionary wavelet spectrum FN False negative FP R False positive rate LSW Local stationary wavelet GABA Gamma Aminob utyric acid MSC Mean square coherence MTLE M esial temporal lobe epilepsy PLAI Phase lag index PMRS Pattern match regularity statistics PNES P sychogenic nonepileptic seizures STL max Short term maximum Lyapunov exponent SVM Support vector machine TN Tr ue negative TP True positive
13 Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy EEG ANALYSIS OF BRAIN DYNA MICAL BEHA VIOR WITH APPLICATIONS IN EPILEPSY By Jui Hong Chien May 201 1 Chair: Panagote M. Pardalos Major: Biomedical Engineering Electroencephalography (EEG) is a technology for measuring the activation of neurons and is used to investigate various pathological conditions of the brain Epilepsy is a common brain disorder that disables a patient with unforeseen seizures One salient characteristic of epileptic EEG data is the synchronous neuronal activation over a n excessive portion of a brain. My studies sought m ethodologies to define epilepsy related interaction between brain regions so that the development of epileptic activity can be monitored and warned or intervened Other findings in my studies may be used to assist epileptic foci localization and epileptic patient classification. T o forecast the coming of seizures features were extracted from EEG data Forecasting has been extensively and successfully studied for stationary time series However, the non stationarity of EEG signals tarnishes many convenient properties of a stationar y time series To overcome this complicat ion and achieve successful seizure warning s I used a signal regularity based dynamic feature and T index to monitor the entrainment process among brain areas utilizing EEG s recorded from ep ileptic patients. The hypothesis was that the regularity entrainment among certain brain sites precedes seizure onsets. An algorithm was proposed and implemented on preprocessed EEG
14 signals. The evaluation included a comparison between the proposed algorit hm and a nave random scheme. The combined p value over 20 cross validation trials showed that the proposed warning algorithm achieved a better performance ( p =0.015) than the proposed random warning scheme. A brain can be viewed as a complex network of ne urons. A brain functional network represent s quantitative interactions amongst EEG ch annels and can be expressed as a graph. Graph theoretical analysis therefore, can be applied to offer a broader scope to inspect the global functional network characteris tics of epileptic brain s and can reveal the existence of small world network structure. I further inspected the inter hemispheric power asymmetry of physiologically and psychologically epileptic brains and found significant differen ce s between the two pati ent groups The degrees of asymmetry of the two patient groups differ ed around the frontal lobe in the delta, theta, alpha and gamma frequency bands.
15 CHAPTER 1 INTRODUCTION Electroencephalogram Electroencephalogram ( EEG ) is a technology measuring the vo ltage of activation from a set of neurons within the brain EEG electrodes can record brai nwaves from either the scalp or, through invasive procedures, deeper layers of the brain tissue. Depending on the location and the type of the electrodes, EEG signal s can reveal different levels of neuronal activities Scalp EEG is the least invasive methodology and is therefore, widely used in many fields of neuroscience including psychology, epilepsy, brain machine interface and sleep research EEG has become a mor e popular means of recording interesting brain activity over the last few decades because improvement s in computing power and storage capacity make possible sophisticated analys e s of very large amount s of EEG data. Like electrocardiogram (ECG), another wi dely used bioelectrical signal in monitoring and diagnosing heart abnormalities, numerous electrodes can be applied on one subject so that the EEG data ha ve both temporal and spatial information of brain activity The location of electrode placement on a s calp usually follows the 10 20 system Distinct landmarks of a head are first identified and then electrodes are placed at 10% or 20% distance intervals along the landmarks. Through analyzing the multi channel EEG data of a short period of time, a temporar y brain network can be identified. Analyzing long term EEG recording data can reveal the dynamics of the brain network as the subject goes through many states of mental activity or pathological manifestation. The EEG data included in this dissertation con tain both scalp and intracranial EEG s All of them are continuous long term recordings from either Allegheny General
16 Hospital (AGH, Pittsburgh, PA) or the Medical University of South Carolina (MUSC, Charleston, SC) Depending on the study design long ter m continuous or extracted EEG recordings were analyzed. All scalp EEG data were recorded using 19 channels and intracranial EEG data used respectively different channels covering the different brain areas depending on the specific case of each epileptic pa tient. All recordings are from human beings who are older than 18 years old and experienced epileptic or psychological non epileptic seizures. Epilepsy E pilepsy is a brain disorder instead of a disease, characterized mostly by recurrent and unforeseen int erruptions of normal brain function. It has been known from ancient times and was at one point, attributed to divine intervention until the great Greek physician Hippocrates realized that it was a disorder of the brain (Iasemidis, 2003) It is also the second most common disorder of the central nervous system and affects about 0.4 1% of the population (Sander, 2003) The sudden interruption of brain function is termed a se izure. Many other physical or psychological sudden and severe events are also called seizures such as a heart seizure. The meaning of a seizure, going back to its origin in Greek, is An epileptic seizure is the manifestation of epilepsy a nd is due to abnormal ly excessive or synchronous neuronal activity in the brain (Fisher et al. 2005) The duration of the manifestation of an epileptic seizure is sometimes called an ictal period. Most epileptic sei zures can be recorded by the electrodes of an EEG and shown in the EEG data. Through reading and analyzing the epileptic EEG data, we gain more information about the onset pattern and type of epilepsy. During the onset of an epileptic seizure (ictal period ), one may observe that the synchronous and rhythmic discharges originate from
17 one part of the brain (partial, focal or localization related seizure s ) or begin simultaneously in both sides of the hemispheres (generalized seizures). F ocal seizure s after the onset may remain localized within one part of the brain or propagate to the other side of the hemisphere and cause a wider range of synchronous neuronal activity (second arily generalized seizures). If a seizure is confined to a limited area of the brain, it usually causes relatively mild and transient cognitive, psychic, sensory, motor or autonomic impairment. However, generalized seizures cause altered consciousness and accompany a variety of motor symptoms from the jerking of limbs to stiff or convulsive movements of the whole body. During the ictal period of focal seizures, consciousness may be unaffected (simple partial seizures) or be impaired (complex partial seizures). Few p atients with focal seizures feel unusual sensations (auras) when experiencing the beginning stage of an ictal period. Most seizures occur in a nearly unpredictable manner, stress ing the patients as well as the people a round them. Since an epileptic seizure is characterized by excessively synchronous neuronal activity, researchers hypothesize that the dysfunction of the inhibition mechanism of neurons favors the development of seizures (Buhl et al. 1996; Cossart et al. 2001) Gamma aminobutyric acid (GABA) is the most pr evalent inhibitory neurotransmitter in human brains. GABAergic inhibition is one of the main control mechanisms influencing the excitability of neurons. Many compounds that increase the GABAergic inhibition are used as anti epileptic drugs (AED). Although many studies confirm the relationship between GABAergic inhibition and the development of epilepsy, about 20 30% of epileptic patient s failed to benefit from the administration of AED (Kwan and Br odie, 2000; Shorvon, 1996) These latter patients expressed poor quality of life due to both
18 depression and neurotoxicity from AED (Gilliam and Kanner, 2002) Although patients have d ifferent seizure onset frequenci es and severities seizure types and frequencies did not correlate with the score of quality of life (Kanner, 2003) These study results suggest that the frequency or severity of seizures is not the primary cause of patient depress ion, but rather the sense of uncertainty regarding seizure onsets and the constant fear of unpredictable onset. The presence of correct seizure onset warnings should alleviate the degree of depression and increase quality of life by decr easi ng the amount of uncertainty and offering patients the ability to prepare for seizure onset. Objective and Contributions of the Dissertation My objective in pursuing a Ph.D. degree in biomedical engineering is to is to be able to offer engineering techniq ues to assist with clinical problems. I have focused my research on EEG signal feature extraction and epilepsy. One of the most intriguing problems in epileptic computational neuroscience is the predictability of epileptic seizures. For more than two decad es engineering techniques have been applied to predict seizure s Many brilliant ideas have been proposed and many promising results have been reported Although important gains in seizure prediction see m near there are in fact still many difficulties and questions that need to be solved. In addition criticism of previous over ly promising results have a risen forc ing us to re think the status quo challenges regarding seizure prediction. In 2007, M orman n et al proposed a more careful study design to corr ectly estimate the sensitivity and false positive rate of an automatic seizure prediction algorithm (Mormann et al. 2007) In my work a n automatic online warning algorithm using a new statistical feature was applie d to long term EEG scalp recordings. The performance was compared with a random predictor and achieved a p value of 0.1.
19 About 20% of epileptic patients cannot have complete seizure control through AED s and may consider brain surgery to achieve seizure con trol. A n interesting clinical challenge in epileptic brain surgery is to decide the size and location of the areas to dissect. There are many methods to find the potential resection zone including imaging techniques and long term video EEG recording. Most of t hese pre surgical evaluations are qualitative but not quantitative. A pilot study was designed to try identifying the epileptic brain network and quantifying the network features quantitatively. This study can be followed up by another study identify in g the critical nodes or regions of an epileptic network through more sophisticated network algorithms. The critical regions are those nodes that are the most important for maintaining the connectivity of the epileptic network. These quantitative results ma y be a quantitative reference for deciding the resection area of an epileptic surgery. Among epileptic patients, there is a subset of patients whose causes of seizures are not physiological but psychological. Patients having psychogenic non epileptic seizu res (PNES) need to be distinguished from other patients having epileptic seizures so that the treatment can be effective Patients have PNES onsets that can be identified by monitoring the ictal pattern shown in EEG and video. The classif ication using icta l symptoms requires several seizing onsets from a patient. The number or frequency of an EEG recording fac ility cannot be precisely anticipated in advance. I deally, one c ould identify PNES patients using only awake and interictal scalp EEG recordings. This study researched the connectivity features of awake and relaxed interictal EEG signals recorded from five patients having PNES and another five patients having complex partial seizures (CPS) with fixed foci.
20 Organ ization of Chapters The research studies presented in this dissertation are organized as follows : Chapter 1 provides an overview of the dissertation and introdu ces the background of my research. Chapter 2 begins with the introduction of non stationary time series and the engineering methodologies used to analy ze them Chapter 3 propose s a possible automatic scheme to warn epileptic patients to reduce the risk of having a seizure as well as lessen the mental stress of living with uncertainty and fear Chapte r 4 proposes another scheme to potentially categorize the PNES and CPS patients using only interictal scalp EEG recordings. Chapter 5 describes a study revealing an epileptic brain network. Chapter 6 summarize s the findings of this dissertation and provide s the direction for continued research
21 CHAPTER 2 THE PROPERTIES AND ANALYSIS OF NON STATIONARY TIME SERIES Foreword Most phenomen a happening with time can be measured and recorded as time series, such as wind speed, heart rate, stock exchange price and EEG signal s By analyzing a time series, we can gain insight into an array of interesting and complex phenomena A human brain is a complex system that responds simultaneously to many external inputs as well as numerous internal processes. Analysis and qu antification of the characteristics of EEG signals offer a window to observe the complex interactions amongst neurons in a brain. For example, in awake and asleep EEG signals one can recognize that the signals usually have swiftly changing backgr ounds Th erefore, EEG signals are usually regarded as non stationary time series consisting of the many distinct background activities that occur when a subject experiences different phases of psychological event s or physiological states Analyzing and forecasting non stationary time serie s have been challenging for researchers. Many possible methodologies have been developed to alleviate the difficulties of analyzing a non stationary process C hapter 2 introduces the basic properties and analysis of a non stationar y times series. Since the rapid growth in the power of computational machines, we now have the ability to investigate very data rich EEG signals with sophisticated and computation demanding techniques For an interesting time series, i t is convenient to a ssume that there are exact mathematical models generating the time series. Understanding the underlying mechanism of these processes is very crucial in engineering and science. It is also intuitive to presume that there are rules governing the route of a t ime series and once
22 we understand those rules we can better forecast how the time series may develop. However, real life phenomena are rarely that simple and usually involve many parameters noise and unobserved variables that can influence the time seri es with different time lags to a different extent. Still, with careful control of the recording environment, it is feasible to reduce many unobserved covariates and try to model the time series as a stochastic process The more precise we model the time se ries, the smaller forecasting error we can achieve To understand all the principles and rules controlling these black box systems at once is nearly impossible Therefore, I would like to introduce simple non stationary paradigms before moving to the more complicated time series. As aforementioned, a straightforward assumption is that a process possesses certain fixed properties in itself s o that it evolves with principle s One of the most convenient properties in a stochastic process is stationarity. Stati onarity in the strict sense means th os e statistics of a stochastic process do not change with time as described in Equation 2 1. (2 1 ) In Equation 2 1, denotes the probability mass function and denotes a time lag of any size. When a process fails to fulfill Eq uation 2 1, we claim that the process is not stationary in the strict sense. Stationary processes possess relatively stable stati s tical properties compared to non stati onary ones ; such as a fixed mean and variance. Nevertheless fixed mean and variance by themselves do not imply stationarity in the strict sense. W ith those fixed statistics we can begin predicting the foregoing process with more confidence. The first sect ion of C hapter 2 focuses on some necessary background for methods that will be
23 introduced later. In the second section of C hapter 2 we will present two well known practical methods to forecast non stationary processes. Properties of Time Series Stationa rity Stationarity in the strict sense is almost impossible in real life. Even a parametric simulated process can have stationarity only theoretically not to mention a time series recorded from real world signals which have much more complex mechanisms. Al though we can hardly have a stationary process in the strict sense, we want to utilize the concept of stationarity when studying non stationary time series because of its convenient statistical properties. W e can still claim that a process has stationarity in a wide sense if it has not only a fixed mean and variance but also time invariant auto covariance. One reason why we want to keep the concept of stationarity in a stochastic process is that our parametric model can hardly include all variables that aff ect the process, not to mention some factors that are beyond control. To take all effects into account, we should model a stochastic process with several deterministic and random terms. One way to implement this is to use the concept of pi ecewise stationarity. Even though the whole process is non stationary, we can still treat each segment of a reasonable length as if they were stationary, after disintegrating it into segments. Based on a piecewise stationary concept in the time domain, one can make an analogue to the time and frequency domain which is the idea of local stationary wavelet analysis which will be introduced in the non pa rametric method section. Local Stationarity One of the crucial parameters for using the concept of piecewis e stationarity is how to choose the length of a segment that is neither too short for calculating statistical
24 estimates nor too long so that the process is no longer stationary (Appel and Brandt, 1983) It should be noted that the criterion for selecting the segment length should always be dependent on the method that we are using in order to model the non stationary process. For example, for a quasi stationary process, the underlying assumption is that the parametri c properties might change their state between different time segments. As a result, an ideal segment should not be too long to overlap segments having different parametric values. Compared to quasi stationary processes, parameters of a time varying autoreg ressive model do not suddenly change but gradually evolve along with time. Proper selection of a method to model a process depends on some basic properties of the raw data. For instance, a seismic wave recording or an electroencephalograph recording conta ining a seizure shows a sudden state change. As a result, it is reasonable to use a quasi stationary method for modeling. For other processes such as the population of some species or the price of a stock, one might only be interested in their smooth evol ution property instead of their sudden change due to special or dramatic events, so it would be reasonable to use a time varying autoregressive model that will be introduced later. For more details of segmentation implementation, readers are referred to th e work of Appel U, Brandt A.V., 1983 (Appel and Brandt, 1983) and Adak S., 1998 (Adak, 1998) Those methods of modeling a non stationary process are called parametric methods In a ddition one can use non parametric methods which do not assume the existence of a model or distribution family in a process. Compared to parametric methods, non parametric methods still involve parameter choosing to some extent but not as much as parametr ic methods. If one intends to observe a non stationary process as a piecewise, quasi, or
25 evolutionary stationary process, a very well known and widely used non parametric tool called time dependent spectrum in the time frequency analysis field can be used. For a more detailed introduction of classes of non stationary process, readers are referred to the work of Adak S. ( 1998 ) (Adak, 1998) Parametric Methods for Analyzing Non Stationary Time Series The benefits of using a parametric method are its simplicity and lucid structure. The main feature of parametric methods is the use of parameters as structural elements to introduce interpretation and interaction so that one can reduce ambiguity within the data itself. This, howe ver, might be achieved at the risk of misinterpretation (Harrison and Stevens, 1976) Some other parametric approaches are based on specific assumptions that process outputs belong to the family or families of the distributions used in the model or can be described as the output of a stochastic process (Pereda et al. 2005) whereas non parametric approaches do not impose specific constraints with regards to the distribution fam ily. One has to first carefully scrutinize the statistical properties of the observed data and the necessary assumption for the parametric method being applied for estimation Schindler et al. 2007) In sum, parametric methods are beneficial only if there is a strong reason to look at the process in a particular structural way or as a member of a certain distribution family. As one can imagi ne, a parametric approach might not perform well on a process having a huge future random term that does not constantly fall in to any distribution family and has a weak link between current and future states of the process. Time Varying Autoregressive (AR) model s (West et al. 1999) just like a stationary process, can use parametric or non parametric means to start the analysis on a non stationary process. One of the most commonly used model s for simulating or fitting no n stationary time
26 series is the time varying AR model. A order linear time varying AR model can be written in a general form as Equation 2 2. ( 2 2) In Equation 2 2, are time varyin g parameters that differentiate a time varying AR model from an AR model that has a fixed coefficient instead of The main difficulty in using a parametric model for analysis lies in choosing proper values for th e parameters. As Eq uation 2 2 shows, and are parameters which need to be chosen based on the training data while is a zero mean white noise process. One advantage of a parametric method is that it renders a lucid dynamic relationship between data points and gives us some sense of how data evolve over time. However, the model may not be unique because of the uncertainty of parameters. One should bear in mind that we can fit models of differe nt orders into the same datasets. In this case the time varying AR model parameter values ( ) might be significantly different especially for a complex dataset. Several methods have been proposed to find the optimal order of an AR mode l such as Akaike information criterion (Akaike, 1973) Bayesian information criterion (Hannan and Quinn, 1979) and the likelihood ratio test (Vuong, 1989) Once the order of the AR model is decided, a re calculated from Yule Walker equations. O ther methods exist for estimating regression coefficients such as least squares approach, modified covariance m ethod covariance method, parametric spectral estimation method, Prony's method etc. However, choosing the order of a time varying AR model is more complex than that of an AR
27 model. A recently proposed method for time varying AR model order uncertainty is the generalized likelihood ratio test (Abramovich et al. 2005) Autoregressive Integrated Moving Average (ARIMA) model : W e model an interesting time series as a non stationary processes whose derivative, of pro per order, i s a stationary process. A n a utoregressive m oving a verage m odel (ARMA) can be represented in Eq u ation 2 3. ( 2 3) In Equation 2 3, is a Gaussian white noise process. A n ARIMA can exhibit not only ARMA but also non stationarity having neither fixed mean nor fixed variance. A very important property of ARIMA is that if one takes the derivative of the process it will result i n a stationary process th at can be represented by another ARMA model (Box and Jenkins, 1994) This relationship can be represented in Equation 2 4. ( 2 4) This relat ionship renders a powerful linkage for changing a non stationary process into a stationary process which can be more easily analyzed. However, one should notice that the high frequency components of the process would be amplified after differentiating (Huang and Chalabi, 1995) Considering an ARMA ( ), with it is possible that the AR part has more terms than the MA part. The realization of an ARM A state space model is g iven in Equation 2 5.
28 ( 2 5) ( 2 6) In Equation 2 5, and are both Gaussian white noise processes. Tilda indicates a vector. If we denote the transition matrix as the noise coupling matrix as and the observation matrix as = Using Eq uation 2 5 and Eq uation 2 6, Eq uation 2 3 and E q uation 2 4 become Equation 2 7 and Equation 2 8. ( 2 7) ( 2 8) We can make a comparison between ARIMA and time varying AR which has the dynamic matrix form as in Equation s 2 9 and 2 10. ( 2 9) (2 1 0) In Equation 2 9, is a zero mean innovation with time varying variance and is described in Equation 2 11. = ( 2 11)
29 Non Parametric Methods for Analyzing Non Stationary Time Series Non parametric methods assume no specific distribution or fixed structure existing for a process. Some parameter settings still exist in the non parametric methods but the parameters are mainly for the sensitivity aspects of the analysis rather than relating to whole process. As a result, non parametric methods are less affected by the choice of assum ptions or parameter value compared to a parametric method. In sum, non parametric methods do not explain the mechanism behind a time series and are more data driven. Wavelet Transforms Wavelet analysis is one of the most utilized tool s for time frequency a nalysis. Most of the stationary processes are generated and recorded in the time domain. However, a weak stationary process can be characterized by its frequency components. That is why Fourier transform is applied in almost every field of science and eng ineering. The Fourier transform represents a process formulating a su perposition of frequencies. ( 2 1 2 ) (2 13) In Equation 2 1 2 denotes an imaginary unit. Although the representation through Fourier transform can clearly show those weighted components of a process in terms of frequencies, it does not include any time component. One can hardly perceive when those frequency components participate in a process. Especially for real data such as seismic waves or biomedical signals, a functional tr ansform showing results in both the time and frequency domain can be much more informative and intuitive. A very widely used time frequency analysis is the wavelet transform which renders a perspective with
30 both time and frequency resolutions through impos ing a mask func tion before the original signal. ( 2 14 ) In Equation 2 14, is a window function. When is a rectangular window the above equation is also called short time Fourier transfo rm, windowed Fourier transform or time dependent Fourier transform. One may observe the characteristic of the wavelet from Eq uation 2 14 that it has a parameter that can characterize the signal at a certain range on the time domain. As we can imagine a real life signal such as seismic waves show characteristics of a short impulsive feature that can hardly be described through Fourier transform decomposing the seismic wave into a combination of sinusoidal functions ranging from negat ive infinity to infinity (Qian and Chen, 1999) Because of the dual property of a wavelet in both the time and frequency domain, one cannot use a wavelet for time frequency analysis only but also for de nosing a process that ha s noise with a certain unique frequency band and that occurs sporadic ally along the time frame (Taswell et al. 2000) Local Stationary Wavelets Assuming that a process possesses stationa rity, we could estimate second order statistics such as variance or auto covariance through only a segment of observation and have confidence in the estimated statistics because its stationarity has rendered characteristic quantitative homogeneity. As we h ave longer term observation, we can estimate statistics better through elongating our observational period. In contrast, for a non stationary process, if we assume that the variance of a process changes slowly over time, saying that variance is a continuou s and differentiable
31 function of time, the variance could still be estimated with certain precision by gathering information around the time point of interested variance. This demand exactly fits the property of wavelet transform which preserves the local characteristics of a process to some extent on both the time and frequency domain respectively A statistically slowly changing process can be called a local stationary process (Nason, 1999) For a covariance stationa ry process, t here is a Cramr representation as Equation 2 15. ( 2 1 5 ) In Equation 2 15, exp means exponential and is a stochastic process with ort honormal increments. Non stationary processes can be written almost the same as Eq uation 2 1 5 but is replaced by a function of time. Following this rationale, a non stationary process can be represented using local stationary wavelet and non decimated discrete wavelets where indicates the scale and indicates time location. To express a local stationary process, one can write it as Equation 2 16. ( 2 1 6 ) In Equation 2 16, is the mutually orthogonal zero mean random innovation. From Eq uation 2 1 5 and Eq uation 2 1 6 one should notice that in Eq uation 2 1 5 corresponds to in Eq uation 2 1 5 and both of them indicate the amplitude of each analysis component. One should also notice that is a replacement of the Fourier harmonics tem in Eq uation 2 15 The subscript indices and in represent scale and time location respectively just like those in The reason for
32 using a non decimated discrete wavelet is that it can be shifted to any time point and not only confined in shifts of compared with a discrete wavelet. Please note that non decimated wavelets do not have orthogonalities and are an overcomplete collections of shifted vectors (Nason, 1999) An example of a discrete non decimated wavelet is the Haar wavelet. Nason also proposed an evolutionary wavelet spectrum (EWS) to quantify how the size of changes over time as in Equation 2 17 (Nason, 1999) ( 2 1 7 ) Please note that in EWS, still has resolution on both the frequency and time domain but with a different time scale which rescales the whole observation time into R escaled time, is calculated by dividing k by the whole observation time length, and Our goal is t o estimate the variance of a local stationary process. To accomplish that, one way is to put together the localized information around the time point of interest, with the concept of local stationary wavelet and EWS, let represent loc alized autocovariance and be defined as Equation 2 18. ( 2 1 8 ) In Equation 2 18, is the autocorrelation function of and denotes the i nteger part of the real nu mber. ( 2 1 9 ) Once localized autocovariance is defined, localized variance at time can be defined as described in Equation 2 20.
33 ( 2 20 ) Further evidence supportin g that a time frequency domain autocovariance in terms of wavelets preserves the intuitive sense of an autocovariance in the time domain can be seen through the relationship between Allan variances (Allan, 1966) whic h is defined as Equation 2 21 and Haar wavelet variances (Nason, 1999) ( 2 2 1 ) It is worth noting that Allan variance involves merely localized time domain representation. However, a Haar wavelet involves both localized time and f requency domain representations. In this case, an Allan variance can be expressed as in Equation 2 22. ( 2 2 2 ) A Haar coefficient on the finest scale or smallest dilation, at translation or location can be expressed as Equation 2 23. ( 2 2 3 ) Then we can see that they have the relationship as Equation 2 24. ( 2 2 4 ) Forecasting Non Stationary Time Series Based on the background mentioned above, we are ready to present prediction algorithms for non stationary processes. A Kalman filter can be used for prediction on a more parametric basis. Compared to a Kalman filter, local stationary wavelets do not necessarily need to fit data with a parametric model and may concern only
34 measurement per se In contrast to a parametric analysis, a non parametric analysis forsakes the idea that there is a model governing the evolution of data and that in any case the observation should be allowed to characterize the data by i tself. Kalman Filter The Kalman filter is a recursive filter that can estimate linear dynamic parameters given some noisy measurements (Kalman, 1960; McGonigal and Ionescu, 1995) Here we shou ld not take the word, filter, as only a passive data processing algorithm but think of it as a n active computer program that gather s information from inputs and then optimally estimate s the transition of a system In the previous subsection introducing AR IMA we know that a non stationary process can be treated as a stationary ARMA after taking appropriate orders of derivative on the ARIMA (McGonigal and Ionescu, 1995) As a result, we can first analyze only ARMA wi thout loss of generality and then transform the result back to the ARIMA model. If an ARIMA fails to be differentiated properly into an ARMA, one can still implement a Kalman filter onto an ARIMA process even if the ARIMA has missing observations (Gomez and Maravall, 1994) Since the Gaussianity and linearity possessed in a Kalman filter are independent from the non stationarity of a process, we can assume linearity and Gaussianity in the Kalman filter and use it to anal yze a non stationarity process. Readers are referred to Gordon et al 1993 (Gordon et al. 1993) and Kitagawa, 1987 (Kitagawa et al. 1987) for more detailed methods for mode ling in the non Gaussian state space. According to work by McGonigal and Ionescu in 1995 a state space model must be found first before applying a Kalman filter for prediction (McGonigal and Ionescu, 1995) They pr method that can model an ARMA state space (Aoki and Havenner, 1991) In Huang et
35 al 1995 (Huang and Chalabi, 1995) t he authors implemented a Kalman filter to predict a time varying AR (2) process which is a non stationary and Gaussian Markovian process and can be described in a general form as Equation 2 25. ( 2 2 5 ) If one compares Eq uation 2 2 5 to Eq uation 2 2, one should notice that the time varying model mentioned in this section has an extra term, Huang added a time varying term, to better track the moving trend in the data of wind speed which i s one of the application s in his work. To describe this time varying AR (2) model in a matrix form, we first define a vector space containing parameters that we need to decide later and then represent the differences of parameters as they evolve with time u sing their dynamic relationship. ( 2 2 6 ) In Equation 2 26, is the identity matrix, is a zero mean white noise vector. If we define as an identity matrix, is a random walk process and then ends up being an integrated random walk process. If is a zero matrix, then should be a white Gaussian n oise process and ends up being a random walk process. Please note that the role of is the first order difference of and it describes the changes between and as time moves on. We can predict our non stationary process which is a time varying AR (2) model once the crucial step of choosing parameters in our model is done and this can be accomplished by Kalman algorithm. A general way of
36 implement ation of a Kalman filter starts with setting up the paramete r model in the state space form. ( 2 2 7 ) (2 28) In Equation 2 28, = is the state vector and is the observation vector of dimension meaning Given the estimation of as one can write the estimation dynamics in terms of a current estimation, a past estimation and a dynamic transform matrix A ( 2 2 9 ) T hen we define the covariance matrix of estimation errors state d isturbance matrix where denotes expectation, and the noise variance matrix of measurement The updated state vector estimation, recursive algorithm ( 2 30 ) ( 2 3 1 ) ( 2 3 2 )
37 Local S tationary Wavelet (LSW) Given a series of measurements, in which we want to use those measurements to predict the process step after we can write our interested process as where and th en define a linear predictor using Equation 2 33. ( 2 3 3 ) The next task is to optimize so that the mean square prediction error (MSPE), is minimized. To utilize the properties of LSW, we represented approximat e MSPE in a wavelet matrix form. ( 2 3 4 ) In Equation 2 34, and is the covariance matrix consisting of localiz ed autocovariance whose size is and th element is which is the localized autocovariance involving a wavelet form as mentioned in the background LSW section. The whole procedure of forecasting a non stationary process with long enough observations involves two training sets. First, divide the long term observations into two training sets and use the first one to generate the matrix. Second, use the other trai ning set to optimize vector such that MSPE is minimized. In practice, the first part involves selection of two more parameters and readers are referred to Fryzlewicz, et al 2003 (Fryzlewicz et al. 2003) for specific details. Once those parameters are selected, the mode l would be ready for prediction In C hapter 2 I introduced the properties of non stationary stochastic time series and some common methods for modeling or analyzing them. Foreca sting a non
38 stationary stochastic process can be achieved if the parametric model constructed is proper and robust or the non parametric methods are precisely updated For the remainder of my dissertation I will introduce a number of methods that are more applicable in the EEG analysis The application s of parametric models, such as AR models, and wavelet transform on EEG signals have provided some insights into seizure prediction (Chisci et al. 2010; Indiradevi et al. 2008; Le Van Quyen et al. 2001a)
39 CHAPTER 3 SIGNAL REGULARITY BASED AUTOMATED SEIZ URE PREDICTION ALGOR ITHM ON SCALP EEG Background and Significance An epileptic seizure is a transient occurrence of s igns and/or symptoms due to abnormal excessive and synchronous neuronal activity in the brain (Fisher et al. 2005) The worldwide prevalence of epilepsy ranges from 0.4% to 1% (Sander, 2003) Some epileptic patients experience a prodrome or an aura (Gupta et al. 1983) which can serve as a warning before the signs of seizure onset. Rare patients learn to abort a seizure without external i ntervention (Lee and No, 2005) However, most patients cannot predict or arrest their seizures. In industrialized countries, where antiepileptic drugs and seizure control devices are readily available, about 70% of epil epsy patients are able to gain satisfactory control of their seizures (Kwan and Brodie, 2000) For patients whose seizures do not respond to antiepileptic medications, less than 50% are candidates for epilepsy surgery (Engel Jr and Shewmon, 1993) Therefore, approximately 15 20% of epilepsy patients have no choice but to live their lives with unforeseen and uncontrolled seizure attacks, which cause considerable stress for these pa tients and their care givers and limit the range of daily activities available due to safety concerns. These lifestyle limitations decrease quality of life and may contribute to the increased prevalence of depression in patients with uncontrolled seizures (Kanner, 2003) If a device could be developed that could warn an epilepsy patient of an impending seizure, it could lessen the psychological stress of epilepsy and improve patient safety. Numerous studies on intracr anial EEG (Chaovalitwongse et al. 2005; Chavez et a l. 2003; Jacobs et al. 2009; Le Van Quyen et al. 2001b; Le Van Quyen et al. 2005;
40 Litt and Lehnertz, 2002; Mormann et al. 2003a; Mormann et al. 2003b) and fMRI (Federico et al. 2005) data in epilepsy patients suggest the existence of a preictal transition between an interictal and ictal state, and many attempts have been made to identify the most consistent patterns of preictal transitions In the last two decades, signal processing methods in nonlinear dynam ics have been popularly applied to EEG analyses b ecause of the non linear nature of neuron al activity and the paroxysmal character of epileptic seizure s (Hughes, 2008; Iasemidis et al. 1990; Lehnertz, 2008; Osterhage and Lehnertz, 2007; Osterhage et al. 2007) Due to the complexity of the implement of prediction, it is not enough to just utilize features showing differe nces between preictal and ictal states. Decisions about issuing a warning based on features need to be made continuously as the EEG signals shift from state to state in the real life environment. Therefore, many other sophisticated engineering techniques d erived from statistical learning and control theory were also employed on EEG data to facilitate seizure prediction performance (Chisci et al. 2010; Netoff et al. 2009; Tsak alis et al. 2006) S tudies utilizing a variety of engineering techniques and also the features on EEG data will be introduced It is not surprising that most seizure prediction studies were on intracranial EEG recordings because the signals are less infl uenced by movement, muscle, and other normal physiologic activities. Chisci, Mavino and coworkers applied classic parametric models and machine learning techniques to accomplish seizure prediction on intracranial EEG from nine patients in the Freiburg data base (Chisci et al. 2010) The a uthors extracted the parameters of autoregressive models for each EEG channel as dynamic features. The feature values were then fed into a classifier combining the support vector machi ne along with Kalman filter
41 techniques. Training and testing procedure s were implemented for 10 runs for each patient. Authors achieved 100% sensitivity and zero false positive rates on data from two patients while the remaining seven patients had false po sitive rate s among 0.12 to 1.03 false warnings per hour. In this study, the definition of a false positive is when a classifier calls an epoch as preictal and then classifies any of the following epochs as interictal before the coming of the next seizure. In contrast, the definition of a true positive is more intuitive as any epoch classified as preictal before a seizure in a segment containing a seizure. The average time difference between warnings and seizure onsets ranged from 6 to 92 minutes among patie nts. Netoff, Park and Parhi reported a study that utilized a cost sensitive support vector machine fed with the frequency power of EEG signal s to achieve seizure prediction (Netoff et al. 2009) They applied their me thod on 219 hour long EEG data for nine out of the 21 patients in the Freiburg dataset. The performance was evaluated using a prediction horizon of five minutes and achieved sensitivity of 77.8% and zero false positive rates. Mirowski, Madhavan and coworke rs applied several bivariate features and classifiers on the Freiburg database to recognize preictal patterns (Mirowski et al. 2009) Cross correlation, nonlinear interdependence, dynamical convergence and wavelet synchrony were estimated in five second windows. Features were then aggregated into spatio temporal patterns to feed into classifiers. Three classifiers were used and compared. The convolutional neural network outperformed logistic regression and the supp ort vector machine. The best result achieved by using wavelet coherence and a convolutional neural network was 71% sensitivity with zero false positive rates for 15 out of 21 patients. For those studies mentioned above they all applied classifiers on the
42 Freiburg data and report promising results. It would be intriguing to see if the classifier could work as well on scalp EEG data. In research conducted by Feldwisch Drentrup, Schelter and colleagues, two dynamic features (mean phase coherence and dynamic similarity index) were applied on 1,456 hours of long term continuous intracranial EEG data from eight patients (Feldwisch Drentrup et al. 2010) The prediction sensitivity increased (from 25% to 43.2%) w benefits of both features. However, the results of using individual methods and logical combination methods were all separately optimized within the same dataset. Therefore, the increase of sensitivity o f combined methods is not surprising. More other Studies regarding earlier seizure prediction methodologies and development using intracranial EEG have been extensively reviewed (Iasemidis, 2003; Litt and Lehnertz, 2002; Mormann et al. 2005; Mormann et al. 2006; Mormann et al. 2007) Scalp EEG has been used to investigate normal and pathological brain function for over 80 years, when Hans Berger published his work and coined the term (Swartz and Goldensohn, 1998) Early attempts to predict seizures from EEG signals began in the 1970s and flourished during the 19 90s (Lehnertz and Litt, 2005) Because of the non linear nature of neuronal function and the paroxysmal character of a seizure, non linear signal processing methods were popularly applied to EEG for predicting seizu res (Hughes, 2008; Iasemidis et al. 1990; Lehnertz, 2008; Osterhage and Lehnertz, 2007; Osterhage et al. 2007) By a pplying a non linear similarity index, Le Van Quyen, Martinerie and colleagues studied preictal EEG dynamics in 26 scalp EEG recordings obtained from 23 patients with temporal lobe epilepsy. In five patients with simultaneous scalp and intracranial recor dings, changes in
43 the similarity index values was observed during preictal period in both types of EEG recordings, and 25 out of the 26 EEG segments showed changes prior to the occurrence seizures (mean 7 minutes). Although this study did not evaluate the specificity of the similarity index change in long term EEG recordings (contained only 50 minutes before seizure onsets), it rendered an encouraging result using only scalp EEG to achieve seizure prediction (Le V an Quyen et al. 2001b) Another study by Hively and Protopopescu used L 1 distance and 2 statistic to estimate the dissimilarity in density functions between the base windows and the test windows in 20 scalp EEG recordings (Hively et al. 2000) The result showed pre seizure changes in all datasets with the forewarning times ranged from 10 to 13660 seconds. However, this study only examined one selected channel in each data set, and similar to (Le Van Quyen et al. 2001b) specificity of the method was not reported. In t wo follow up studies by the same group, one us ing all available recording channels and the other us ing a fixed channel, showed similar results (Hively and Protopopescu, 2003; Protopopescu et al. 2001) While the results seemed promising, no validation study has been reported for this method. In the study reported by Corsini, Shoker and c oworkers, 20 sets of simultaneous scalp and intracranial EEG recordings were analyzed using blind source separation and a short term Lyapunov exponent (STLmax) (Corsini et al. 2006; Iasemidis et al. 1990) This study reported changes of STLmax values before seizure onsets and noted that the scalp EEG may give better predictive power over intracranial EEG when the intracranial electrodes did not record the electrical activity in the epileptic focus. However, the practicality of this method on long term scalp EEG recordings may be limited due to the lack of an automatic procedure to select the most relevant source
44 component. Schad et al investigated seizure detection and prediction in 423 hour s of long term simultaneous scalp and intracranial EEG recordings from six epileptic patients (Schad et al. 2008) The method used techniques based on simulated leaky integrate and fire neurons. The study reported that 59%/50% of the 22 seizures were predicted using scalp/invasive EEGs given a maximum number of 0.15 false predictions per hour. In the study by Bruzzo et al ., a small sample of scalp EEG recordings (115 hours from 3 epileptic patients) was analyzed using p ermutation entropy (PE) (Bruzzo et al. 2008) By examining the area under the receiver operating characteristic (ROC) curve, they reported that the decrease of PE values was correlated with the occurrence of seizure s. However, the authors also concluded that the dependency of PE changes on the vigilance state may restrict its possible application for seizure prediction. More recently, Zandi, Dumont and colleagues reported a prediction method based on the positive ze ro crossing interval series (Zandi et al. 2009) The method was applied on a 21.5 hour scalp EEG dataset recorded from 4 patients with temporal lobe epilepsy. They reported a training result of 87.5% sensitivity (16 seizures) with a false prediction rate of 0.28 per hour, where the average prediction time was approximately 25 minutes. James and Gupta analyzed long term continuous scalp EEG recordings from nine patients (5 in training set and the other 4 in test datase t) (James and Gupta, 2009) The data were processed by a sequence of techniques consisting of independent component analysis, phase locking value, neuroscale, and Gaussian mixture model. The prediction performance of this method achieved a sensitivity of 65 100% and specificity of 65 80% as the pr ediction horizon ranged from 35 65 minutes in the test dataset.
45 Among these studies, one of the common features used in many of the seizure warning algorithms is the chan ge of synchronization of EEG signals recorded from different ways such as coherence (Zaveri et al. 1999) phase synchronizati on (Mormann et al. 2003b) (Iasemidis et al. 2005) In this study, we have attempted to id entify preictal transitions by detecting dynamic entrainment based on the convergence of PMRS among multiple electrode sites PMRS is a probabilistic statistic that quantifies the regularity of a time series It is especially useful when the moment statist ics (e.g., mean, variance, etc.) or frequency cannot detect changes in a signal By further applying a paired t convergence of two PMRS time series over time, we constructed an automated seizur e prediction algorithm that monitors the change of T index and issues a warning of an impending seizure when the T index curve exhibits the pattern defined by the algorithm. The general hypothesis is that seizures are preceded by PMRS entrainment; this hyp othesis was based upon findings reported previously using a different measure of signal order, STLmax (Iasemidis et al. 1990; Iasemidis et al. 2004; Sackella res et al. 2006) The prediction parameters and the specific EEG channels to be monitored were determined by the use of a training dataset Algorithm performance was then assessed using an independent test dataset. The performance was further validated b y comparison with that from a random warning scheme that did not use any information from the EEG signals
46 Methods Seizure Warning Mechanism Overview EEG signals were first filtered by a fifth order Butterworth filter with a band passing frequency between 1 to 20 Hz (the bandwidth within which most ictal epileptiform patterns occur) (Tsakalis et al. 2006) After the filtering process, for each channel, PMRS was calculated for each non overlapping 5.12 se cond epoch Based on the PMRS values T indices were then calculated for each of the selected channel groups To increase the sensitivity of seizure warning, the proposed algorithm independently monitored four T index curves (i.e., from four channel groups ). A warning was issued when any of the monitored T index curves met convergence criteria Pattern Match Regularity Statistic (PMRS) PMRS is a probabilistic statistic quantifying signal regularity (Kelly et al. 2010; Shiau, 2001) One of the characteristic features of EEG signals during a seizure is the rhythmic and regular discharges over a wide range of the brain. Therefore, the first step of the warning algorithm presented in this study calculated PMRS sequentially for each EEG signal analyzed. The rationale of applying this pattern match method (instead of value match) is due to its robustness over scalp signal values, which are usually more unstable than their up an d down trends The pr ocedure for calculating PMRS is described below: Given a time series with standard deviation a tolerance coefficient and a fixed integer the two segments in ( ) are considered pat tern matched to each other when Equation 3 1 is fulfilled.
47 (3 1) In Equation 3 1, t he first two criteria require value match to some extent at both the beginning and ending points of two segments, where was set to be 0.2 empirically. The third criterion requires pattern match between x i and x j with in a range of (set as 3 in this study). To calculate PMRS, we first d efine a conditional probability, ( 3 2) Given can be estimated as as in Equation 3 3. ( 3 (3 3 ) In Equation 3 3, Finally a PMRS can be estimated. (3 4) As the time series develops into a more regular state, s become larger and decreases as a result. PMRS convergence (T index) T index is basically the paired t statistic function used to quantify the degree of convergence between tw o PMRS time series. For two time series and (the PMRS
48 value time series), their values in a calculation window with a size of n data points are prese nted as and in Equation 3 5 and Equation 3 6. ( 3 5) ( 3 6) T hen the pair wise differences between and can be expressed in Equation 3 7. ( 3 7) The T index over the calculation window between the two time series is calculated using Equation 3 8. ( 3 8) In Equation 3 8, and are the sample mean and the sample standard deviation of respectively. PMRS convergence is a process during which one EEG signal is influenced by or coupled with another with respect to the signal regularity This phenomenon was used as a dynamical pattern for warnin g of an impending seizure Figure 3 1 shows the PMRS traces derived from three EEG channels (F8, T4, and T6) (upper panel) and their average T index values (bottom panel) over a 350 minute interval containing a seizure As show n in the PMRS plot, all PMRS values of the three channels dropped seconds after the seizure onset (indicated by a black vertical line), which was due to the extreme signal singularity during the ictal period. More importantly, the PMRS values became conver gent about 60 minutes before the seizure onset, which caused the decrease of T index values (shown in Figure 3 1 ).
49 Selection of channel groups for seizure warning monitoring Channel groups were selected such that they were bilaterally symmetric along the m idsagittal line. In order to avoid frequent eye movement artifact, electrodes Fp1 and Fp2 were excluded In addition, to avoid the regular alpha rhythm pattern that might give potential false positive warnings because of the nature of the regularity statis tics used in our prediction algorithm, electrodes O1 and O2 were also excluded. The frontal electrodes F3 and F4 were included in the analysis, whereas electrodes F7 and F8 were not included in order to minimize muscle artifact due to chewing. The four ch annel groups selected in this study resemble the well known and popular anterior These four channel groups were: (F7, T3, T5), (F3, C3, P3), (F4, C4, P4), and (F8, T4, T6). The main rationale for this selection was that most seizures occurring in patients with temporal lobe epilepsy start from a channel group on one side (hemisphere) of the brain, and therefore, intuitively, these channels were more likely to be entrained with each other during the preictal trans ition period. By monitoring four channel groups that cover both hemispheres, the algorithm was enabled to detect a PMRS convergence that preceded an impendi ng seizure initiated from either hemisphere For each channel group k k = 1 4, the group T index over the calculation window is denoted as: ( 3 9) The window size for calculating one group T index is 60 PMRS data points, with 59 point s overlapping from t to t+1 (i.e., sliding window). The warning algorithm only monitored the four group T index curves instead of individual pair wise T indices.
50 Detection of PMRS c onvergence The proposed prediction algorithm, instead of attempting to dete ct a signal threshold crossing, was set to detect a certain pattern of T index dynamics that gradually descends from a baseline value. Therefore, the first step to detect such a pattern was to determine an upper threshold (i.e., baseline value, see an illustration in Fig ure 3 2 ). A decrease of a group T index from U T was considered as a necessary condition for a potential PMRS convergence and a convergence was identified when the group T index values fel l below a lower threshold It was further defined that a period of convergence should be maintained for at least several minutes. For each group T index, its was set to be the asymptotic 95 th percentile in the preceding 12 minutes, as described in equation 10. The duration of 12 minutes was decided by training across multiple patients in the training dataset. ( 3 10) However, if half of the followi ng 16 T index values (representing 81.92 seconds of data) were greater than U T U T was updated as the median of these 16 values. This operation is described in equations 11 and 12. (3 11) (3 12) In Equation 3 12, is an indicator function, is the criterion of updating and is the group T index value at time If is f alse, then is kept as it is. denotes the Heaviside step function. is not updated with a maximum value because the existence of artifact in raw EEG recordings could affect the T index causing an abrupt
51 surge and t herefore produce an abnormally high Once is determined, is equal to minus where is a parameter that would determine the sensitivity of the algorithm In general, the bigger is, the less sensitive the al gorithm will be, but the less susceptible the algorithm will be to false warnings. In addition to the above detection rules, the period of a descending T index pattern from U T to L T must continue for more than minutes to be regarded as a PMRS convergence Therefore, the algorithm would issue a warning of an impending seizure only when any of the monitored T index curves traveled from U T to L T with the traveling time longer than to ensure a gradual descendent pattern. Once a conver gence was identified, other convergence s that followed (identified within the seizure warning horizon (SWH) by the warning algorithm using the same parameter settings) were silenced due to the possibility that a convergence preceding an onset may be much l onger than the value and cause several convergence identifications. The work flow is depicted in Figure 3 3 and gives an overview o f our seizure warning algorithm. EEG Data Characteristics Subjects and EEG recording specifications All subject s were 18 years of age or older admitted to either Allegheny General Hospital (AGH, Pittsburgh, PA) or the Medical University of South Carolina (MUSC, Charleston, SC) for inpatient seizure monitoring for diagnostic purposes or presurgical evaluation Data collection procedure w as approved by the Investigational Review Boards of AGH and MUSC and the Western Investigational Review Board (WIRB). The EEG recordings at MUSC were obtained using XLTEK monitoring systems (Oakville,
52 Ontario, Canada) with a sampling rate of 256 Hz and th e EEG recordings at AGH used 128 channel Nicolet BMSI 6000 systems (Viasys, Madison, WI, USA) with a 400 Hz sampling rate. The EEGs recorded at both institutions used a referential montage and the 19 electrode international 10 20 system of electrode placem ent The exact locations of referential electrodes placed in our dataset were decided on site and usually followed the recommended location of the Cz and Pz electrodes as suggested by the American Clinical Neurophysiology Society (Peters et al. 2001) All segments were reviewed by the collaborative clinical sites (AGH and MUSC) and all seizure events were verified by KK and JH, respectively Data s election Because a scalp EEG might be severely contaminated with arti fact typical muscle contraction and movement artifacts include blinking, chewing, and talking, only EEG segments with a tolerable level of artifact i.e., not present for more than 50% of the recording and not involving more than 50% of the recording cha nnels, were included in this study There were a total of 71 EEG recording segments from 71 patients selected for this study. All EEG segments were long term recordings (mean = 25 16 hours) containing at least one seizure. EEG segment s containing more than one seizure within any two hour interval were excluded to avoid potential overlap between preictal and post ictal periods of consecutive seizure s. Thus, the warning algorithm ha d a sufficiently long period of observation to detect the transition from inter ictal to ictal states The dataset was randomly divided into complementary training ( n=35 ) and test ( n=36 ) sets The total EEG recordings used in this study was 1786 hours and contained 98 seizures. The individual recording durations of each subject in the whole dataset (both training dataset and test dataset) are shown in Figure 3 4
53 Statistical Evaluation Estim ation of performance statistics The performance of the seizure warning algorithm varies with the length of the seizure warning horizon (SWH) the longer the SWH, the lower the sensitivity and false positive rate (FPR). In this work, SWH was defined as the time window following a seizure warning during which the patient was likely to have a seizure (Sacke llares et al. 2006) It is worth noting that Winterhalder et al (Winterhalder et al. 2003) defined the SWH (same as seizure prediction horizon, SPH) differently as a short intervention preparation period. Th one seizure onset occurred within the SWH following the warning. Otherwise, the warning was considered a f alse positive. After determining the outcome (i.e., true or false) of each warning, sensitivity was estimated by dividing the number of correct warnings by the total number of seizures in the EEG segments, i.e. the proportion of seizures correctly predict ed (within SWH). FPR (per hour) was calculated by dividing the total number of false positives by the total recording hours outside the SWH before each seizure. The SWH period before each seizure was excluded because it was impossible for a false warning t o occur within the SHW before each seizure. Performance s tatistical validation comparison with a random seizure warning scheme Several methods have been proposed for validating the performance of a seizure warning algorithm (Andrzejak et al. 2003; Kreuz et al. 2004; Schelter et al. 2006) Each study tested a specific statistical hypothesis. In this study, we compared the
54 prediction sensitivity estimates with a random prediction scheme that issued random seizure warnings (in time) during the recording. Under this random scheme, the interval between two consecutive warnings followed an exponential distribution, with a condition that the length of any interval coul d not be smaller than the SWH. This additional condition ensured that the two compared algorithms (test algorithm and prediction scheme) were consistent in terms of the management of nearby warnings (within the SWH). Although this condition diminished some rando mness from the random prediction scheme, it was necessary in order to have a meaningful comparison of prediction performance between the two prediction methods. In order to have an unbiased comparison for sensitivity, two parameters were controlled: 1) len gth of the SWH, and 2) total number of warnings For each of the test EEG segments, the total number of warnings allowed by the random predictor was set to be the same as that issued by the test algorithm. For example, if the test algorithm issued two true positive and two false warnings in a specific segment, the random predictor would only randomly i ssue a total of four warnings. Under the same number of warnings allowed, the study compared overall prediction sensitivity across all test patients A distri bution of overall sensitivity by the random prediction scheme was generated from 1000 simulations, and was used to assess how significant the prediction performance by the test algorithm was over the random scheme. It is worth mentioning that it was neces sary to impose the above conditions in the random warning algorithm compared in this study to make it comparable with the proposed algorithm. Since this was a conditioned random process, the specific null ensitivity (or FPR) of the test
55 algorithm is the same as that of a random prediction scheme with the same number of Cross v alidation on p erformance c omparison The comparison of performance between the proposed warnin g algorithm and the random scheme was further implemented through cross validation so that the comparison result is more independent from the constituent s in the training dataset and test dataset. T o execute the cross validate, w e randomly divided all 71 s egments into training and test datasets for 20 times to generate 20 individual trials and then repeated the training and test procedures on each trial In each trial, the parameters ( and ) were selected as those causing the performance to exceed 70% sensitivity and resulting in the least FPR in the training dataset. The parameters were input to the test algorithm run on test dataset to generate the final performance results in the trial. For each trial, the final performance result of the test dat aset was compared with the random seizure warning scheme mentioned above Please note that each test segment in each trial can have a different number of warnings and even the same segment in different trial would have a different number of warnings if the parameters used in the two trials were different. Using fixing parameters for a group of epileptic patients with different seizure types and physiological backgrounds can hardly outperform that using best setting parameters for each individual patient. Ho wever, these results can offer a perspective on the robustness of using the test algorithm among several patients when the previous individual patient data is not available or sufficient for training. Combining p values In each trial, the comparison result between the test algorithm and the random scheme is quantified as a p value, and indicates how significant the prediction
56 performance by the test algorithm was over the random scheme. After obtaining all the p values from the 20 cross validation trials, t he overall statistical significance of comparison between the two seizure prediction schemes was estimated by mapping the p values to a z distribution. The z transform method was chosen for testing combined p value s because it is not differentially sensiti ve to data that support or refute a common null hypothesis another widely used method for testing combined p values, is more sensitive to data refuting a common null hypothesis (Rice, 1990) Al l p values of each trial were transformed to a z value, in a N (0.1) distribution such that Pr[N(0,1) ] = p value. The overall z value, is carried out by calculating Last ly is transformed back to a combined p value. Results Training R e sults The proposed algorithm was applied on the training dataset (n = 35) to optimize settings for the parameters and Parameter D requests the minimal decrease that a group T index value needs to fall from to be considered as a PMRS con vergence, and gives the constraint on the minimal time length for a group T index curve traveling from to One of the training result s is presented in Figure 3 5 which shows the receiver operating characteristic (ROC) curves (sensi tivity vs. false positive rate) under different ranging from 10 to 3 0 min utes with 10 minute increments. Ea ch ROC curve contains eight settings of parameter each representing the performance under a certain value ranging from 1 to 8 (in T index units ) with increments of 1 As increased, each ROC curve first reac hed its peak sensitivity and then declined as kept increasing This was due to the constraint of : when was set to a small value
57 the t ime interval during which a group T index curve drop s from to was likely to be short er than the minimally required and therefore the drop of the group T index curve could not be considered as a PMRS convergence As D be came bigger, the sensitivity started increasing until the value became too large to start affecting the number of group T index drops (i.e., sensitivity started to reduce) For the training processes of all trials, the best performance was de fined as when the proposed algorithm achieves sensitivity of at least 70% with the lowest possible false positive rate The best parameter configuration was optimized within each trial according to the definition As a result, the optimized prediction para meter values were decided for each trail after review ing all performances using a full range of possible parameter settings. These optimized parameters were fixed for the performance evaluation in the paired test dataset within the same trail Test R esult s T he test dataset consisted of 36 long term EEG segments from 36 patients in each trial With the parameter settings determined in the training dataset of the same trial, the test algorithm gave sensitivities among 0.6 and 0.72 (mean = 0.65) while the fal se positive rate was among 0.22 and 0.26 (mean = 0.24) The random warning scheme was designed to issue random warnings with the same number as that issued by the test algorithm for each segment in the test dataset A fter the random scheme ran through all of the patients in the test dataset an overall sensitivity was calculated The random scheme repeated 1000 times on the test dataset and thus generated 1000 overall sensitivity estimates. The overall sensitivities of all 20 trials are listed in the T able 3 1. The overall sensitivity achieved by the test algorithm
58 was then compared with the sensitivity distribution of the random warning scheme. For example, the performance comparison in one of t he trials is shown in Figure 3 6 T he combined p value over 20 cross validation trials showed that the proposed warning algorithm achieved a better performance ( p =0.015) than the compared random warning scheme. Discussion Dataset Requirement for Performance Evaluation Amongst the early seizure prediction studies, shor t EEG recording containing seizure onsets were often used. The claim of sensitivity can be very satisfactory especially when the result was optimized for all the data without separating them into training and test dataset. Although these studies were absol utely valuable as pioneers in the field simply interpreting those results as successfully predicting seizures would be too optimistic. C oncerns of overestimation of results soon developed and drew much attention amongst many research groups. Using interic tal EEG recordings and a separate test dataset for performance evaluation was advocated to fairly estimate the false positive rate and avoid misleading conclusions (Lehnert z et al. 2007; Mormann et al. 2006; Mormann et al. 2007) In sum, potential misleading factors that could overrate the performance should be rule d out. These concerns should not be regarded as unduly skeptic al or critical of the rigor of these early st udy design s ; rather they represent the positive intention of expanding the previous promising results to more challenging test scenarios closer to reality. In this study, for the interictal recording requirement, we included 1786 hour long EEG segments in total, which included 98 seizures in 71 segments. Each seizure is at least three hours away from each other. If one defines a preictal period as three hours before the onset of a seizure, the dataset contains 1492
59 hours of interictal period. These interic tal EEG signals should be sufficient for a reasonable estimate of the false positive rate. For the separate testing dataset requirement, we randomly divided the database into training and test dataset for 20 trials. Every EEG segment was recorded from a di fferent patient and belongs either to the training or test dataset in a trial. This scheme of dataset treatment avoids both in sample and patient specific optimization. Patient Specific Optimization If one has a long enough EEG recording containing severa l seizures for each patient in the dataset, patient specific training and a test scheme can be implemented and sensitivity as well as specificity can still be fairly estimated. For many studies using classifiers to predict seizures, the classifier is usua lly optimized using a portion of the EEG segment of one patient and then test ed on the rest of EEG segment s from the same patient (Chisci et al. 2010; Mirowski et al. 2009; Netoff et al. 2009) Training parameters vary from method to method. As long as all parameters were fixed once after the training procedure ends, no more optimization or changing of parameter values should be allowed in the test process; otherwise, the result should not be viewed as an independent test result but a result of further in sample optimized procedure. Although the dataset is divided into two parts, the further optimized result should not be viewed as an independent test result since some para meters have fixed before the further optimization. In a recent study by Feldwisch Drentrup et al ., they divided their dataset into two sets (Feldwisch Drentrup et al. 2010) They executed channel selectio n in the first part of the data and optimized prediction intervention time and thresholds for two individual methods as well as two logical combinations of two methods in the remaining part of the dataset. The results of the four prediction methods were all in sample
60 optimized. In our study, we test ed our algorithm using the same parameters across many patients for one trial and predefine d a fixed prediction horizon for all patients. These test results are valid as long as the test dataset s have long enou gh intracranial segment s for estimating a false positive rate and enough seizure onsets for estimating sensitivity. The key point is to interpret the results with the design and the hypothesis of the study combined Our study is intended to test our algori thm across patients against a random predictor issuing the same number of warnings uniformly. Vigilance S tate s as a Possible C onfounding factor In addition to using interictal segments and separate test dataset s for performance evaluation, there are other concerns related to study design. One possible issue is the control of other confounding factors. For example, the scalp EEG signals respond to vigilance states more obviously than intracranial EEG signals. Chewing and talking muscle artifacts on the temp oral region electrodes during awake state and prevailing alpha waves on occipital region electrodes during eyes closed but awake state are very dominant in the scalp EEG signals. If the interictal segments were all recorded during the awake state while s eizures occurred during sleep, an algorithm detecting drowsy state can achieve a high sensitivity predicting seizures. One way to rule out of this erroneous relationship is selecting many interictal and preictal recordings under the same vigilance state fo r both training and test. Whether our algorithm is more sensitive to a vigilance state or an epileptic state change can be further clarified by conducting another continuous study. However, the existence of many different vigilance states in the data is ne cessary if one wants to test the robustness of a proposed algorithm under many vigilance states. In this study, all of the EEG segments are continuous long term EEG segments and most of them are long enough to go through several vigilance
61 states. A seizure could happen during the sleep or awake state s in a continuous recording. By including many continuous long term segments from many different patients, we scattered seizures among several vigilance states instead of one prevailing condition. Significance Test The result of this study should not be interpreted to mean that the performance of the proposed algorithm outperformed any random prediction scheme T h is study tested a specific null hypothesis based on the random prediction scheme designed In fact, a comparison of performance with a completely random prediction scheme may have very little meaning, both scientifically and practically. Certain conditions should be imposed on the random scheme to prevent potential bias in the comparison, especially for the assessment of a sophisticated EEG based prediction algorithm Nevertheless, there are several studies that have applied different random warning schemes or surrogate data to test different hypotheses (Andrzejak et al. 2003; Kreuz et al. 2004; Schad et al. 2008; Schelter et al. 2006) The comparison of the proposed algorithm with other random warning schemes should be followed Furthermore, studi es on a larger sample size and the assessment of the proposed algorithm with different SWHs should be conducted in the future.
62 Table 3 1. The sensitivities, false positive rates and p values achieved by the proposed algorithm on training and test dataset. The sensitivities and FPRs of the training datasets presented in this table are the best performances defined as reaching sensitivity over 70% meanwhile, result ed in the least FPR. The sensitivities and FPRs of the test datasets resulted from applying t he proposed algorithm using the best parameter configuration found using the training dataset. The p values in the test datasets w ere estimated by comparing the test algorithm to the random warning scheme (1000 simulations in each trial) that issued the sa me number of warnings. Training Dataset Test Dataset Trial Number Sensitivity FPR (times/hour) Sensitivity FPR (times/hour) 1 0.71 0.22 0.62 0.24 2 0.71 0.27 0.66 0.26 3 0.73 0.25 0.7 0.24 4 0.72 0.23 0.71 0.25 5 0.76 0.25 0.62 0.23 6 0.7 0.23 0.6 2 0.24 7 0.74 0.28 0.67 0.26 8 0.71 0.24 0.61 0.22 9 0.72 0.24 0.6 0.23 10 0.76 0.25 0.65 0.25 11 0.76 0.24 0.67 0.25 12 0.75 0.24 0.57 0.22 13 0.77 0.24 0.68 0.22 14 0.71 0.24 0.6 0.23 15 0.7 0.24 0.72 0.25 16 0.72 0.24 0.65 0.26 17 0.7 0.26 0. 7 0.24 18 0.72 0.23 0.61 0.24 19 0.72 0.24 0.67 0.23 20 0.72 0.24 0.7 0.25 Average 0.73 0.24 0.65 0.24
63 Fig ure 3 1. Dynamic features of three EEG electrode signals. Top: PMRS traces of F8, T4, and T6 electrodes. There is a sudden drop of PMRS values in all three channels right after the seizure onset (denoted by vertical dashed black lines in both panels at the 200 minute time point). Bottom: averaged T index among the three channels. Approaching the time of seizure occurrence, a gradual decrea se of the T index ( convergence ) from approximately the 120 min point to the 160 min point can be observed. The T index values remain small before the seizure.
64 Figure 3 2 The group T index plot with warning algorithm parameter indications. The seizure warning sensitivity relates to parameters and As the algorithm finds both and the difference of time and the group T index between them, as indicated by and in the figure, are compared with and respectively. Only the decrease d of a T index curve was larger than and also larger than would trigger a seizure warning. For example, with the proper setting of and the T index drop at 100 min was not considered a warning event because a large occurred during a short However, the T index drop from the 120 min point to the 160 min point shows a gradual and persistent decrease and is considered a proper warning event. The dashed vertical line at 200 minute denoted the seizure onset time point.
65 Figure 3 3 Flow chart of the seizure warning mechanism Calculating T indices Calculating PMRS Filter EEG signal Find next in each TSPG Calculating TSPG T indices Are 50% of the following 16 T indices bigger than ? No Find the time when T indices cross Yes Update Time that T indices fall from to Longer than ? No No Issue a warning Ye s Performance evaluation Recor d ing end? Yes
66 Figure 3 4. Recording duration of all 71 subjects in the dataset. The range of recording duration was between 6.18 to 70.25 hours. The mean was 25.16 hours and the standard deviation was 10. 98 hours.
67 Fig ure 3 5 One of the t raining result s from 20 trials Three ROC curves are shown under different values. Each point in the ROC curve corresponds to an overall result over the training dataset with a specific value. The inset indicates the parameter used in each ROC curve. For example, tt10 denotes the ROC curve using param eter = 10 minutes. The best pa rameter configuration of this training dataset was observed at the value equal to 20 min utes and equal to 3, which achieved sensitivity above 0.7 with a false positive rate of 0.232 per hour
68 Fig ure 3 6 The performance comparisons of the 20 th trial The figure shows the compari sons between the performances achieved by the random scheme and that achieved by the proposed algorithm. The bar plot shows the histogram of the sensitivities a chieved by the random warning scheme (generated from 1000 repetitions) T he sensitivity of the proposed algorithm applied on the same test dataset is denoted by the dashed vertical line which was better than 91.6% of the sensitivities of 1000 repetitions achieved by the random scheme.
69 CHAPTER 4 PSYCOGENIC NON EPILEPTIC SEIZURE AND COMPLEX PARTIAL SEIZURE PATIENT S CLASSIFICATION Preface About 25 % of the patients visiting an epilepsy monitoring unit are diagnosed as psychogenic nonepileptic seizure (PNES) p atients and need to be distinguished from other patients having epileptic seizures in order to be correctly treated. Patients have PNES onsets that can be identified by monitoring the ictal pattern shown either in scalp electroencephalogram (EEG) or video. The classif ication using ictal symptoms requires several seizing onsets from a patient. The number or frequency of seizure onsets during convenient for EMU schedule and cost e fficien cy if one c ould identify PNES patients using only awake and interictal scalp EEG recordings. This study researched the connectivity and features of awake and relaxed interictal EEG signals. The subjects include seven patients having PNES and another seven patients having complex partial seizures (CPS) with fixed foci. The hypothesis is that the inter hemispheric connectivity and network dynamic features from the EEG of a CPS patient are different from that of a PNES patient because the repetitive foc al discharges impair more around the epileptogenic zone rather than the contralateral hemisphere. Conversely, patients with PNES onsets should not have this abnormal physioelectrical consequence and are supposed to have a sounder functional network or mor e similar values of features between both hemispheres.
70 Small World Network in EEG Data A brain is a structurally and functionally complex network of neurons. The functional network reflect s the connectedness among brain regions in terms of neuronal activit y Graph theoretical analysis is a mathematical tool to reveal topological characteristics of a network. Applying graph theoretical analysis on the EEG data reveals the brain functional network features. One of the network structures called world n independence and global integration in the network. The balance can be evaluated by quantifying two graph features, called the local clustering coefficient and the characteristic path length. A small world network has a relatively high cluster coefficient and a small characteristic path length. Small world networks are usually compared to a network with a lattice like configuration or to a random network ; the two extreme cases. A regular latti ce network is characterized by a relatively high cluster coefficient and a long average path length. On the other hand, a random network has a relatively lower cluster coefficient and a shorter average path length (Boersma et al. ) Small world networks are efficient at information processing, cost effective, and are relatively resilient to network damage and as a result, can be regarded as the ideal model for a normally functioning brain network (de Haan et al. 2009) Graph Theoretical Analysis o n E EG Signals Preprocess While reviewing the EEG data, one can easily observe that EEG data has one salient characteristic, the rhythmic oscillation. Not surprisingly, therefore, frequency an alysis has been broadly used to preprocess and analyze EEG signals. It is widely accepted that the different strengths of the frequency components can reveal different
71 brain states. (de Haan et al. 2009; Dimitriadis et al. 2010; Gal et al. 2010; Palva et al. 2010; Ponten et al. 2007b) In analyzing EEG signals, filtering is usually implemented to remove artifacts or to help focus on the frequency bands of interest. Before doing a network analysis of EEG data it is important to remov e the artifacts that contribute to the multiplication of channels Without doing this first, the shared input from the artifact s would likely res ult in the existence of a spurious association amongst channels having the same artifact. Therefore, applying a power line frequency notch filter on EEG data is usually necessary and highly recommended. Any other artifact that would affect more than one ch annel, such as movement of the reference electrode, should be removed; otherwise, the result of graph should be carefully inspected to rule out the spurious coupling strength amongst vertices with common artifact sources. Some montages can be used to elimi nate possible common artifact s between two channels (Rubinov et al. 2009) Beside using montage, average reference is also commonly used when dealing with EEG data (Dimitriadis et al. 2010; Gal et al. 2010) In addition to using conventional fast Fourier transform based filtering techniques to obtain EEG signals within a range of frequencies, some studies implemented wavelet transform to extract c omponents with more specific time frequency resolution. For example, Palva et al (2010) used Morlet wavelets to filter EEG data into 36 frequency bands before construct the oscillatory phase synchronized network from EEG data (Palva et al. 2010) ). Revealing network structures in EEG data The brain functional network can be represented as a graph. G raph theoretical analysis is then applied after the graph has been established. To do so, I first had to
72 define vertices a nd edges in the EEG data. If the EEG channels are designated as the vertices of a graph, an edge between two vertices signifies a functional connection between two EEG channels. One would expect a larger correlation between two EEG channels when there is a n edge between them. Edges can also be values quantifying how well the two vertices correlate in weighted graphs. Most of the studies using graph theoretic analysis on EEG data assume that statistical interdependencies between EEG time series reflect funct ional interactions between neurons in the brain regions (Dimitriadis et al. 2009) There are many statistical metrics computing the degree of association between time series. Some of the most commonly used met rics in EEG functional connectivity graph studies will be discussed below Phase locking value (PLI ) is a statistic measuring the frequency specific synchronization between two neuroelectric signals (Lachaux et al. 1999) This is a method focusing on the phase information of time series and is different from coherence, which gives the interrelation of both amplitude an d phase between signals. The PLI is calculated by first extract ing the instantaneous phase informa tion of signals through either wavelet transform (Lachaux et al. 1999) or Hilbert transform (Mormann et al. 2000) Both methods lead to a similar result in real world EEG d ata (Quian Quiroga et al. 2002) When doing the wavelet transform on a signal x(t) a wavelet function, is first chosen. Gabor and Morlet wavelet function s have both been applied on EEG data (Lachaux et al. 1999; Le Van Quyen et al. 2001a) Then the wavelet coefficient time series, can be computed by convolute x(t) and (4 1) Then the phase time series, (t), can be computed.
73 (4 2) When doing the Hilbert transform, the phase information, (t), of a signal s (t) is obtained thr ough Equation 4 3. (4 3) (4 4) In Equation 4 4, denotes the Cauchy principal value in the equation. Then the relative phase, 1,1 (t), can be calculated as Equat ion 4 5. (4 5) The subscript n and m signify the relative phase relationship between channel e and f Finally, we can compute the PLI, PLV using Equation 4 6. ( 4 6) There are other computation al details including widowing and unwrapping the instantaneous phase being considered Interested readers can refer to those references mentioned above. PLI was used in a study by Dimitriadis et al to construct the functio nal connectivity graphs from 30 electrode EEG data. They utilized surrogate data to generate a baseline distribution of random PLIs and then determined the functional connections (edges) if there was a significantly different (p<0.001) PLI for a specific p air (Dimitriadis et al. 2010) The surrogate data is generated by permuting the order of trials of one signal repeatedly (Le Van Quyen et al. 2001a) Synchronizati on likelihood (SL) is a statistic measuring the non linear similarity between time series. SL offers an extended perspective of correlation that is not limited
74 in terms of linear relationship ; while coherence, has the limitation of rendering only the linea r correlation as a function of frequency (Stam and Van Dijk, 2002) The computation of SL first requires the interested time series, ( k M as channel number and i N as time indices), to be reconst ructed using the method of delays (Takens, 1981) (4 7) In Equation 4 7, denotes the lag and m is the embedding dimension. A number is then defined to denote the number of channels and that are closer than a crucial distance, (4 8) (4 9) (4 10) We can see that SL, is a measure describing how well channel k at time i is correlated to all o ther M 1 channels. Ponten et al used Synchronization Likelihood (SL) as the connectivity metric in their study about the relationship between epilepsy and small world network s (Ponten et al. 2007b) The edges betw een vertices were determined by applying a threshold on the SL value between two time series recorded at the two vertices within an epoch. In their study, for the threshold it was decided to keep the graph density as a constant. In other words, the thresh old increase s slightly
75 until the average number of edges of each vertex in a graph is equal to a certain number. Non linear I ndependence Measure (NIM) estimates the strength of functional coupling of two time series embedded in their respective phase space s (Arnhold et al. 1999; Quian Quiroga et al. 2002) Let and j w denote the time indices of the w nearest neighbors of two reconstructed vectors and in a phase space respectively For each the mean squared Euclidean distance to its w neighbors can be defined as Equat ion 4 11. (4 11) We can also define Y conditioned mean squared Euclidean distance as Equation 4 12 but replace the nearest neighbors by the equal time partners of the closest neighbors of (4 12) The NIM, N(X|Y) is computed as Equation 4 13. (4 13) The asymmetry of NI M (N(X|Y) ~=N(Y|X)) is the main advantage of this nonlinear measure because it can deliver directional information between vertices. In the study done by Dimitriadis et al ., they used both SL and NIM to quantify the extent of association between vertices (Dimitriadis et al. 2009) In their study, they found that the correlation dimension was around six from the sleeping EEG data of 10 healthy
76 subjects. Furthermore SL and NIM are weakly correlated and considered as revealing complementary information regarding the functional connectivity. Phase lag index (PLAI) quantifies the asymmetry of the distribution of instantaneous phase differences between two time series (Stam et al 2007) Suppose we have calculated the relative phase time series between channel e and f e f (t) then the PLAI can be defined as ( signifies the expectation value operator) PLAI ranges from 0 to 1. PLAI values greater than 0 suggest the existence of phase locking to some extent and value s equal to 0 signify n o coupling or no coupling with a phase difference centered around 0 radians. PLAI is supposed to rule out the synchronization due to instantaneous volume conduction or a common source that is the main cause giving spurious synchronization. Results of S t udies Using Small World Network Analysis to EEG Data The existence of small world structure in the cortex of a human brain has been observed through many measuring methods including fMRI (Achard et al. 2006; Salvador et al. 2005) magnetoencephalography (MEG) (Bassett et al. 2006; Stam, 2004) and EEG (Ferri et al. 2007) Many studies use grap h theoretical analysis to distinguish pathological consequence s or identify the state of human brains. The following paragraphs introduce summaries of some interesting studies related to small word network s in EEG data. Readers interested in more detail ab out the application of small world network s in EEG will find more information following the references. An ageing study found out that the clustering coefficient and the value of path length were both lower in the elderly compared to the young subject grou p. This implies
77 that the brains of the elderly subject group appear to be closer to random networks (Gal et al. 2010) Small world network s are supposed to be very efficient for data transfer. Epilepsy as a disease having excessive synchronization between neurons is assumed to have a relationship with the small world architecture in the functional brain network. Ponten et al conducted a study doing graph analysis on intracerebral EEG recordings from patients having drug resistant mesial temporal lobe epilepsy and found an increase in the clustering coefficient in the lower frequency band (1 13 Hz), and an increase in the path length in the alpha and theta bands during and after a seizure compared to interictal recor dings (Ponten et al. 2007b) This implies that the functional brain network changes from a more random organization to a small world structure. The efficiency of information transmission within a small world networ k implies that some pathological damage to a brain may result in losing the traits of a small world et al ., EEG data from patients with AD were compared with control subjects (de Haan et al. 2009) All subjects showed small world network traits in their functional graphs in all frequency bands, but in the beta band alone they observed significantly less small worldness in the AD group compared to the controls. Palva et al did a study using magnetoencephalogram (MEG) and EEG to investigate the network properties when the subjects are engaged in visual working memory (VWM) tasks (Palva et al. 2010) The study recor ded both MEG and EEG from 13 healthy subjects and applied PLV as the metric of association between brain regions. Their results implied that small world network structures appeared dynamically during
78 the VWM task execution within the alpha and beta band. M oreover, the small worldness was dependent on the VWM memory load. G raph theoretic analysis has also been applied to the quantitative EEG data of epilepsy patients. In a study by van Dellen, interictal electrocorticography (EcoG) of 27 patients with pharma co resistant temporal lobe epilepsy (TLE) were analyzed (Van Dellen et al. 2009) Because an epileptic seizure is a manifestation of overly hyper synchronous activity between neurons, a metric measuring synchroni zation, PLAI, was implemented on those electrode s located on the temporal cortex and showed lower values on those patients with a longer history of TLE. In addition, the cluster coefficient and small world index were both negatively correlated with TLE dur ation in the broad frequency band (0.5 48 Hz). This may have result ed from the accumulative damage on the brain tissue due to intermittent excessive epileptic discharges over a long period of time The optimal balanced small world structure had been impair ed to form a more randomized one as the TLE duration continued In an interesting model undertaken by Rothkegel and Lehnertz (Rothkegel and Lehnertz, 2009) they observed the co occurrence of local wave patterns an d global collective firing in a two dimentional small world network. Boersma et al conducted a study on resting state EEG in developing young brains (Boersma et al. ) They recorded resting state eyes closed EEG (14 chann els) from 227 childern when they were 5 and 7 years of age and found out that the clustering coefficient increased in the alpha band with age. Path lengths increased in all frequency bands with age. This suggests that a brain shifts from random towards mor e ordered, small world like configurations during
79 maturation. Girls showed higher mean clustering coefficients in the alpha and beta bands compared with boys. S chizophrenia has been suspected as the result of a more disconnected brain network among certai n crucial areas in the brain (Peled, 1999) Rubinov et al did a study investigating the disconnection hypothesis. They recorded resting state scalp EEG from 40 subjects with a recent first episode of schizophrenia a nd another 40 healthy matched controls. Nonlinear interdependences were identified from bipolar derivations of EEG data and weighted graphs were generated. Graphs of both groups showed features consistent with a small world topology but graphs in the schi zophrenia group displayed lower clustering and shorter path lengths. This result can be interpreted as a pathological process that the small world network transformed to a more randomized small world network in a schizophrenia brain. Th is randomization may be the reason why schizophrenia evidences cognitive and behavioral disturbances (Rubinov et al. 2009) In another similar study (Pachou et al. 2008) scalp EEG data from 20 s chizophrenics and 20 controls (age and sex matched) collected when they were performing working memory tasks were analyzed using conventional coherence. After applying a threshold on the values of coherence so that the number of average edges is five in a graph, binary graphs were obtained. The results showed that schizophrenics bear a more random neuronal network organization especially in the alpha band. Inter Hemispheric Power Asymmetry The symmetric parts of a human body such as limbs and sense organs often have the same function and structure. The brain of a human also has symmetric shape although some functions happen on a dominant side. The EEG signals from a pair of
80 symmetric channels usually have similar morphologies. The asymmetry of EEG signals is regarded as a pathological consequence or unusual phenomenon if the reason is not found. Many studies have used asymmetry as an index to quantify the morbidity of brain disease or abnormal states (Debener et al. 2000; Hagemann et al. 1999; Thatcher et al. 2001) In this study, the hypothesis is that the relative frequency powers of symmetric pairs of EEG channels are more different in CPS patients than that in PNE S patients. The powers of EEG signals in narrow frequency bands associate with different brain functions or motifs (Dimitriadis et al. 2010; Miltner et al. 1999) The degree of asymm etry is quantified through the relative frequency power of several frequency bands and the T index, a statistic measuring the degree of divergence between two groups. Methods EEG Data Characteristics Subjects and EEG r ecording s pecifications In this study I included seven PNES and seven CPS patients. All subject s were 18 years of age or older who were admitted to the Medical University of South Carolina (MUSC, Charleston, SC) for inpatient seizure monitoring for diagnostic purposes or presurgical evaluati on The data collection procedure was approved by the Investigational Review Boards MUSC and the Western Investigational Review Board (WIRB). The EEG recordings at MUSC were obtained using XLTEK monitoring systems (Oakville, Ontario, Canada) with a sampli ng rate of 256 Hz and the EEG recordings at AGH used 128 channel Nicolet BMSI 6000 systems (Viasys, Madison, WI, USA) with a 400 Hz sampling rate. The EEGs recorded at both institutions used a referential montage and the 19 electrode international 10 20 sy stem of electrode placement The
81 exact locations of referential electrodes placed in our dataset were decided on site and usually followed the recommended location of the Cz and Pz electrodes as suggested by the American Clinical Neurophysiology Society (Peters et al. 2001) Awake and relaxed state EEG data s election The EEG signals of each subject were reviewed to select out the awake and relaxed state sections containing the dominant alpha wave s over the occipita l regions and no eye blinking event s around the frontal electrodes. To classify patients using only interictal information, the awake and relaxed state EEG signals containing epileptiform discharges or other suspicious epileptic activities were also exclud ed. All selected awake and relaxed state sections were at least five hours before the first CPS or PNES happened. During the awake and relaxed state, the brain is supposed to be in a resting state and not actively involved in any goal oriented events. As a result, the awake and relaxed state s offer controlled background for the brain network to be compared between subjects. On the other hand, an awake and alert section of long term continuous EEG from different patients could contain a variety of ongoing ps ychological activit y, raising concerns of confounding if used for comparison analysis. Functional Network Graph Since oscillation frequency is a main characteristic of the brain and changes when the brain undergo es different psychological conditions or exe cute s different cognitive tasks, all selected epochs were filtered to specific frequency bands including delta, theta, alpha, and beta. All signals were filtered using filters that do not distort the phase information of the filtered signals so that the in stantaneous phase time series of the original and filtered signals are the same. This is crucial for this analysis because the connection strength was evaluated using phase information solely
82 PLAI was estimated in every five second epoch of all selected awake and relaxed state section s A weighted graph representing the functional network can be generated after all pair wise PLAI are calculated amongst all electrodes. The weighted graph can also be presented as a weighted adj acency matrix as in Figure 4 3 To convert the weighted graph to a binary one, a threshold can be chosen and applied to all PLAI values such that the edges between vertices is either connected or disconnected. The choice of threshold is done by controlling the density of a graph so tha t the non trivial structure of the network can be reveal ed I selected to present the network in a graph with density around 0.75 so that only one quarter of the strongest connections (larger PLAI) remained in the functional network graph and the other wea ker edges (smaller PLAI) are ignored. The threshold in each epoch can be different from each other and is increased slight ly from a small value until the desired density of a graph is reached. The threshold changes from epoch to epoch because the brain may undergo many phases of signal process ing and show different connection strength between functional regions while the structure of the information flow should persist in a small world configuration so that the information is efficiently shared and processe d among functional regions. For different subjects, it is reasonable to have individual thresholds for each epoch so that the dominant structure of the functional network can be revealed and compared across subjects. A network measure is a value quantifyin g a characteristic of the topology of a network. Several measures have been proposed and used in analysis of network structures. For determination of the existence of a small world configuration in a network, a clustering coefficient and minimum path lengt h are two crucial measures to
83 evaluate how small world a network is (Watts and Strogatz, 1998) A clustering coefficient quantifies how locally entangled a network is and a minimum path length reflects how globally integrated a network is. These two measures of a functional network graph are compared with those generated from 50 randomized graphs keeping the same number of vertices and edges. Given a clustering coefficient as Cp and a minimum path length as Lp for an epoch of selected awake and relaxed state EEG; the randomized graphs generate Cp s and Lp s respectively. Given is the number of edges between neighbors of vertex and is t he number of neighbors of vertex the local clustering coeffi cient can be calculated as Equation 4 14. In Equation 4 15 t he clustering coefficient, is the average of local clustering coefficients of every vertices in the graph as V is the number of vertices in the graph. (4 14) (4 15) (4 16) In Equation 4 16, is the shortest path length from vertex to vertex .The ratio of Cp/ < Cp s > to Lp/ < Lp s > is called the small world network index, (Humphries et al. 2006) If world network because the investigated network possesses a higher clustering coefficient or l ower minimum path length comparing to 100 randomized networks possessing the same number of edge s and vertices.
84 Inter Hemispheric Power Asymmetry Frequency power density, can be estimated by applying discrete Fourier transform on the interesting ti me series, (4 17) However, discrete Fourier transform is not a consistent estimator and modification is needed to better estimate frequency power density. For each five second epoch, w as divided into eight windows that overlapped half of the neighbor windows. For each window, hamming window technique was applied to reduce the noise due to truncation of signal and power density was estimated. All power density functions from eight window s were averaged as the estimate of the power density of the epoch. Relative powers of delta (1 4Hz), theta (4 7Hz), alpha (8 12Hz), beta (13 30Hz), and gamma (30 58Hz) band were computed from the power density function of frequency. For example, the alpha band relative power would be the ratio of the sum of powers in the alpha band to the sum of powers from 1 58Hz. These ratios, were further transformed to a variable, as described in Equation 4 18 so that has a distribution close to the n ormal distribution (Gasser et al. 1982; John et al. 1980) (4 18) The relative powers in each frequency band were compared to those of the contralateral hemispheric EE G signals by the computation of the T index. T index is a function to compute the degree of divergence between paired samples from two groups. In this case, the sample groups were the relative powers of a certain frequency band from an anatomically symmetr ic, left and right, pair. Left and right relative powers in each epoch
85 should be paired up because they both reflected the state of the brain during the same period of time. (4 19) In Equation 4 16, is the mean of the differences of relative powers, of each awake and relaxed state EEG epoch and is the standard deviation of the diffe rences. Due to the fact that each subject has a different length of awake and relaxed state EEG recording s, and the number of samples (degree of freedom), n for calculating T index must be fixed so that the comparison is meaningful, random sampling was pe rformed for each subject so that each subject has 36 (9 random samples from 4 continuous awake and relaxed state EEG segments) epochs input to the T index function. Results Network Measures For each five second epoch, one functional network graph was gene rated after the PLAIs were estimated pair by pair. For each functional network graph, the clustering coefficient, minimum path length, and small world index was calculated from the adjacency matrix. All network measures of every epoch from a patient were a veraged into one value for each network measure. As a result, every subject has one final value for each network measure. Totally, the PNES or CPS patients group has seven subjects and therefore, seven values for each network measure. The small world index is larger than 1 in all frequency bands for both patient groups support ing the existence of a small world network structure in both groups in all frequency bands. A S tudent T test
86 was performed to test if the in CPS or PNES patient group was larger than one. The p values showed significance in all frequency bands for both groups (see Table 4 1) Wilcoxon signed rank test was also performed and all subjects show ed p values as 0.0156. The Wilcoxon signed rank test results were the same for all frequen cy bands since the sample number is always seven and there were seven patients in both groups and all samples were larger than one. Following this I tested if these network measures show ed a difference between the CPS and PNES patient groups. A Mann Whit ney U test was performed to assess if the network measures were different enough between both groups and the p values were shown in the Table 4 2. Only the minimum path length ratio ( Lp/Lp s ) in the delta band showed a significant p value. Inter Hemispheri c Power Asymmetry The inter hemispheric power asymmetry was quantified by In the delta, theta, and alpha bands, the means of of every pair (besides T5 T6) in the CPS group were larger than that of the PNES group. In the beta band, the means of of every pair (besides P3 P4) in th e CPS group were larger than that of the PNES group. In the gamma band, the means of of every pair in the CPS group were larger than that of the PNES group. Almost every anatomically symmetric pair in the CPS patient group showed more powe r asymmetry than that in the PNES patient group. A Mann divergences (T indexes) of individually symmetric pairs were the same between the two groups and the p values are listed in Table 4 3 and Table 4 4 respectively. The
87 on the C3 C4 pair in the CPS group had significantly larger divergence than that in the PNES group. Analysis of variance (ANOVA) p artitions an observed variance into several components of some possible factors and provides a test of whether the means of groups are equal. I hypothesized that the variance of the T indexes, which quantify the inter hemispheric asymmetry, can be partitio ned into components explained by patient groups, pairs and frequency bands. Two way ANOVA considers two factors in a linear model to explain the interesting dependent variable and tests if the means of factor groups are the same. Amongst three factors, th e main interest is to test if the patient group is a strong factor explaining the variance of the T index. I first chose the patient groups and anatomically symmetric pairs as potential factors for explaining the variance of T indexes in the two way ANOVA for each frequency band. The results showed that patient groups were the significant factor in the T indexes and the T indexes were significantly different between the two patient groups under all frequency bands except the beta band. The ANOVA p values a re presented in Table 4 5. I later used the patient groups and frequency bands as factors and did the two way ANOVA again for each anatomically symmetric pair. The results are presented in Table 4 6 and show that the patient group was a significant factor explaining the variance of T index for those anatomically symmetric pairs in the frontal brain area. Discussion The global network structure and specific symmetric pair connections were investigated in both PNES and CPS EEG recordings. The small world net work index indicates how small world a graph is comparing to relatively random graphs having the same number of vertices and edges. In Table 4 1, the small world index was always
88 test showed significance in all frequency b ands for both patient groups. The non parametric Wilcoxon signed rank test showed the same results of being significant in all frequency bands. Both patient groups showed small world network structure in their functional graphs and the small world indexes had no significant difference between two groups. Only the minimum path length in the delta band showed significant difference between the two patient groups. Minimum path lengths in the delta band of the CPS group were mostly shorter than that of the PNES group ( p value= 0.037 ). This implied that the global integration is more efficient in CPS patients than that in PNES patients during the awake and relaxed state s This is not surprising because seizure itself is a synchronous activity amongst many brain ar eas. Furthermore about half of partial seizures happen during sleep. Moreover, temporal lobe complex partial seizures were more likely to secondarily generalize during sleep than during wakefulness (Herman et al. 2 001) My result could explain these interactions between partial seizures and sleep. Before entering the sleep state, the epileptic brain functional network showed higher global integration and could possibly facilitate the generation of seizures or the s econdary generalization during sleep. This study did not include any sleep EEG data so the hypothesis should be further pursued. Other network measures did not show that significant difference s could result from several reasons and the possible combination of these reasons. First, the pathological network structure may not manifest during the interictal awake and relaxed state of a patient. Second, the metric of association (PLAI) may not be sensitive enough to differentiate the pathological nuance of neuro nal interaction. Third, the network structure of the brain itself may be attack tolerant. Even though there have been some lesions
89 due to epileptic discharge, the brain could be organized in a fashion such that tissue damage does not affect much of the exi stence of the small world feature in the functional network. These hypotheses can be confirmed through further study For example, a study comparing the interictal and preictal network structure could determine whether or not the pathological network appea rs during the interictal state. Applying other interdependence measures could show a different perspective of the functional network in the awake and relaxed state. It would also be intriguing to design an integration index from different interdependence m easures so that the integration index values correspond to the synchronous activities of the greatest interest. I nter hemispheric power asymmetry is a more specific and local measure to investigate the quality of connection between hemispheres. For CPS pa tients, the brain tissue around the foci could be damaged by the recurrent onsets of partial seizures (So, 2000; Spooner et al. 2006; Worrell et al. 2000) The lesion of the brain tissue should affect the ensemble neuronal activity and cause the EEG signal to deviate from the EEG signal of the anatomically symmetric channel. On the other ha nd, PNES patients have more similarity between EEG signals from anatomically symmet ric channels due to symmetrically commensurate tissue integrity. Other than the lesion, the fixed focal hemisphere of partial seizures could result from substantially different pathological structure s within the hemisphere. If both of the abovementioned ca uses are valid and additive, the asymmetry could exacerbate
90 Table 4 1. Small world network index, of functional networks of CPS and PNES patients during awake and relaxed state. Delta Theta Alpha Beta CPS Mean 1.182 1.119 1.124 1.108 SD 0.052 0.063 0.036 0.063 T test p value 8.72E 05* 0.002510* 0.000103* 0.003772* PNES Mean 1.150 1.162 1.112 1.129 SD 0.054 0.049 0.038 0.037 T test p value 0.000315* 0.000124* 0.000232* 9.96E 05* *significant result when the hypothesis was that the is equal to one. Table 4 2. Mann Whitney U test results of network measures of CPS and PNES patient groups during awake and relaxed state. CPS PNES U test Mean SD Mean SD p value Delta(1 4Hz) 1.182 0.052 1.150 0.054 0.259 Theta(4 7Hz) 1.119 0.063 1.162 0.049 0.209 Alpha(8 12Hz) 1.124 0.036 1.112 0.038 0.456 Beta(13 30 Hz) 1.108 0.063 1.129 0.0 37 0.383 Lp/ < L p s > Delta(1 4Hz) 0.997 0.015 1.004 0.003 0.037* Theta(4 7Hz) 1.003 0.009 1.000 0.010 0.902 Alpha(8 12 Hz) 1.004 0.009 1.003 0.003 0.805 Beta(13 30 Hz) 1.001 0.008 1.003 0.006 0.383 Cp/ < C p s > Delta(1 4Hz) 1.177 0.041 1.155 0.056 0.383 Theta(4 7Hz) 1.121 0.064 1.162 0.043 0.383 Alpha(8 12 Hz) 1.128 0.035 1.115 0.040 0.535 Beta(13 30 Hz) 1.109 0.062 1.132 0.039 0.318 *significant result.
91 Table 4 3. Mann Whitney U test results of relative power asymmetry of CPS and PNES p atient groups during awake and relaxed state. F3 F4 C3 C4 P3 P4 O1 O2 F7 F8 T3 T4 T5 T6 Delta 0.209 0.017* 0.805 0.209 0.209 0.535 0.053 Theta 0.097 0.073 0.805 0.535 0.165 0.038* 0.805 Alpha 0.318 0.053 0.097 0.318 0.073 0.620 0.710 Beta 0.259 1.000 0.259 1.000 0.053 0.383 0.902 Gamma 0.535 0.710 1.000 0.710 0.128 0.535 1.000 *significant result. Table 4 sample T test results of relative power asymmetry of CPS and PNES patient groups during awake and relaxed state. F3 F4 C3 C4 P3 P4 O1 O2 F7 F8 T3 T4 T5 T6 Delta 0.184 0.029* 0.293 0.201 0.326 0.257 0.049* Theta 0.415 0.047* 0.633 0.412 0.236 0.054 0.944 Alpha 0.262 0.050* 0.057 0.364 0.075 0.623 0.532 Beta 0.244 0.424 0.450 0.644 0.110 0.170 0.752 Gamma 0.222 0.405 0.568 0.529 0.147 0.394 0.476 *significant result. Table 4 5. P values of two way ANOVA (patient groups as factor) in each frequency band. Frequency bands Delta Theta Alpha Beta Gamma P value 0.0036 0.0034 0.0019 0.0584 0.0254 *significant result. Table 4 6 P values of two way ANOVA (patient groups as factor) for each anatomically symmetric pair. Symmetric pair F3 F4 C3 C4 P3 P4 O1 O2 F7 F8 T3 T4 T5 T6 P value 0.011 0.001 0.109 0.070 0.002 0.017 0.905 *significant result.
92 Figure 4 1 E ye closed awa ke and relaxed state EEG signals of a patient having CPS. Figure 4 2. E ye closed awake and relaxed state EEG signals of a patient having PNES
93 Figure 4 3. A weighted adjacency matrix. The matrix has 18 rows and 19 columns. The adjacency matrix is suppo sed to be a square symmetric matrix. To eliminate the redundant information, the last row is eliminated and the diagonal as well as the lower triangular part of the original adjacency matrix are forced to be zero.
94 Figure 4 4. An adjacency matrix after a pplying a threshold on the weighted adjacency matrix in Figure 4 3. The matrix has 18 rows and 19 columns. The adjacency matrix is supposed to be a square symmetric matrix. To eliminate the redundant information, the last row is eliminated and the diagonal as well as the lower triangular part of the original adjacency matrix are forced to be zero. The threshold is applied and the entries in the matrix are either one or zero.
95 Figure 4 5. of individual anatomically symmetric pairs in the delta frequency band. Figure 4 6. of individual anatomically symmetric pairs in the theta frequency band.
96 Figure 4 7. of individual anatomically symmetric pairs in the alpha frequency band. Figure 4 8 of individual anatomically symmetric pairs in the beta frequency band.
97 Figure 4 9. of individual anatomically symmetric pairs in the gamma frequency band.
98 CHAPTER 5 CONNECTIVITY TRANSIT ION FROM INTERICATL TO ICTAL STATES ON NEOCORTICA L EPILEPSY Preamble Neocortical seizures originate in the neocortex the external surface part of the cerebral hemispheres. Neocortical epilepsy differs from mesial temporal epilepsy in that it is difficult to clearly define a single area from which the sei zures originate. Since seizures associated with neocortical epilepsy generally do not respond well to medication, epilepsy surgery is often one of the few options that patients with neocortical epilepsy have. However, surgery for neocortical epilepsy has a significantly lower success rate than other kinds of epilepsy. In patients with other types of epilepsy, such as mesial temporal lobe epilepsy (MTLE), surgeries can provide seizur e freedom in more than 70 to 90 % of cases (Benbadis and Tatum, 2000) Possible reasons that could result in the low success rate of neocortical epilepsy are: (1) there may be multiple seizure onset zones, (2) the initiation site may vary from seizure to seizure, and (3) with currently ava ilable technologies, it is difficult to precisely identify the duration and extent of seizure onset zones. A recent study by Lee et al in 2005 (Lee et al. 2005) on nonlesional neocortical epilepsy showed that, based on the epileptogenic focus locations, only 47 out of 89 patients were seizure free (Engel Class I) after the surgery, and an additional 7 experienced significant reduction in seizure frequency (Engel Class II). Therefore, developing a more reliable and eff ective method for identifying suitable parts and critical regions for neocortical epilepsy surgery would be a major contribution to improve the service of medicare and quality of life for those patients.
99 Many epilepsy researchers have postulated that there exists a network of anatomical regions in the cerebral cortex and its connections that initiates the transitions between normal (interictal) and seizure (ictal) states, resulting in recurrent seizures (Bertram et al. 1998; Monto et al. 2007; Spencer, 2002) In this case, identifying the configuration of the network and understanding the characteristics and dynamics of the network could lead to a greater understanding o f neocortical epilepsy. One of the most obvious hallmarks of the epileptogenic brain is the presence of interictal spikes in the EEG. However, areas that generate interictal spikes do not consistently correspond to the seizure onset zones (brain site where the seizure discharge first appears). Nevertheless, it is very likely that areas with interictal spikes and the seizure onset zones are both within the same epileptic brain network Furthermore, there may be other critical brain regions that are related to seizure occurrences but are not exhibited through interictal spikes or ictal discharges. My research, by analyzing intracranial EEG s seeks to explore quantitatively the existence and structure of a neocortical epilepsy brain network. The method s and re sult s section indicates how to identify the nodes and edges in the network quantitatively and how the characteristics of the brain network change from an interictal state to a seizure (ictal) state. Another important diagnostic tool for neocortical epilep sy surgery is high resolution magnetic resonance imaging (MRI). Often this offers a high predictive value for surgical outcomes (Bronen, 1992; Cascino et al. 1993) However, MRI does not wor k well in 29 % of patients with partial epilepsy (Semah et al. 1998) Many patients referred to epilepsy centers for surgery have normal MRI results (no leision) and previous studies
100 report that surgical outcomes were poor for patients with neocortical epilepsy with a normal MRI result (Cascino, 1990; Cascino et al. 1992; Cascino, 2004) Intracranial EEG recording is also indispe nsable for localizing seizure onset zones. However, sampling error can lead to false or missed localization during intracranial EEG recording. This problem is particularly common in cases with extratemporal seizure origin (Williamson et al. 1992a; Williamson et al. 1992b; Williamson, 1992) Nevertheless, several studies have demonstrated the usefulness of focus localization in neocortical epilepsy with quanti tative analysis on intracranial EEG. For example, Andrzejak et al in 1999 (Andrzejak et al. 1999) applied methods derived from nonlinear dynamics on intracranial interictal EEGs (N = 8) to lateralize the primar y epileptogenic area in two patients with neocortical epilepsy. They found that a nonlinear deterministic measure can contribute to an interictal localization of the primary epileptogenic area in patients with neocortical epilepsy. Based on a focus inde x measure, Roopun et al (Roopun et al. 2009) used a few examples to demonstrate that basic dynamic changes in focal epilepsy of neocortical origin may be useful in localizing the origin of seizures. The findings of these studies and many others support the hypothesis that the dynamics of intracranial EEGs are associated with the behavior of epileptogenic focus in neocortical epilepsy. Recent observations in humans with MTLE and in the animal models for this condit ion (Bertram, 1997; Spencer and Spencer, 1994) helped us concept ualize multifocal seizure onset (i.e. seizures that begin focally within different limbic structures with each seizure) as wel l as synchronized regional seizure onset (i.e. presumed simultaneous seizure initiation). These observations, together with multifocal physiological and anatomical changes in the animal models (Ben Ari et al. 1980;
101 Bertram et al. 1990; Bertram and Cornett, 1994; Bertram and Lothman, 1993; Cavalheiro et al. 1991; Du et al. 1993; Wuarin and Dudek, 1996) have raised the possibility of a widely distributed neural network (e.g., specific cortical and subcortical networks) in the gensis and expression of partial onset seizures Spencer (2002) (Spencer, 2002) defined a network to be a functionally and anatomically connected, bilaterally represented, set of cortical and subcortical brain structures and regions in servation on which this definition is based is that vulnerability to seizure activity in any one part of the network is influenced by activity everywhere else in the network, and that the network as a whole is responsible for the clinical and electrographi c phenomenon that we associate with human seizures. She further distinguishes between the network and the the ictal EEG. It is the area where the seizure discharge is first seen on ECoG (subdural electrodes). It may involve one or many electrode sites. In some patients, seizures (based on clinical manifestations) can start with EEG discharges beginning in different izure progresses. Spencer (2002) suggested that the existence of an epileptic network is supported also by the responses to invasive therapy (Spencer, 2002) If human epilepsy is the expression of specific, abnorma lly active, intrinsically defined and connected cortical/subcortical/bilateral networks, then one could theoretically alter seizure expression by intervening in any part of the specific network. Operations involving anterior temporal lobe, medial structure s only, lateral structures only, or more or less extensive lateral temporal resection can cure this disorder. Procedures with no
102 anatomic overlap are similarly successful (Spencer, 1996) This cannot be explained unless the multiple areas are all critical in the production of the intractable seizures of this disorder Then interruption of the network in any one of those areas would be (and apparently is) sufficient to alter the seizures. The similarly excellent re sponse with cessation of seizures after temporal lobe resection in well selected patients who have bilateral independent medial temporal lobe origin of seizures is another example of the existence of a network, interference with which at any site alters th e expression of the intractable seizures (Hirsch et al. 1991a; Hirsch et al. 1991b) Based on the above observations, it may be pertinent to consider study of other kinds of phenomena in ind ividual patients, which may define the network in better terms than we have sought in the past because of our single minded attention to defining regions of so called seizure onset For example, quantitative intracranial EEG analysis, background patterns, sleep effects on interictal and ictal activity, and other types of functional assessments may contribute considerably to our understanding of the role of networks in the expression of the epilepsies. Studying broad regions of brain structures related by th e presence of such networks, using quantitative EEG analyses and sophisticated approaches, may detect alterations in the behavior of the network before the more traditional manifestation clinically or on traditional EEG. Material Intracranial EEGs provides the most convincing evidence to support the network hypothesis (Spencer, 2002) Because the entire network participates in the expression of the seizure activity and can be entrained from any of its various parts, initial electrical
103 occurrence within the network. The initial area of apparent seizure involvement is no t really an onset area, because at any place in the network, and might even vary from seizure to seizure in a given patient. This locational variability performed in only one part of the network (King and Spencer, 1995; Spencer et al. 1992; Spencer and Spencer, 1994; Spencer, 1998) This may be the main reason for surgical failure in the neocortical epilepsy. A reliable quantitative method based on intracranial EEGs for determining the spatial distribution of the nodes of the epileptic network may lead to a better understanding of the mechanisms that lead to the generation of a seizure, and provide insights into more effective approaches to seizure control. We selected one subject with neocortical epilepsy and used the intracranial EEG of the subject to construct a representative network which is capa ble of showing the critical region with recording electrodes as vertices in the network. The patient is a 21 year old male, with a history of intractable seizures. Continuous VEEG monitoring was performed with XLTEK monitoring systems (Oakville, Ontario, C anada) equipment with a sampling rate of 499.707 Hz. Totally, 44 electrodes were applied on the left hemisphere cortex of subject. Data were obtained, stored, and interpreted according to ACNS guidelines. During the recording, frequent interictal spikes di scharges were present at the G4, G5, G9, G14, G15, IH2, and IH3 electrodes. Two EEG segments were cut and used to construct a brain network from a continuous EEG recording containing one neocortical epilepsy seizure from the patient. The seizure events we re identified by the recording facility based on both clinical
104 observation and a review of the EEG. In order to avoid the potential overlapping between ictal, postictal and preictal stages, only seizures that were at least 2 hours apart from the previous o r next one were used. The two segments analyzed were both 150 min long. One contain ed a seizure starting at 120 min and the other was seizure free. For the later, the previous or the next seizure happened at least 5 hours away from the beginning or ending of the seizure free segment. Both s egments were preproce ssed with a band pass filter (1 220 Hz) and several notch filters (60 Hz, 120 Hz, and 180 Hz) to remove the alias, power line artifacts and DC component. Methods This paper sought to identify the cri tical regions in the brain that are closely related to the development of neocortical epileptic seizure. In intracranial EEG recording, each electrode was treated (channel) as a node (vertex) in the network and cross correlation functions w ere used to calc ulate the weighting on the edge of each pair of nodes. To convert a weighted network into a binary one, I kept an edge stay in the network if and only if the two nodes adjacent to it were sufficiently connected. The brain network was constructed as follows Both 150 minute segments were further chopped into non overlapping calculation windows of 6 seconds to construct transient brain networks. In every calculation window, each electrode was treated as a node in the network. Pair wire cross correlation coeff icients were calculated among all 44 nodes. A cross correlation coefficient, with delay was calculated as Equation 5 1. (5 1)
105 (5 2) where and are EEG signal value s of a pair of electrodes at time and is the number of sampling points in each 6 second calculation window. The cross correlation coeffic ients with various delay (range from 1 to +1 second) were calculated. Then, the maximum of the absolute value of the cross correlation coefficient was obtained by selecting delay For each pair of nodes and the weight of edge is d efined by this maximum value. In this way, one obtains a weighted network. To convert it to a binary network, only the edges with a weight higher than a threshold (0.9) stay in the network. (5 3) (5 4) In each six second epoch, a graph representing the instant brain network was constructed and then the degree of electrode was calculated in Equation 5 5. (5 5) In Equation 5 5, is the total number of nodes implanted under the scalp and equals 44 for this neocortical epilepsy patient. The five minute mean occurrences of degree of node was calculated and then an average degree over all nodes, was calculated by averaging over all channels as Equation 5 6.
106 (5 6) Results and Hypothesis Testing After the coupling strength between each pair of nodes was quanti fied as the maximum cross correlation, the functional epileptic network was constructed There are many features of network structures and degree is one of the most widely used basic properties of network analysis. Based on the degree results, the connect ivity of an instant network was defined as the average degree of all nodes during a five minute period and quantified as One of the main purposes of this study was to find if there was a significant change in the instant network structures before and after a seizure onset Therefore, one of those proper hypotheses to test this idea was that the connectivity of a n e pileptic brain network increases when a brain approach es a seizure onset and then decreases after the onset. The results are shown as follows and the hypothesis test will be done in the discussion section. To have a better time resolution and globally spat ial layout of connectivity, was used to present the connectivity. After the calculation, the trajectories of were drawn The slope of the average degree over all channels was estimated between the beginning and 120 minute time point (i.e., onset time) of th e segment The slope from the 0 to 120 minute time point in the seizure segment was 0.02 / 5 minutes (p = 0.028), and was 0.0036 (p = 0.447) during the same time period of the interictal s egment. The degree s of nodes in the network were averaged over 300 calculation windows that span ned over 30 minutes and then plotted in Figure 5 4.
107 Discussion In Figure 5 2 and Figure 5 3, the regression slope of the segment with an onset is larger than zero while the segment without any onset is less than zero. This su ggested that the connectivity increases before an onset. Although the interictal segment had high values of around the 95 minute, the regression slope was not larger than zero. Based on this observation, the following can be concluded If the reco rding is done during interictal state, the connectivity should fluctuate within a range and have a regression slope around zero. If we have a preictal recording, we should more likely detect a positive regression slope within a period right before the seiz ure onset. The p value provides information about how unlikely the regression slope is compared to the value of zero. The regression slope equal to zero means that the connectivity does not correlate with the time within the segment. This is also what was presumed to happen during an interictal state. Additional support for this hypothesis is the result that the p value is equal to 0.0028 (<0.05) in the segment preceding an onset and 0.447 (>0.05) in the interictal segment. An epileptic seizure can be inte rpreted as an activity during which the neurons from different region overly correlate with each other. Therefore, one should observe a decrease of connectivity after a seizure. Network connectivity in terms of 30 minutes before and 30 minutes after the seizure onset were compared. The average degree before the seizure was 0.33 versus 0.09 after the seizure (p = 0.027). Applying the same comparison of average degrees in the interictal segment, the differ ence was not significant (p = 0.366). The results of this pilot study suggest that transitions in the brain networks may exist and are related to the underlying dynamics of seizures caused by neocortical epilepsy. An increase of connectivity can be observe d before seizure
108 onset. The onset of a seizure can cause the reset of the connectivity closer to a baseline level. Prokopyev et al (2003) conducted additional analys e s on the brain and found significant changes in the network characteristic during a seizu re (Boginski et al. 2005; Prokopyev et al. 2003) In this paper, a decrease of connection s was observed in results after the onset (in Figure 5 4 ) ; this is consistent with the reset pheno menon after seizure (Iasemidis et al. 2004; Sabesan et al. 2009) Monto et al in 2007 (Monto et al. 2007) studied epileptic brain networks by quanti fying the long range temporal correlation in subdural human EEG recordings. Their observations on the spread of abnormally large LRTCs suggested that the epileptic focus is associated with significant changes in network behavior even in the cortical areas immediately surrounding the clinically determined focus. By applying synchronization and graph analysis to intracranial EEG recordings, Ponten et al (2007) investigated the hypothesis that functional neuronal networks during temporal lobe seizures change in configuration before and during seizures (Ponten et al. 2007a) The study found that the functional brain networks change from a more random configuration during an interictal recording period to a more ordered configuration during seizures, especially during seizure spreading. The authors further suggested that the findings supported the theory that a random network (during interictal periods) even had a stronger tendency to synchronize (Chavez et al. 2006; Netoff et al. 2004) which could cause seizures. More recently, Zaveri et al (2009) (Zaveri et al. 2009) calculated a magnitude squared coherence on back ground intracranial EEG as a measure of functional connectivity to investigate the network effects within and outside the seizure onset area. The analysis demonstrated
109 an inverse relationship between the connectivity strength and the distance from the seiz ure onset area (especially during the frequency band). Another similar work by Ortega et al in 2008 (Ortega et al. 2008) showed the presence of clusters of increased synchronization in different locations on the l ateral temporal cortex in patients with temporal lobe epilepsy. As yet unproven it is possible that increased connectivity (as can be defined by coherence, synchronization clusters, or nonzero functional connectivity) helps initiate seizures. If increased brain connectivity or alteration in network topology helps initiate seizures, then network nodes and pathways could serve as targets for resective or disconnective surgery, implantable devices, and investigations of seizure anticipation. Applying classifi cation techniques for c ircumscribing some communities within the network could be done after the network is constructed and these communities may serve as good entities for observing the entrainment preceding a seizure. This is due, in part, b ecause they a re supposed to be less correlated during interictal periods and then entrain as evolving to ictal state. Several methods of identifying communities in network were introduced in Po r in 2009 (Porter et al. 2009)
110 Figure 5 1. Electrode placement of the neocortical epilepsy patient. There were totally 44 electrodes placed on the left frontal cortex of the subject.
111 Figure 5 2. of the segment with a seizure onset at the 120 minute. The blue line denotes the trajectory of the connectivity states for every five minutes. Least square linear regression was implemented and the result is indicated as the red line.
112 Figure 5 3 of the segment without any seizure activity. The blue line denotes the trajectory of the connectivity states for every five minutes. Least square linear regression was implemented and the result is i ndicated as the red line.
113 Fig ure 5 4. A ver aged degree over 300 calculation windows of nodes in the network. It includes five periods: 120~90 min before; 90~60 min before; 60~30 min before; 30~0 min (seizure onset) before; 0~30 min after seizure onset. Red lines separate different periods. In each period, all the 44 nodes (electrodes) in the network were plotted in bars. The height of each bar denotes the averaged degree of the corresponding node.
114 Figure 5 5. The average of over 30 minutes before and after a seizure onset. The error bars indicate the standard error of the mean. The p value is 0.027(<0.05) if we assume that the values of should be the same within 30 minutes before and after the seizure onset.
115 Figure 5 6. The average of over 30 minutes before and after an imaginary seizure onset that is located at 120 minute point. The error bars indicate the standard error of the mean. The p value is 0.366(>0.05) if we assume that the values of should be the same within 30 minutes
116 CHAPTER 6 CONCLUSION My work investigated several time series features applied to several applications relating to epileptic EEG signals. The first part of my research focused on a possible scheme for seizure p rediction. While there have been numerous studies that show promising results in anticipating seizure onsets using intracranial EEG (Chaovalitwongse et al. 2005; Le Van Quyen et al. 2005; Mormann et al. 2003b; Pacia and Ebersole, 2005; Rabbi et al. 2010; Schad et al. 2008) transfering these successes to seizure prediction based on scalp EEG is st ill a challenging task. Many studies so far have investigated the seizure prediction problem only under very well controlled conditions. Developing a robust scheme for seizure prediction outside the laboratory environment requires a variety of signal proce ssing techniques. The first part of my research explored the performance of a single feature based algorithm on long term continuous scalp EEG. Th e focus here was on evaluat ing an automated seizure prediction algorithm that issues seizure warnings by monit oring the convergence of signal regularity among EEG channels in continuous long term scalp EEG recordings. For each trial, t he algorithm was optimized with a training dataset and its performance was evaluated using a separate test dataset In the test da taset, the algorithm achieved an average sensitivity of 65% (i.e., ~ 2 out of every 3 seizures) with an average false positive rate of 0.24 per hour (~ 1 per 4 hours) The algorithm performance in the test dataset was nearly the same as that in the trainin g dataset This implies that the algorithm generated a stable performance i n epileptic patients Statistically, the overall sensitivity achieved by the test algorithm was better than a random warning scheme ( p value = 0.0145) While these results are encou raging, the performance of the algorithm
117 may not yet be sufficient for some clinical applications. For example, the performance would be of limited utility for inpatient monitoring applications, such as the EMU setting, due to the high false positive rate. However, the performance may be useful for driving seizure control devices, such as the vagus nerve stimulator. Furthermore, in this type of application, fine tuning the prediction algorithm to optimize performance for each individual patient could yield even better performance. Further research will be required to evaluate this type of closed loop seizure control application. Although the performance was not yet sufficient for its usage in a clinical environment, the results demonstrated the supremacy of the algorithm compared to a random scheme. This implied that the method found some events preceding the onsets of seizures. Further analysis needs to follow to be more specific about the positive findings in this study. It is also informative to investiga te the common characteristics of these negative results. For instance, the false positives may result mostly from a certain type of event and could be muted if the triggering event is found. This could bring down the false positive rate and possibly furthe r increase the sensitivity of the algorithm when the parameters are chosen to be more sensitive and the false positive rate is well suppressed to an acceptable level. Incorporating engineering techniques could help the prediction scheme become more applica ble in reality. For a complex system such as the brain, it may be too ambitious to hope that only one dynamic feature or method would be sufficient to accomplish effective seizure prediction. Aggregating the advantages from several time series dynamic feat ures can possibly join the benefit of each feature (Feldwisch Drentrup et al. 2010) Until now, most studies have compared one kind of method to another, such as linear vs. nonlinear (McSharry et al. 2003) and
118 univariate vs. bivariate features (Mirowski et al. 2009; Mormann et al. 2005; Mormann et al. 2007) These results demo nstrated that the bivariate features outperformed univarite features. Nonlinear features have been popular since 1999 (Andrzejak et al. 1999; Andrzejak et al. 2006; Hughes, 2008; Iasemidis et al. 2004; Lehnertz et al. 2001; Lehnertz, 2008; Maiwald et al. 2 004; Osterhage and Lehnertz, 2007; Pereda et al. 2005; Stam, 2005) but progress since then has not been substantial (Hughes, 2008) Until recently, linear methods such as the AR model or coherence have still been used and have generated satisfactory results in many EEG research fields including seizure prediction (Chisci et al. 2010) Different features explain the same signal from different perspectives and should be clearly correlated to the raw signal so that one can better utilize a feature at the proper time for a suitable purpose. Besides signal features, additional engineering methods could benefit the prediction scheme in the preparation and post feature stages. In the preparation stage, filtering, inverse source locating, and network analysis could help extract more seizure relevant information from the EEG signals before the dynamic features are applied to the signals. During the post feature stages, a decision to iss ue a warning or an index denoting the susceptibility of having a seizure should be output. A decider should be designed to consider situational information about the different features from the patient and integrate them in an optimal way. Optimization ski lls could be involved to assist in the decision process. Additionally, a prediction scheme with an updating threshold or a learning classifier would be more adaptive to the possibly changing background of EEG activities if one considers applying the predic tion scheme in an ambulant environment. There are still many aspects that can be improved and deserve greater analysis or hypothesis testing in
119 seizure prediction. The understanding of applicable features on the raw data and the integration of information are fundamental and necessary to accomplish this challenging task. In my initial study, the focus was on the temporal transition from the interictal to ictal states of a group of epileptic patients. Following this, I investigated the possible differences o f interictal network features between psychologically and physiologically epileptic patient groups. Graph theoretical analysis was applied to reveal the functional network structure of CPS and PNES patients when they were in the awake and relaxed intericta l state. The results showed the existence of a small world network configuration in both patient groups. Furthermore, the delta band minimum path length was significantly smaller in the CPS patient group. This implied that the low frequency global integrat ion was more efficient in CPS than in PNES patients. Both groups had the same small world network indexes and clustering coefficients. These results implied that the brain either somehow maintained the structure of a small world network or the original sma ll world network was attack tolerant so that the pathological influence of epileptic discharges did not bring down the small world configuration. The results above viewed the connection in a comprehensive scope over the whole brain. However, if the patholo gical effect was only locally misleading, the large scale analysis could fail to track down the difference. To take a closer look at the connection strength of several specific pairs, I retrieved the frequency power information and compared the inter hemis pheric asymmetry for both groups. The difference of power asymmetry was revealed in all frequency bands on most of the frontal pairs. The results showed that CPS patients had more power discrepancies between contralateral pairs than PNES
120 patients. These di fferences were not significant for those pairs around the rear part of the brain which may result from the epoch selection criterion. The awake and relaxed state EEG segments were selected based on the dominant alpha band oscillations appearing around the occipital channels. This selection criterion forced the presence of high alpha power in all segments and may have reduced the power of the statistical test. Some EEG studies in other fields also focus on power asymmetry only for the frontal region (Benca et al. 1999; Coan and Allen, 2003; Coan and Allen, 2004; Knott et al. 2001; Miskovic and Schmidt, 2010; Sutt on and Davidson, 2000) Future research for this study could involve applying a classifier to distinguish if a subject is a CPS or PNES patient by analyzing only a short period of the awake and relaxed state EEG. This would greatly decrease the recording time and resources of an EMU. Due to the requirement that a CPS patient must receive anti epileptic treatment, the classifier should be designed to preclude misdiagnosing a CPS patient as a PNES patient, but could be permitted to incorrectly identify a fe w PNES patients as CPS patients. The results of this study also suggest that future hypotheses may need to include more sensitive scales for detecting these important distinctions. For the big scale network analysis, these differences were not obvious or s ignificant. However, a more precise comparison of specific pairs disclosed more locally detailed information and showed the expected difference as hypothesized. In the third part of my dissertation, I examined the local functional network with a higher sp atial resolution around the onset region of a neocortical epilepsy patient. I used cross correlation to quantify the connection strength of every pair and then revealed the weighted graph of the functional network. The averaged degree of each vertex was th en
121 computed as the network feature unit. The overall averaged degree increased in the preictal segment while the interictal segment had a significantly smaller increase. The neocortical epilepsy patient that was the subject for this part of the research be came seizure free after a focus removal surgery. It would be intriguing to see if this case has a special network configuration compared to other neocortical epilepsy patients who were unable to benefit from focus removal surgery. For this particular study I had access to only one neocortical epilepsy case; but I look forward to continuing this research when more cases become available. The focus of my dissertation research was to discuss the potential usage of EEG to obtain clinical efficiency and provide relief for epileptic patients. It concludes by noting that the increased integration of engineering techniques and the utilization of detailed medical information are expected to improve the accuracy of the final decisions and provide real diagnostic and psychological relief to many current and future epileptic patients.
122 LIST OF REFERENCES Abramovich YI, Spencer NK, Turley MDE. Time varying autoregressive (TVAR) adaptive order and spectrum estimation. Signals, Systems and Computers, 2005. Con ference Record of the Thirty Ninth Asilomar Conference on 2005: 89 93. Achard S, Salvador R, Whitcher B, Suckling J, Bullmore E. A resilient, low frequency, small world human brain functional network with highly connected association cortical hubs. Journa l of Neuroscience 2006; 26: 63 72. Adak S. Time dependent spectral analysis of nonstationary time series. Journal of the American Statistical Association 1998; 93: 1488 9. Akaike H. Information theory and an extension of the maximum likelihood principle. Proceedings of the 2nd International Symposium on Information Theory, 1973: 267 281. Allan DW. Statistics of atomic frequency standards. Proceedings of the IEEE 1966; 54: 221 30. Andrzejak RG, Widman G, Lehnertz K, David P, Elger CE. Nonlinear determini sm in intracranial EEG recordings allows focus localization in neocortical lesional epilepsy. Epilepsia 1999; 40: 171 2. Andrzejak RG, Mormann F, Widman G, Kreuz T, Elger CE, Lehnertz K. Improved spatial characterization of the epileptic brain by focusing on nonlinearity. Epilepsy Research 2006; 69: 30 44. Andrzejak RG, Mormann F, Kreuz T, Rieke C, Kraskov A, Elger CE, et al. Testing the null hypothesis of the nonexistence of a preseizure state. Physical Review E 2003; 67: 010901(R). Aoki M, Havenner A. State space modeling of multiple time series. Econometric Reviews 1991; 10: 1 59. Appel U, Brandt AV. Adaptive sequential segmentation of piecewise stationary time series. Information Sciences 1983; 29: 27 56. Arnhold J, Grassberger P, Lehnertz K, Elger CE. A robust method for detecting interdependences: Application to intracranially recorded EEG. Physica D: Nonlinear Phenomena 1999; 134: 419 30. Bassett DS, Meyer Lindenberg A, Achard S, Duke T, Bullmore E. Adaptive reconfiguration of fractal small world human brain functional networks. Proceedings of the National Academy of Sciences 2006; 103: 19518 23.
123 Ben Ari Y, Tremblay E, Ottersen OP, Meldrum BS. The role of epileptic activity in hippocampal and "remote" cerebral lesions induced by kainic acid. Brai n Research 1980; 191: 79 97. Benbadis SR, Tatum WO. When drugs don't work: An algorithmic approach to medically intractable epilepsy. Neurology 2000; 55: 1780 4. Benca RM, Obermeyer WH, Larson CL, Yun B, Dolski I, Kleist KD, et al. EEG alpha power and al pha power asymmetry in sleep and wakefulness. Psychophysiology 1999; 36: 430 6. Bertram EH. Functional anatomy of spontaneous seizures in a rat model of limbic epilepsy. Epilepsia 1997; 38: 95 105. Bertram EH, Cornett JF. The evolution of a rat model of chronic spontaneous limbic seizures. Brain Research 1994; 661: 157 62. Bertram EH, Lothman EW, Lenn NJ. The hippocampus in experimental chronic epilepsy: A morphometric analysis. Annals of Neurology 1990; 27: 43 8. Bertram EH, Zhang DX, Mangan P, Fountai n N, Rempe D. Functional anatomy of limbic epilepsy: A proposal for central synchronization of a diffusely hyperexcitable network. Epilepsy Research 1998; 32: 194 205. Bertram EH,3rd, Lothman EW. Morphometric effects of intermittent kindled seizures and l imbic status epilepticus in the dentate gyrus of the rat. Brain Research 1993; 603: 25 31. Boersma M, Smit DJA, de Bie H, Van Baal GCM, Boomsma DI, de Geus EJC, et al. Network analysis of resting state EEG in the developing young brain: Structure comes wi th maturation. Human brain mapping Boginski V, Pardalos PM, Vazacopoulos A. Network based models and algorithms in data mining and knowledge discovery. Handbook of combinatorial optimization. Springer, Boston, MA, USA 2005: 217 58. Box GEP, Jenkins G. Time series analysis, forecasting and control. Prentice Hall PTR, Upper Saddle River, NJ, USA 1994. Bronen RA. Epilepsy: The role of MR imaging. American Journal of Roentgenology 1992; 159: 1165 74.
124 Bruzzo AA, Gesierich B, Santi M, Tassinari CA, Birbau mer N, Rubboli G. Permutation entropy to detect vigilance changes and preictal states from scalp EEG in epileptic patients. A preliminary study. Neurological sciences : official journal of the Italian Neurological Society and of the Italian Society of Clin ical Neurophysiology 2008; 29: 3 9. Buhl EH, Otis TS, Mody I. Zinc induced collapse of augmented inhibition by GABA in a temporal lobe epilepsy model. Science 1996; 271: 369 73. Cascino GD. Surgical treatment for epilepsy. Epilepsy Research 2004; 60: 179 86. Cascino GD. Epilepsy and brain tumors: Implications for treatment. Epilepsia 1990; 31: S37 44. Cascino GD, Jack Jr CR, Sharbrough FW, Kelly PJ, Marsh WR. MRI assessments of hippocampal pathology in extratemporal lesional epilepsy. Neurology 1993; 43 : 2380. Cascino GD, Jack CR, Parisi JE, Marsh WR, Kelly PJ, Sharbrough FW, et al. MRI in the presurgical evaluation of patients with frontal lobe epilespy and children with temporal lobe epilepsy: Pathologic correlation and prognostic importance. Epilepsy Research 1992; 11: 51 9. Cavalheiro EA, Leite JP, Bortolotto ZA, Turski WA, Ikonomidou C, Turski L. Long term effects of pilocarpine in rats: Structural damage of the brain triggers kindling and spontaneous recurrent seizures. Epilepsia 1991; 32: 778 82. Chaovalitwongse W, Iasemidis LD, Pardalos PM, Carney PR, Shiau DS, Sackellares JC. Performance of a seizure warning algorithm based on the dynamics of intracranial EEG. Epilepsy Research 2005; 64: 93 113. Chavez M, Hwang DU, Amann A, Boccaletti S. Synch ronizing weighted complex networks. Chaos: An Interdisciplinary Journal of Nonlinear Science 2006; 16: 015106. Chavez M, Le Van Quyen M, Navarro V, Baulac M, Martinerie J. Spatio temporal dynamics prior to neocortical seizures: Amplitude versus phase coup lings. Biomedical Engineering, IEEE Transactions on 2003; 50: 571 83. Chisci L, Mavino A, Perferi G, Sciandrone M, Anile C, Colicchio G, et al. Real time epileptic seizure prediction using AR models and support vector machines. Biomedical Engineering, IEE E Transactions on 2010; 57: 1124 32. Coan JA, Allen JJB. Frontal EEG asymmetry as a moderator and mediator of emotion. Biological Psychology 2004; 67: 7 50. Coan JA, Allen JJB. Frontal EEG asymmetry and the behavioral activation and inhibition systems. P sychophysiology 2003; 40: 106 14.
125 Corsini J, Shoker L, Sanei S, Alarcon G. Epileptic seizure predictability from scalp EEG incorporating constrained blind source separation. IEEE Transactions on Bio Medical Engineering 2006; 53: 790 9. Cossart R, Dinocou rt C, Hirsch JC, Merchan Perez A, De Felipe J, Ben Ari Y, et al. Dendritic but not somatic GABAergic inhibition is decreased in experimental epilepsy. Nature Neuroscience 2001; 4: 52 62. de Haan W, Pijnenburg YA, Strijers RL, van der Made Y, van der Flier WM, Scheltens P, et al. Functional neural network analysis in frontotemporal dementia and alzheimer's disease using EEG and graph theory. BMC Neuroscience 2009; 10: 1 12. Debener S, Beauducel A, Nessler D, Brocke B, Heilemann H, Kayser J. Is resting anterior EEG alpha asymmetry a trait marker for depression? Neuropsychobiology 2000; 41: 31 7. Dimitriadis SI, Laskaris NA, Del Rio Portilla Y, Koudounis GC. Characterizing dyn amic functional connectivity across sleep stages from EEG. Brain Topography 2009; 22: 119 33. Dimitriadis SI, Laskaris NA, Tsirka V, Vourkas M, Micheloyannis S. What does delta band tell us about cognitive processes: A mental calculation study. Neuroscien ce Letters 2010; 483: 11 5. Du F, Williamson J, Bertram E, Lothman E, Okuno E, Schwarcz R. Kynurenine pathway enzymes in a rat model of chronic epilepsy: Immunohistochemical study of activated glial cells. Neuroscience 1993; 55: 975 89. Engel Jr J, Shewm on DA. Overview: Who should be considered a surgical candidate. Surgical Treatment of the Epilepsies, 2nd edition 1993: 23 34. Federico P, Abbott DF, Briellmann RS, Harvey AS, Jackson GD. Functional MRI of the pre ictal state. Brain 2005; 128: 1811 7. Fe ldwisch Drentrup H, Schelter B, Jachan M, Nawrath J, Timmer J, Schulze Bonhage A. Joining the benefits: Combining epileptic seizure prediction methods. Epilepsia 2010; 51: 1598 606. Ferri R, Rundo F, Bruni O, Terzano MG, Stam CJ. Small world network organ ization of functional connectivity of EEG slow wave activity during sleep. Clinical Neurophysiology 2007; 118: 449 56. Fisher RS, Boas WE, Blume W, Elger C, Genton P, Lee P, et al. Epileptic seizures and epilepsy: Definitions proposed by the international league against epilepsy (ILAE) and the international bureau for epilepsy (IBE). Epilepsia 2005; 46: 470 2.
126 Fryzlewicz P, Van Bellegem S, von Sachs R. Forecasting non stationary time series by wavelet process modelling. Annals of the Institute of Statisti cal Mathematics 2003; 55: 737 64. Gal ZA, Boha R, Stam CJ, Molnr M. Age dependent features of EEG reactivity spectral, complexity, and network characteristics. Neuroscience Letters 2010; 479: 79 84. Gasser T, Bcher P, Mcks J. Transformations towards the normal distribution of broad band spectral parameters of the EEG. Electroencephalography and clinical Neurophysiology 1982; 53: 119 24. Gilliam F, Kanner AM. Treatment of depressive disorders in epilepsy patients. Epilepsy & Behavior 2002; 3: 2 9. Go mez V, Maravall A. Estimation, prediction, and interpolation for nonstationary series with the kalman filter. Journal of the American Statistical Association 1994; 89: 611 24. Gordon NJ, Salmond DJ, Smith AFM. Novel approach to nonlinear/non gaussian baye sian state estimation. Proceedings of the IEEE 1993; 140: 107 13. Gupta AK, Jeavons PM, Hughes RC, Covanis A. Aura in temporal lobe epilepsy: Clinical and electroencephalographic correlation. British Medical Journal 1983; 46: 1079 83. Hagemann D, Naumann E, Lrken A, Becker G, Maier S, Bartussek D. EEG asymmetry, dispositional mood and personality. Personality and Individual Differences 1999; 27: 541 68. Hannan EJ, Quinn BG. The determination of the order of an autoregression. Journal of the Royal Statis tical Society 1979; 41: 190 5. Harrison PJ, Stevens CF. Bayesian forecasting. Journal of the Royal Statistical Society.Series B (Methodological) 1976; 38: 205 47. Herman ST, Walczak TS, Bazil CW. Distribution of partial seizures during the sleep wake cyc le: Differences by seizure onset site. Neurology 2001; 56: 1453 9. Hirsch LJ, Spencer SS, Spencer DD, Williamson PD, Mattson RH. Temporal lobectomy in patients with bitemporal epilepsy defined by depth electroencephalography. Annals of Neurology 1991a; 30 : 347 56. Hirsch LJ, Spencer SS, Williamson PD, Spencer DD, Mattson RH. Comparison of bitemporal and unitemporal epilepsy defined by depth electroencephalography. Annals of Neurology 1991b; 30: 340 6.
127 Hively LM, Protopopescu VA. Channel consistent forewa rning of epileptic events from scalp EEG. IEEE Transactions on Biomedical Engineering 2003; 50: 584 93. Hively LM, Protopopescu VA, Gailey PC. Timely detection of dynamical change in scalp EEG signals. Chaos: An Interdisciplinary Journal of Nonlinear Scie nce 2000; 10: 864 75. based on information theoretic approaches in time series analysis. Physics Reports 2007; 441: 1 46. Huang Z, Chalabi ZS. Use of time series analysis to model and forecast wind speed. Journal of Wind Engineering & Industrial Aerodynamics 1995; 56: 311 22. Hughes JR. Progress in predicting seizure episodes with nonlinear methods. Epilepsy and Behavior 2008; 12: 128 35. Humphries MD, Gurney K, Prescott TJ The brainstem reticular formation is a small world, not scale free, network. Proceedings of the Royal Society B: Biological Sciences 2006; 273: 503 11. Iasemidis LD. Epileptic seizure prediction and control. IEEE Transactions on Biomedical Engineering 2 003; 50: 549 58. Iasemidis LD, Chris Sackellares J, Zaveri HP, Williams WJ. Phase space topography and the lyapunov exponent of electrocorticograms in partial seizures. Brain Topography 1990; 2: 187 201. Iasemidis LD, Shiau DS, Sackellares JC, Pardalos P M, Prasad A. Dynamical resetting of the human brain at epileptic seizures: Application of nonlinear dynamics and global optimization techniques. IEEE Transactions on Biomedical Engineering 2004; 51: 493 506. Iasemidis LD, Shiau DS, Pardalos PM, Chaovalitw ongse W, Narayanan K, Prasad A, et al. Long term prospective on line real time seizure prediction. Clinical Neurophysiology 2005; 116: 532 44. Indiradevi KP, Elias E, Sathidevi PS, Dinesh Nayak S, Radhakrishnan K. A multi level wavelet approach for automa tic detection of epileptic spikes in the electroencephalogram. Computers in Biology and Medicine 2008; 38: 805 16. Jacobs J, Zelmann R, Jirsch J, Chander R, Dubeau CCF, Gotman J. High frequency oscillations (80 500 hz) in the preictal period in patients with focal seizures. Epilepsia 2009; 50: 1780 92.
128 James CJ, Gupta D. Seizure prediction for epilepsy using a multi stage phase synchrony based system. Conference proceedings : Annual International Conference of the IEEE Engineering in Medicine and Biology Society.IEEE Engineering in Medicine and Biology Society.Conference 2009; 2009: 25 8. John ER, Ahn H, Prichep L, Trepetin M, Brown D, Kaye H. Developmental equations for the electroencephalogram. Science 1980; 210: 1255 8. Kalman RE. A new approach to l inear filtering and prediction problems. Journal of Basic Engineering 1960; 82: 35 45. Kanner AM. Depression in epilepsy: Prevalence, clinical semiology, pathogenic mechanisms, and treatment. Biological Psychiatry 2003; 54: 388 98. Kelly KM, Shiau DS, Ke rn RT, Chien JH, Yang MCK, Yandora KA, et al. Assessment of a scalp EEG based automated seizure detection system. Clinical Neurophysiology 2010; 121: 1832 43. King D, Spencer S. Invasive electroencephalography in mesial temporal lobe epilepsy. Journal of clinical neurophysiology : official publication of the American Electroencephalographic Society 1995; 12: 32 45. Kitagawa G, Kohn R, Ansley CF, Martin RD, Raftery AE, O'sullivan F, et al. Non gaussian state space modeling of nonstationary time series. com ment. rejoinder. Journal of the American Statistical Association 1987; 82: 1032 41. Knott V, Mahoney C, Kennedy S, Evans K. EEG power, frequency, asymmetry and coherence in male depression. Psychiatry Research: Neuroimaging 2001; 106: 123 40. Kreuz T, An drzejak RG, Mormann F, Kraskov A, Stgbauer H, Elger CE, et al. Measure profile surrogates: A method to validate the performance of epileptic seizure prediction algorithms. Physical Review E 2004; 69: 61915. Kwan P, Brodie MJ. Early identification of refr actory epilepsy. New England Journal of Medicine 2000; 342: 314 9. Lachaux JP, Rodriguez E, Martinerie J, Varela FJ. Measuring phase synchrony in brain signals. Human Brain Mapping 1999; 8: 194 208. Le Van Quyen M, Soss J, Navarro V, Robertson R, Chavez M, Baulac M, et al. Preictal state identification by synchronization changes in long term intracranial EEG recordings. Clinical Neurophysiology 2005; 116: 559 68.
129 Le Van Quyen M, Foucher J, Lachaux JP, Rodriguez E, Lutz A, Martinerie J, et al. Comparison of hilbert transform and wavelet methods for the analysis of neuronal synchrony. Journal of Neuroscience Methods 2001a; 111: 83 98. Le Van Quyen M, Martinerie J, Navarro V, Boon P, D'Hav M, Adam C, et al. Anticipation of epileptic seizures from standard EEG recordings. The Lancet 2001b; 357: 183 8. Lee SA, No YJ. Perceived self control of seizures in patients with uncontrolled partial epilepsy. Seizure: European Journal of Epilepsy 2005; 14: 100 5. Lee SK, Lee SY, Kim KK, Hong KS, Lee DS, Chung CK. Surg ical outcome and prognostic factors of cryptogenic neocortical epilepsy. Annals of Neurology 2005; 58: 525 32. Lehnertz K. Epilepsy and nonlinear dynamics. Journal of Biological Physics 2008; 34: 253 66. Lehnertz K, Litt B. The first international collab orative workshop on seizure prediction: Summary and data description. Clinical Neurophysiology 2005; 116: 493 505. Lehnertz K, Andrzejak RG, Arnhold J, Kreuz T, Mormann F, Rieke C, et al. Nonlinear EEG analysis in epilepsy:: Its possible use for intericta l focus localization, seizure anticipation, and prevention. Journal of Clinical Neurophysiology 2001; 18: 209. Lehnertz K, Mormann F, Osterhage H, Mller A, Prusseit J, Chernihovskyi A, et al. State of the art of seizure prediction. Journal of Clinical Ne urophysiology 2007; 24: 147 53. Litt B, Lehnertz K. Seizure prediction and the preseizure period. Current Opinion in Neurology 2002; 15: 173 7. Maiwald T, Winterhalder M, Aschenbrenner Scheibe R, Voss HU, Schulze Bonhage A, Timmer J. Comparison of three nonlinear seizure prediction methods by means of the seizure prediction characteristic. Physica D: Nonlinear Phenomena 2004; 194: 357 68. McGonigal D, Ionescu D. An outline for a kalman filter and r ecursive parameter estimationapproach applied to stock ma rket forecasting Canadian Conference on Electrical and Computer Engineering, 1995; 2: .1148 1151 McSharry PE, Smith LA, Tarassenko L. Comparison of predictability of epileptic seizures by a linear and a nonlinear method. Biomedical Engineering, IEEE Tran sactions on 2003; 50: 628 33.
130 Miltner WHR, Braun C, Arnold M, Witte H, Taub E. Coherence of gamma band EEG activity as a basis for associative learning. Nature 1999; 397: 434 5. Mirowski P, Madhavan D, LeCun Y, Kuzniecky R. Classification of patterns of EEG synchronization for seizure prediction. Clinical Neurophysiology 2009; 120: 1927 40. Miskovic V, Schmidt LA. Frontal brain electrical asymmetry and cardiac vagal tone predict biased attention to social threat. International Journal of Psychophysiology 2010; 75: 332 8. Monto S, Vanhatalo S, Holmes MD, Palva JM. Epileptogenic neocortical networks are revealed by abnormal temporal dynamics in seizure free subdural EEG. Cerebral Cortex 2007; 17: 1386 93. Mormann F, Elger CE, Lehnertz K. Seizure anticipat ion: From algorithms to clinical practice. Current Opinion in Neurology 2006; 19: 187 93. Mormann F, Lehnertz K, David P. Mean phase coherence as a measure for phase synchronization and its application to the EEG of epilepsy patients. Physica D: Nonlinear Phenomena 2000; 144: 358 69. Mormann F, Andrzejak RG, Elger CE, Lehnertz K. Seizure prediction: The long and winding road. Brain 2007; 130: 314 33. Mormann F, Kreuz T, Andrzejak RG, David P, Lehnertz K, Elger CE. Epileptic seizures are preceded by a dec rease in synchronization. Epilepsy Research 2003a; 53: 173 85. Mormann F, Andrzejak RG, Kreuz T, Rieke C, David P, Elger CE, et al. Automated detection of a preseizure state based on a decrease in synchronization in intracranial electroencephalogram recor dings from epilepsy patients. Physical Review E 2003b; 67: 21912. Mormann F, Kreuz T, Rieke C, Andrzejak RG, Kraskov A, David P, et al. On the predictability of epileptic seizures. Clinical Neurophysiology 2005; 116: 569 87. Nason GP. Wavelets in time se ries analysis. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 1999; 357: 2511 26. Netoff T, Park Y, Parhi K. Seizure prediction using cost sensitive support vector machine. Annual International Conferenc e of the IEEE Engineering in Medicine and Biology Society, 2009: 3322 2 5. Netoff TI, Clewley R, Arno S, Keck T, White JA. Epilepsy in small world networks. Journal of Neuroscience 2004; 24: 8075 83.
131 Ortega GJ, de la Prida LM, Sola RG, Pastor J. Synchroni zation clusters of interictal activity in the lateral temporal cortex of epileptic patients: Intraoperative electrocorticographic analysis. Epilepsia 2008; 49: 269 81. Osterhage H, Lehnertz K. Nonlinear time series analysis in epilepsy. International Jour nal of Bifurcation and Chaos 2007; 17: 3305 23. Osterhage H, Mormann F, Staniek MAU, Lehnertz K. Measuring synchronization in the epileptic brain: A comparison of different approaches. International Journal of Bifurcation and Chaos 2007; 17: 3539 44. Pac hou E, Vourkas M, Simos P, Smit D, Stam CJ, Tsirka V, et al. Working memory in schizophrenia: An EEG study using power spectrum and coherence analysis to estimate cortical activation and network behavior. Brain Topography 2008; 21: 128 37. Pacia SV, Ebers ole JS. Intracranial EEG substrates of scalp ictal patterns from temporal lobe foci. Epilepsia 2005; 38: 642 54. Palva S, Monto S, Palva JM. Graph properties of synchronized cortical networks during visual working memory maintenance. NeuroImage 2010; 49: 3257 68. Peled A. Multiple constraint organization in the brain: A theory for schizophrenia. Brain Research Bulletin 1999; 49: 245 50. Pereda E, Quiroga RQ, Bhattacharya J. Nonlinear multivariate analysis of neurophysiological signals. Progress in Neurob iology 2005; 77: 1 37. Peters TE, Bhavaraju NC, Frei MG, Osorio I. Network system for automated seizure detection and contingent delivery of therapy. Journal of Clinical Neurophysiology 2001; 18: 545 9. Ponten SC, Bartolomei F, Stam CJ. Small world netwo rks and epilepsy: Graph theoretical analysis of intracerebrally recorded mesial temporal lobe seizures. Clinical Neurophysiology 2007a; 118: 918 27. Ponten SC, Bartolomei F, Stam CJ. Small world networks and epilepsy: Graph theoretical analysis of intrace rebrally recorded mesial temporal lobe seizures. Clinical Neurophysiology 2007b; 118: 918 27. Porter MA, Onnela JP, Mucha PJ. Communities in networks. Notices of the American Mathematical Society 2009; 56: 1082 97. Prokopyev OA, Boginski V, Chaovalitwong se W, Pardalos PM, Sackellares JC, Carney PR. Network based techniques in EEG data analysis and epileptic brain modeling. Data Mining in Biomedicine 2003; 7: 559 73.
132 Protopopescu VA, Hively LM, Gailey PC. Epileptic event forewarning from scalp EEG. Journa l of Clinical Neurophysiology 2001; 18: 223 45. Qian S, Chen D. Joint time frequency analysis. Signal Processing Magazine, IEEE 1999; 16: 52 67. Quian Quiroga R, Kraskov A, Kreuz T, Grassberger P. Performance of different synchronization measures in real data: A case study on electroencephalographic signals. Physical Review E 2002; 65: 041903. Rabbi AF, Aarabi A, Fazel Rezai R. Fuzzy rule based seizure prediction based on correlation dimension changes in intracranial EEG. Annual International Conference of the IEEE Engineering in Medicine and Biology Society 2010; 1: 3301 4. Rice WR. A consensus combined P value test and the family wide significance of component tests. Biometrics 1990; 46: 303 8. Roopun AK, Traub RD, Baldeweg T, Cunningham MO, Whittaker RG, Trevelyan A, et al. Detecting seizure origin using basic, multiscale population dynamic measures: Preliminary findings. Epilepsy and Behavior 2009; 14: 39 46. Rothkegel A, Lehnertz K. Multistability, local pattern formation, and global collective fir ing in a small world network of nonleaky integrate and fire neurons. Chaos: An Interdisciplinary Journal of Nonlinear Science 2009; 19: 015109. Rubinov M, Knock SA, Stam CJ, Micheloyannis S, Harris AWF, Williams LM, et al. Small world properties of nonlin ear brain activity in schizophrenia. Human Brain Mapping 2009; 30: 403 16. Sabesan S, Chakravarthy N, Tsakalis K, Pardalos P, Iasemidis L. Measuring resetting of brain dynamics at epileptic seizures: Application of global optimization and spatial synchron ization techniques. Journal of Combinatorial Optimization 2009; 17: 74 97. Sackellares JC, Shiau DS, Principe JC, Yang MCK, Dance LK, Suharitdamrong W, et al. Predictability analysis for an automated seizure prediction algorithm. Journal of Clinical Neuro physiology 2006; 23: 509 20. Salvador R, Suckling J, Coleman MR, Pickard JD, Menon D, Bullmore E. Neurophysiological architecture of functional magnetic resonance images of human brain. Cerebral Cortex 2005; 15: 1332 42. Sander JW. The epidemiology of ep ilepsy revisited. Current Opinion in Neurology 2003; 16: 165 70.
133 Schad A, Schindler K, Schelter B, Maiwald T, Brandt A, Timmer J, et al. Application of a multivariate seizure detection and prediction method to non invasive and intracranial long term EEG recordings. Clinical neurophysiology : official journal of the International Federation of Clinical Neurophysiology 2008; 119: 197 211. Schelter B, Winterhalder M, Maiwald T, Brandt A, Schad A, Schulze Bonhage A, et al. Testing statistical significance of multivariate time series analysis techniques for epileptic seizure prediction. Chaos: An Interdisciplinary Journal of Nonlinear Science 2006; 16: 013108. Semah F, Picot MC, Adam C, Broglin D, Arzimanoglou A, Bazin B, et al. Is the underlying cause of epi lepsy a major prognostic factor for recurrence? Neurology 1998; 51: 1256 62. Shiau D. Signal identification and forecasting in nonstationary time series data. 2001. Shorvon SD. The epidemiology and treatment of chronic and refractory epilepsy. Epilepsia 1996; 37: S1 3. So EL. Integration of EEG, MRI, and SPECT in localizing the seizure focus for epilepsy surgery. Epilepsia 2000; 41: S48 54. Spencer SS. Neural networks in human epilepsy: Evidence of and implications for treatment. Epilepsia 2002; 43: 219 27. Spencer SS. Substrates of localization related epilepsies: Biologic implications of localizing findings in humans. Epilepsia(Copenhagen) 1998; 39: 114 23. Spencer SS. Long term outcome after epilepsy surgery. Epilepsia 1996; 37: 807 13. Spencer SS, Spencer DD. Entorhinal hippocampal interactions in medial temporal lobe epilepsy. Epilepsia 1994; 35: 721 7. Spencer SS, Guimaraes P, Katz A, Kim J, Spencer D. Morphological patterns of seizures recorded intracranially. Epilepsia 1992; 33: 537 45. Spoon er CG, Berkovic SF, Mitchell LA, Wrennall JA, Harvey AS. New onset temporal lobe epilepsy in children: Lesion on MRI predicts poor seizure outcome. Neurology 2006; 67: 2147 53. Stam CJ. Nonlinear dynamical analysis of EEG and MEG: Review of an emerging fi eld. Clinical Neurophysiology 2005; 116: 2266 301. Stam CJ. Functional connectivity patterns of human magnetoencephalographic recordings: A  small world'network? Neuroscience Letters 2004; 355: 25 8.
134 Stam CJ, Van Dijk BW. Synchronization likelihood: An unbiased measure of generalized synchronization in multivariate data sets. Physica D: Nonlinear Phenomena 2002; 163: 236 51. Stam CJ, Nolte G, Daffertshofer A. Phase lag index: Assessment of functional connectivity from multi channel EEG and MEG with dim inished bias from common sources. Human Brain Mapping 2007; 28: 1178 93. Sutton SK, Davidson RJ. Prefrontal brain electrical asymmetry predicts the evaluation of affective stimuli. Neuropsychologia 2000; 38: 1723 33. Swartz BE, Goldensohn ES. Timeline of the history of EEG and associated fields. Electroencephalography and Clinical Neurophysiology 1998; 106: 173 6. Takens F. Detecting strange attractors in turbulence. Dynamical systems and turbulence 1981; 898: 366 81 in Rand, David and Lai Sang Young, ed s., Dynamical Systems and Turbulence, Warwick 1980 Taswell C, Toolsmiths C, Stanford CA. The what, how, and why of wavelet shrinkage denoising. Computing in Science & Engineering [see also IEEE Computational Science and Engineering] 2000; 2: 12 9. Thatc her RW, North DM, Curtin RT, Walker RA, Biver CJ, Gomez JF, et al. An EEG severity index of traumatic brain injury. Journal of Neuropsychiatry and Clinical Neurosciences 2001; 13: 77 87. Tsakalis K, Chakravarthy N, Sabesan S, Iasemidis LD, Pardalos PM. A feedback control systems view of epileptic seizures. Cybernetics and Systems Analysis 2006; 42: 483 95. Van Dellen E, Douw L, Baayen JC, Heimans JJ, Ponten SC, Vandertop WP, et al. Long term effects of temporal lobe epilepsy on local neural networks: A gr aph theoretical analysis of corticography recordings. PLoS One 2009; 4: e8081. Vuong QH. Likelihood ratio tests for model selection and non nested hypotheses. Econometrica 1989; 57: 307 33. orks. Nature 1998; 393: 440 2. West M, Prado R, Krystal AD. Evaluation and comparison of EEG traces: Latent structure in nonstationary time series. Journal of the American Statistical Association 1999; 94: 375 6.
135 Williamson PD, Thadani VM, Darcey TM, Spe ncer DD, Spencer SS, Mattson RH. Occipital lobe epilepsy: Clinical characteristics, seizure spread patterns, and results of surgery. Annals of Neurology 1992a; 31: 3 13. Williamson PD, Boon PA, Thadani VM, Darcey TM, Spencer DD, Spencer SS, et al. Parieta l lobe epilepsy: Diagnostic considerations and results of surgery. Annals of Neurology 1992b; 31: 193 201. Williamson PD. Frontal lobe seizures. problems of diagnosis and classification. Advances in Neurology 1992; 57: 289 309. Winterhalder M, Maiwald T, Voss HU, Aschenbrenner Scheibe R, Timmer J, Schulze Bonhage A. The seizure prediction characteristic: A general framework to assess and compare seizure prediction methods. Epilepsy and Behavior 2003; 4: 318 25. Worrell GA, Lagerlund TD, Sharbrough FW, Br inkmann BH, Busacker NE, Cicora KM, et al. Localization of the epileptic focus by low resolution electromagnetic tomography in patients with a lesion demonstrated by MRI. Brain Topography 2000; 12: 273 82. Wuarin JP, Dudek FE. Electrographic seizures and new recurrent excitatory circuits in the dentate gyrus of hippocampal slices from kainate treated epileptic rats. Journal of Neuroscience 1996; 16: 4438. Zandi AS, Dumont GA, Javidan M, Tafreshi R. An entropy based approach to predict seizures in temporal lobe epilepsy using scalp EEG. Conference proceedings : Annual International Conference of the IEEE Engineering in Medicine and Biology Society.IEEE Engineering in Medicine and Biology Society. Conference 2009; 2009: 228 31. Zaveri HP, Pincus SM, Gonchar ova II, Duckrow RB, Spencer DD, Spencer SS. Localization related epilepsy exhibits significant connectivity away from the seizure onset area. Neuroreport 2009; 20: 891 5. Zaveri HP, Williams WJ, Sackellares JC, Beydoun A, Duckrow RB, Spencer SS. Measuring the coherence of intracranial electroencephalograms. Clinical Neurophysiology 1999; 110: 1717 25.
136 BIOGRAPHICAL SKETCH Jui Hong Chien was born in 1981 in Kaohsiung, Taiwan, R.O.C. He was the first grandson in his furniture manufacturing family. He gre w up with his family in Kaohsiung and graduated from Kaohsiung Senior High School in 1999. He later received his Bachelor of Science degree in civil engineering from National Taiwan University in 2003. He then served as a Second Lieutenant of Engineering A rmy from 2003 to 2005. Due to the death of his father when he was 18 years old, Jui Hong decided to pursue his graduate degree in biomedical engineering. After an honorable discharge ing at the University of Florida in 2006. After earning his Master of Science degree in biomedical engineering in 2008, he continued the Doctor of Philosophy program in biomedical engineering at the University of Florida under the mentorship of Dr. Panos P ardalos in the Center for Applied Optimization Under the supervision of Dr. Pardalos, Jui Hong was introduced to an epilepsy research company, Optima Neuroscience Inc., and started to work as an intern data analyst from 2008 till his completion of graduate study. During that period of time, he worked closely with his adviser, Dr. Pardalos, and his supervisor, Dr. Deng Shen Shiau, in Optima Neuroscience Inc. to apply engi neering techniques to the analysis of epileptic electroencephalographic signal. A fter h e received his doctoral degree from the University o f Flori da in the spring of 2011 h is research interests r emain in a nalysis of bi o logical si g n als