Identification of epilepsy stages from ECoG using genetic programming classifiers

Computers in Biology and Medicine 43 (2013) 1713–1723

Contents lists available at ScienceDirect

Computers in Biology and Medicine journal homepage: www.elsevier.com/locate/cbm

Identification of epilepsy stages from ECoG using genetic programming classifiers Arturo Sotelo a, Enrique Guijarro b, Leonardo Trujillo a,n, Luis N. Coria a, Yuliana Martínez a a Doctorado en Ciencias de la Ingeniería, Departamento de Ingeniería Eléctrica y Electrónica, Instituto Tecnológico de Tijuana, Blvd. Industrial y Av. ITR Tijuana S/N, Mesa Otay C.P. 22500, Tijuana BC, Mexico b Departamento de Ingeniería Electrónica, Universidad Politécnica de Valencia, Spain

art ic l e i nf o

a b s t r a c t

Article history: Received 18 December 2012 Accepted 21 August 2013

Objective: Epilepsy is a common neurological disorder, for which a great deal of research has been devoted to analyze and characterize brain activity during seizures. While this can be done by a human expert, automatic methods still lag behind. This paper analyzes neural activity captured with Electrocorticogram (ECoG), recorded through intracranial implants from Kindling model test subjects. The goal is to automatically identify the main seizure stages: Pre-Ictal, Ictal and Post-Ictal. While visually differentiating each stage can be done by an expert if the complete time-series is available, the goal here is to automatically identify the corresponding stage of short signal segments. Methods and materials: The proposal is to pose the above task as a supervised classification problem and derive a mapping function that classifies each signal segment. Given the complexity of the signal patterns, it is difficult to a priori choose any particular classifier. Therefore, Genetic Programming (GP), a population based meta-heuristic for automatic program induction, is used to automatically search for the mapping functions. Two GP-based classifiers are used and extensively evaluated. The signals from epileptic seizures are obtained using the Kindling model of elicited epilepsy in rodent test subjects, for which a seizure was elicited and recorded on four separate days. Results: Results show that signal segments from a single seizure can be used to derive accurate classifiers that generalize when tested on different signals from the same subject; i.e., GP can automatically produce accurate mapping functions for intra-subject classification. A large number of experiments are performed with the GP classifiers achieving good performance based on standard performance metrics. Moreover, a proof-of-concept real-world prototype is presented, where a GP classifier is transferred and hard-coded on an embedded system using a digital-to-analogue converter and a field programmable gate array, achieving a low average classification error of 14.55%, sensitivity values between 0.65 and 0.97, and specificity values between 0.86 and 0.94. Conclusions: The proposed approach achieves good results for stage identification, particularly when compared with previous works that focus on this task. The results show that the problem of intra-class classification can be solved with a low error, and high sensitivity and specificity. Moreover, the limitations of the approach are identified and good operating configurations can be proposed based on the results. & 2013 Elsevier Ltd. All rights reserved.

Keywords: Epilepsy diagnosis Genetic programming Classification

1. Introduction Epilepsy is a common neurological ailment, characterized by the chronic seizures it causes as part of its symptomatology. According to different studies, epilepsy incidence varies depending on the region and the considered population. For instance, [1] estimates that the number of people with epilepsy is between n

Corresponding author. Tel.: þ 52 6646072942. E-mail addresses: [email protected] (A. Sotelo), [email protected] (E. Guijarro), [email protected] (L. Trujillo), luis.coria@tectijuana. edu.mx (L.N. Coria), [email protected] (Y. Martínez). 0010-4825/$ - see front matter & 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.compbiomed.2013.08.016

11/100,000 and 134/100,000. Other studies concluded that 3% to 5% of the general population experiences one or more seizures during their life-time [2]. A more conservative estimate is that epilepsy affects 1% of the world population [3,4]. According to [5], two thirds of affected individuals develop seizures that can be controlled by anti-epileptic medication, while another 7% or 8% can be cured by surgery. However, the symptomatology of the remaining 25% cannot be controlled by current therapies. Epileptic seizures are sudden disruptive episode of mental functions, that develop over four principal stages [6]: (1) The Basal stage, (2) Pre-Ictal Stage, (3) The Ictal Stage; and (4) The Post-Ictal Stage. At each stage we can appreciate different frequencies and

1714

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

waveforms; i.e., each stage is characterized by a different signal morphology. However, while a human expert has no problem in identifying each stage, an automatic method for stage identification has not been developed. Concretely, the proposed system is designed to automatically determine the stage to which a small signal segment belongs. The task is posed as a supervised learning problem, where the system input is a short signal segment and the output is the corresponding stage of the seizure. However, deriving automatic processing methods for these signals is definitely not a straightforward endeavor, given the complexities of the signals generated during a seizure, which often contain several frequency components [7,8]. Therefore, in this paper the problem is solved using a Genetic Programming (GP) classifier, that analyzes statistical features of each signal segment and derives a non-linear mapping following a symbolic regression strategy [9,10]. GP is a stochastic population based meta-heuristic which is well suited to automatically derive mapping functions. In particular, GP is used when the overall structure of a solution cannot be defined a priori and only a highlevel description of the desired functionality is available. This paper presents a continuation of preliminary research carried out by the authors addressing stage identification [11]. However, the present work presents vastly more experimental tests and validation, two orders of magnitude more in total experimental runs to obtain statistically sound results while considering many different parametrizations of the system. Moreover, experimental tests include more test subjects, two variants of GP classifiers and a prototype for an embedded system that could perform stage identification in a real-world scenario. The remainder of this paper proceeds as follows. First, Section 2 presents a brief introduction to epilepsy models and brain signals recorded through ECoG. Then, Section 3 reviews related work. Afterwards, Section 5 presents the supervised learning problem and the proposed GP-based solution. Section 6 presents the experimental setup and provides a detailed discussion of the results. Finally, concluding comments are given in Section 7.

2. Epilepsy model, seizure stages and signal recording Epileptic seizures are mostly spontaneous, a characteristic that makes them particularly difficult to study. Therefore, for research purposes seizures are often elicited in a controlled manner using animal subjects, commonly rodents. Indeed, one of the most used models is the amygdala Kindling model for temporal lobe epilepsy [12,13], since signal morphology is very similar to those produced by a human brain during a seizure. Using animal models, it is possible to reproduce a chronic brain dysfunction that leads to epilepsy, a strategy that has allowed for research regarding the causes and mechanisms behind epilepsy [14]. The Kindling model is used to study acquired epilepsy that is elicited by electrical impulses delivered to a previously healthy (nonepileptic) animal. Epileptic conditions are achieved in the animal as a result of short duration electrical stimulus in the limbic regions of the brain, such as the amygdala or hippocampus. Following the Kindling model, spontaneous seizures are elicited by an electrical stimulus discharged directly to the brain of an animal, in this case a rodent. The approach has several advantages, such as a precise focal activation and the development of chronic epileptogenesis [15]. Kindling seizures are rated, depending on their symptoms, into a five level scale called the Racine scale [16]. This scale rates seizure intensity from focal to generalized, depending on the symptomatology exhibited by rats from the Wistar breed, where a level-5 is the highest intensity. Symptom for each level are: (0) no seizure response; (1) immobility, eye closure, twitching of vibrissae; (2) head nodding associated with more severe facial clonus; (3) clonus of one

forelimb; (4) bilateral forelimb clonus with rearing; and (5) rearing and falling on the back accompanied by generalize clonic seizures. In the present study, level-5 seizures (generalized motor seizures) are used for the experimental analysis. Despite ongoing discussions regarding the generalization of experimental work based on the Kindling model [17], it is still widely used for the study of new treatments [12,18]. 2.1. Brain signals during epileptic seizures Brain activity produces a highly non-periodical signal with amplitudes in the range of 0:5 μV–100 μV. Such signals can be detected by non-invasive methods when they are recorded at scalp level. In fact, Electroencephalography (EEG) is the most common tool for the diagnosis and treatment of many neurological disorders, including epilepsy. However, EEG signals are normally contaminated by noise or artifacts, mostly produced by muscle activity such as the blinking of an eye or even normal heart activity. It is also possible to obtain less distorted signals through intracranial recording methods, which are far more amendable to precise analysis and clinical evaluations [19], what it is called Electrocorticogram (ECoG). Intracranial detection is accomplished using an array of metal disks that are placed on the surface of the brain, what are called subdural electrodes, or by inserting metal needles within the brain at the required depths, these are called deep electrodes. Moreover, such electrodes are bidirectional, which allows them to be used for both detection and direct stimulation [20]. Of course, the main drawback of such methods is that they require surgical access to a patient's brain. This paper presents an approach to automatically discriminate the three main stages (Pre-Ictal, Ictal, and Post-Ictal) within a brain signal during an epileptic seizure, based on the local dynamics and morphology of the recorded ECoG from elicited epileptic episodes. 2.2. Seizure stages Epileptic seizures are sudden disruptive episodes of mental, motor, sensorial and autonomic functions, caused by a paroxysmal malfunction of brain cells, which is considered an abnormal increase of neural synchrony [21]. Epilepsy may affect the brain of a patient partially or completely, respectively producing partial or generalized seizures [15]. An epileptic seizure normally develops over several stages [6], these are: (1) The Basal stage, (2) PreIctal Stage, (3) The Ictal Stage; and (4) The Post-Ictal Stage. These stages are defined by signal morphology and time domain features, because symptomatology in the patient is only evident for the Ictal stage. The Basal stage corresponds to normal brain functions, in this stage brain signals are characterized by a low amplitude and a relative high frequency. In general these brain signals exhibit normal activity induced by muscle movements or the state of vigilance [22] of the subject; the Basal-stage is not part of an epileptic seizure, it is therefore not considered in the experimental tests for stage identification in this work. In the Pre-Ictal stage, the corresponding EEG or ECoG signals show a considerable amplitude increase relative to the Basal stage, since hyper-excitable cells start to discharge epileptogenic focus and the signal begins to show spikes and transitory activity but no definitive evolution. Nonetheless, during this stage the patient may exhibit clinical manifestations similar to a focal seizure, depending upon the origin of the epileptic discharges. The Ictal stage is precisely when the seizure happens; depending on the affected brain zone it could produce jerky movements, olfactory sensations and even a loss of consciousness if it is a complex partial seizure. If these discharges become sufficiently widespread a convulsive response becomes possible. Brain signals during this stage are identifiable by a high amplitude, a low frequency and a

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

1715

Fig. 1. ECoG signal of a level-5 seizure on Racine scale, (a) Seizure recorded by the cortex electrode, (b) Deep recording through the stimulus electrode.

Fig. 2. Illustration of approximate stereotaxic locations of stimulus, recording, and reference electrodes in a adult male Wistar rat [42].

predominant rhythm. The last stage is the Post-Ictal stage, where signals show general amplitude depression as the patient gradually returns to the Basal stage and the symptomatology begins to cease. Fig. 1(a) depicts the electrical activity of an ECoG signal taken over an entire epileptic episode, where the three main stages of the seizure are clearly marked. If an expert neurologist analyzes the EEG or ECoG signal of a patient that suffered a seizure he can identify the seizure stages as they appear in the time domain. 2.3. Signal recording The electrode implants, the Kindling experiments using live rodents (Wistar rats), and the signal recordings were all carried out at the Centro de Investigación, Hospital General Universitario de Valencia, in Valencia Spain; in accordance with institutional guidelines for human animal care. Stimulation and signal recording were achieved by inserting electrodes within the skull of each rat through symmetric burr holes at stereotaxic locations following [23,18,24]. Fig. 2 shows the approximated stereotaxic location of the electrodes, where the black marks represent stimulus electrodes, orange are for reference, and blue and red represent the cortical frontal and occipital recording electrodes. To induce the seizure, an electrical discharged is applied to the amygdala of each subject as stimulation, through an electrode made of twisted pair of Teflon-coated 0.25 mm diameter stainless steel wires separated by 0.5 mm at the tip and 8 mm in length and implanted through a burr hole. Two stainless steel screws served as cortical recording electrodes, attached to a connector assembly.

After this process was done, the stimulation of rodent subjects began after 7 days. The electric stimulation of the subject's brain tissue and deep recording of the signal is done after the electrodes are implanted and the connector is plugged in. Stimulation and recording of brain activity begins as soon as the rat is connected so the rat does not remove the cable. The applied stimulus consist of a 500 μA at 50 Hz rectangular signal with a 5% duty cycle by 1 s. Electrical manifestations are of variable amplitude, with useful frequency components from 0.5 Hz to 60 Hz [25], and may find useful components up to 100 Hz [26] or 400 Hz [27]. For this work, the signal was bandpass filtered at a 0.5 Hz to 100 Hz bandwidth, sampling rate was 256 Hz to avoid aliasing, using a 12 bit resolution. Fig. 1(a) represents the complete time-series record for a seizure, from the Pre-Ictal stage that begins at second 480 and ends at second 535. Then, the Ictal stage continues up to second 575, and finally the Post-Ictal stage represents the final part of the signal. Most seizures last between 90 and 180 s, with an average duration of roughly 104.6 s for the entire episode. However, this includes some discarded noise at the beginning and end of each recording, as well as some brain activity related to the Basal stage. The average duration of each stage is 30.4, 19.4 and 25.8 s, respectively for the Pre-Ictal, Ictal and Post-Ictal stage. The plot of Fig. 1(b) is the signal from the deep recording electrode, which provides a reference to determine when the seizure is about to start. When the seizure is detected at the cortex level, this means that the stimulus has produced a discharge that stimulated nearby neurons up to the cortex, producing a generalized seizure in the rodent. On the other hand, in some experiments the deep recording electrode shows epileptic activity, but the cortex electrode does not, which represents a local seizure. However, only seizures with a level-5 rating on the Racine scale are considered in this study. In some cases, more than a single stimulus needs to be applied to a subject to induce a level-5 seizure. Fig. 1(a) shows an ECoG for a level-5 seizure.

3. Automatic analysis of epilepsy signals Hitherto, most efforts towards the automatic analysis of epileptic brain signals have centered on the problem of predicting the onset of a seizure or identifying if a signal exhibits traits of an epileptic episode. For instance, [28] analyzes the non-linear

1716

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

dynamics of the signal in an attempt to anticipate a seizure. Similarly, [29,30] use hybrid features and apply computational techniques to recognize the signal dynamics exhibited during the Pre-Ictal stage. Moreover, given the difficulty of the problem, several researchers have turned to heuristic approaches to produce results that are competitive with human analysis. In particular, as done in this work, several authors have proposed GP approaches, where basic signal features are used as decision variables and the task is posed as a classification problem. For instance, [31] uses GP to detect two common patterns found in many brain signals from epileptic seizures (spike-and-slow-wave-complex and spike-or-sharp-wavecomplex), and in this way determine if an EEG signal exhibits seizure features. Classification is performed by GP considering seven basic features, using a relatively small training and test set. In [32], GP is used to detect signal patterns that are precursors to an epileptic episode. The authors analyze signals that are recorded intracranially and consider a set of common signal features from signal processing literature. A final example is the work presented in [33], where GP is used to derive signal features that are then given as input to a conventional classifier, what is known as a wrapper method, to detect the presence of an epileptic seizure in an EEG recording. In fact, progress on the problem of seizure prediction has allowed researchers to move from the experimental domain to clinical applications. Particularly, it has lead to the development of wearable electronic devices that can detect the onset of a seizure [34,35]. This work differs in goal and methodology from the previously reviewed literature. The problem is discriminating between each of the three main stages of an epileptic seizure; i.e., stage identification. It differs from the detection of seizure precursors or the binary-classification of determining the presence of an epileptic episode. In effect, this problem takes a fine grained look at seizure analysis, focusing on discriminating between signal patterns that share common features. Such a system can allow for the automatic study of seizure evolution, with possible applications in the analysis and validation of epileptic treatments over the complete duration of a seizure. This could also lead to more specific treatment methods during seizure evolution, not just to counteract the onset of a seizure. Moreover, by automatically describing an entire epileptic episode, false positives could be avoided more easily. In fact, to the authors knowledge, only a small group of papers have studied the problem of stage identification. For instance, [36,37] used Gabor atom density (GAD) as a measure of ECoG signal complexity, from signal segments, or epochs, of 2 s duration. Considering that ECoG signals, in general, are composed of multiple atoms of different frequency bandwidth and duration, the atoms have a different bandwidth-duration ratio (aspect ratio) [37]. Finally, [38] employs the 0–1 test for chaos to discriminate between each epileptic stage. Their approach discriminates between the Pre-Ictal and Ictal stages, but cannot successfully discriminate the Pre-Ictal and the Post-Ictal stages. Note that most papers cited above present results obtained with the best empirical configuration. While such an analysis is important, it does not consider the fact that different patients exhibit different signal patterns, even if all share a similar underlying structure [8]. In other words, while each stage is identifiable when you analyze the time-series of a seizure as a whole, if only a single 1-s segment is considered, for example, then determining the underlying stage is not trivially done. Moreover, given the stochastic nature of GP, and of all heuristic learning methods, the performance of such systems will greatly depend on the particular system configuration and data preprocessing. Therefore, the aim of the experimental work in this paper is to illustrate the statistical tendencies of the proposed approach, and provide a wide overview

of how the performance of the algorithm varies with, or depends on, several system parameters.

4. Genetic programming Evolutionary computation (EC) as a field, deals with the development and analysis of meta-heuristics (and hyper-heuristics) that are based on a simplified model of Neo-Darwinian evolution. These are population-based search methods, where candidate solutions are stochastically selected and varied to produce new solutions for a specified problem. This process is carried out iteratively until a predefined termination criterion is met. In general, to implement and execute an evolutionary algorithm (EA) the following aspects must be defined based on domain-specific knowledge. 1. An encoding scheme to represent and manipulate candidate solutions for a given problem. Each instantiated solution is referred to as an individual, and the set of solutions manipulated by the EA is called the population. 2. An evaluation function f that measures the quality of each solution based on the high-level goal of the problem, this function assigns a fitness value to each solution generated by the EA. 3. In order to produce new solutions, the EA applies variation operators,1 that take one or more solutions from the population as input, referred to as parents, and produce one or more solutions as output, referred to as offspring. The most common variation operators are mutation and crossover. Mutation stochastically modifies a single solution to produce a single offspring, while crossover takes two solutions as input and randomly swaps information between them to generate two new offspring. 4. Solutions are chosen by the variation operators based on their fitness using a predefined selection mechanism. Selection favors individuals with a better fitness value. 5. Finally, a survival strategy decides which individuals within the population will appear in the following iteration or generation. Among EAs, GP is arguably the most advanced technique, which can be used for automatic program induction [10]. In standard GP each individual is represented by a syntax-tree, because such structures can efficiently express simple computer programs, functions, or mathematical operators. Tree nodes contain a single element from a specified finite set of primitives P ¼ T⋃F. Leaf nodes contain elements from the set of terminals T, which normally correspond to inputs, while internal nodes contain elements from the set of functions F, which are the basic operations used to build more complex expressions. In essence, the set P defines the nature of the underlying search space for the evolutionary search, and even when a maximum depth or size limit for individual trees is enforced, normally the search space is very large but finite. In GP, the variation operators work as follows. During mutation, a single node is selected randomly, and the subtree rooted at that node is deleted and substituted by a randomly generated subtree; this operators is called subtree-mutation. Conversely, subtree-crossover randomly selects a subtree from the two parents and interchanges them. Selection in GP is normally done using a k-tournament, where each parent used by a genetic operator is the best individuals among k individuals chosen uniformly at random from the population. Finally, after 1 These are also referred to as genetic operators, or more generally as search operators.

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

1717

Fig. 3. General evolutionary loop executed by a tree-based GP algorithm. The image depicts how individuals are expressed as trees, how individuals are evaluated and selected as parents, and how the genetic operators are applied to generate new offspring.

each iteration the new population contains all of the offspring produced by the variation operators and the best solution found thus far (what is referred to as elitism). Fig. 3 presents a graphical depiction of the general stages of a GP-based search, illustrating how individuals are represented, evaluated and selected, as well as how the genetic operators are used to generate new candidate solutions.

5. Problem statement and proposed solution In this paper, the goal is to detect the three main stages of an epileptic seizure (Pre-Ictal, Ictal and Post-Ictal) given a short segment of an ECoG signal. This goal is posed as a classification problem, where the signal segment represents a pattern x A Rp , with p the total number of sample points which is dependent on the sampling rate and the signal duration. For instance, since the sampling rate during recording is 256 Hz, if we take a 2 s signal then p ¼512. Then, this can be defined as a supervised learning problem where a training set X of n-dimensional patterns with a known classification are used to derive a mapping function gðxÞ : Rp -M, where M are the three distinct epilepsy stages. The proposal in this work is to solve this problem using GP to derive the mapping function g. GP can be used in different ways to solve a supervised classification tasks such as the one presented here, see for instance [9,39]. However, in this work we test two GP classifiers that have achieved good results in difficult real-world domains [40]. 5.1. Static range selection GP classifier The first approach is called the Static Range Selection GP Classifier or SRS-GPC. In this approach, R is divided into M nonoverlapping regions, one for each class. Then, GP evolves a mapping gðxÞ : Rp -R, such that the region in R where pattern x is mapped to, determines the class to which it belongs. The fitness function is simple, it consists on maximizing the total classification accuracy of g. For the present problem, since there are three classes (the three seizure stages), R is divided into the following three ranges for each stage: Pre-Ictal ðinf ; 1Þ, Ictal [  1,1] and

Table 1 Parameters for the PGPC system used in the experimental tests. Parameter

Description

Population size Generations Initialization Operator probabilities Function set

200 individuals 200 generations Ramped half-and-half, with 6 levels of maximum depth Crossover pc ¼ 0:8; Mutation pμ ¼ 0:2 fþ ; ; n; =; pffi; sin ; cos ; log ; xy ; j  j; if g

Terminal set Bloat control Initial dynamic depth Hard maximum depth Selection

xμ , xm, xs , xmax, xmin, xs and xk Dynamic depth control 6 levels 20 levels Lexicographic parsimony tournament Keep best elitism

Survival

Post-Ictal ð1; inf Þ. This is a very simple and straightforward GP implementation, that is easy to setup and use. However, an obvious shortcoming is that it requires an a priori definition of the order and size of the region boundaries. 5.2. Probabilistic GP classifier We refer to the second approach as the Probabilistic GP Classifier or PGPC [40,41]. In PGPC, it is assumed that the behavior of the mapping function g can be modeled using multiple Gaussian distributions, each corresponding to a single class [40]. The distribution of each class N ðμ; sÞ is derived from the examples provided for it in set X , by computing the mean μ and standard deviation s of the outputs obtained from g on these patterns. Then, from the distributions generated for each class a fitness measure can be derived using Fisher's linear discriminant; for a two class problem it proceeds as follows. After the Gaussian distributions N i and N j are derived for class i and class j respectively, a Zhang and Smart [40] employ a distance measure between both classes given by di;j ¼

jμi μj j ; si þ sj

ð1Þ

where μi and μj are the means of the Gaussian distribution of each class, and si and sj their respective standard deviations. When this

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

0.6 0.5 0.4 0.3 0.2 0.1

2

er

0.4 0.3 0.2 0.1 5

4

3

2

1 0

20

60

Ov

e

80

0.4 0.3 0.2 0.1 3

2

Day

1 0

20

60

Ov

erl

80

0.4 0.3 0.2 0.1 4

3

Day

2

1

0

20

40

60

Ov

er

80

lap

20

Ov

er

0.3 0.2 0.1 4

3

2

1

0

20

80

60

40

Ov

er

0.5 0.4 0.3 0.2 0.1 4

3

2

1

0

20

60

Ov

er

0.4 0.3 0.2 0.1 4

3

Day

2

1

0

20

Ov

0.1 5

4

3

60

er

lap

80

2

1 0

20

80

60

40

Ov

erl

ap

0.7 0.6 0.5 0.4 0.3 0.2 0.1 5

4

3

2

1 0

20

80

60

40

Ov

erl

ap

0.7 0.6 0.5 0.4 0.3 0.2 0.1 4

3

Day

0.5

40

0.2

0 5

lap

0.6

5

0.3

0

80

0.7

0

0.4

Day

0.6

40

0.5

lap

0.7

5

0.6

Day

0.4

0

Classification Error

0.5

0

0.7

0

lap

0.5

Day

0.6

1

40

0.6

ap

0.7

0 5

2

60

80

0.7

5

Classification Error

0.5

4

3

Day

0.6

40

4

0

p rla

0.7

0 5

0.1

Day

0.5

40

0.2

5

Classification Error

Classification Error

Ov

0.6

Day

Classification Error

20

0.3

0

lap

0.7

0

Classification Error

1 0

0.4

Classification Error

3

Day

40

0.5

Classification Error

4

60

0.6

Classification Error

0 5

80

0.7

Classification Error

0.7

Classification Error

Classification Error

1718

2

1 0

20

60

40

80

ap

erl

Ov

0.7 0.6 0.5 0.4 0.3 0.2 0.1 0

5

4

3

Day

2

1

0

20

40

Ov

80

60

er

lap

Fig. 4. Median classification error plotted against day used for training and amount of overlap for the SRS-GPC. Each row corresponds to each subject (top row is S1 and bottom row S4) and each column corresponds to a different segment duration; leftmost row is 1 s, middle row 2 and rightmost 3. (a) S1, 1s, (b) S1, 2s, (c) S1, 3s, (d) S2, 1s, (e) S2, 2s, (f) S2, 3s, (g) S3, 1s, (h) S3, 2s, (i) S3, 3s, (j) S4, 1s, (k) S4, 2s and, (l) S4, 3s.

measure tends to 0, it is the worst case scenario because the mapping of both classes overlaps completely, and when it tends to 1 it represents the optimal case with a maximum separation. Moreover, the class separation measure is normalized by d^i;j ¼

1 ; 1 þ di;j

ð2Þ

which gives an estimation of the separation that the mapping function g induces between classes i and j. For a multiclass problem, with three or more classes, fitness is determined by considering the overlap between all pairs of classes.

Therefore, if c is the total number of classes then there are ð2c Þ overlap regions. Then, the fitness for a mapping g is given by the average d^ computed over all overlap regions. In particular, for the three class problem considered in this work, with the Pre-Ictal (1), Ictal (2) and Post-Ictal (3) classes, fitness is given by   ^ þd^ þd^ Þ 1 : f ðgÞ ¼ ðd1;2 1;3 3;2 3

ð3Þ

Finally, after executing the GP, the best individual found determines the parameters for the Gaussian distribution N i

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

0.5 0.4 0.3 0.2 0.1 4

3

2

Day

1 0

20

60

40

Ov

erl

80

ap

0.2 0.1 4

3

2

1

0

20

60

40

Ov

80

0.1 4

3

2

1

0

20

80

60

40

Ov

erl

e

p rla

0.4 0.3 0.2 0.1 4

3

2

Day

1

0

20

Ov

0.3 0.2 0.1 4

3

Day

2

1 0

20

40

Ov

60

erl

80

0.6 0.5 0.4 0.3 0.2 0.1 0

ap

5

4

3

Day

0.7

2

1 0

20

40

60

Ov

er

80

0.4 0.3 0.2 0.1 4

3

Day

2

1 0

20

40

Ov

60

erl

80

ap

Classification Error

0.5

0.5 0.4 0.3 0.2 0.1 4

3

Day

2

1 0

20

40

Ov

0.1 5

4

3

60

erl

2

ap

80

1 0

20

40

Ov

80

60

erl

ap

0.7 0.6 0.5 0.4 0.3 0.2 0.1 4

3

2

1 0

20

40

Ov

60

er

80

lap

0.7 0.6 0.5 0.4 0.3 0.2 0.1 5

4

3

2

Day

0.6

0 5

0.2

0

lap

0.7

0.6

0.3

Day

Classification Error

0.4

Classification Error

0.5

0.4

0 5

p rla

0.7

0.6

0.5

0

80

60

e

0.6

Day

0.5

40

0.7

ap

0.6

0 5

0.7

Classification Error

0.2

Classification Error

0.3

Classification Error

Classification Error

0.4

Day

Classification Error

0.3

0.7

0.5

0 5

0.4

Day

0.6

0 5

0.5

0 5

0.7

0 5

0.6

Classification Error

0 5

Classification Error

0.7

0.6

Classification Error

Classification Error

0.7

1719

1 0

20

80

60

40

ap

erl

Ov

0.7 0.6 0.5 0.4 0.3 0.2 0.1 0

5

4

3

Day

2

1 0

20

40

Ov

60

erl

80

ap

Fig. 5. Median classification error plotted against day used for training and amount of overlap for the PGPC. Each row corresponds to each subject (top row is S1 and bottom row S4) and each column corresponds to a different segment duration; leftmost row is 1 s, middle row 2 and rightmost 3. (a) S1, 1s, (b) S1, 2s, (c) S1, 3s, (d) S2, 1s, (e) S2, 2s, (f) S2, 3s, (g) S3, 1s, (h) S3, 2s, (i) S3, 3s, (j) S4, 1s, (k) S4, 2s and, (l) S4, 3s.

associated to each class. Then, a new test pattern x is assigned to class i when N i gives the maximum probability.

6. Experiments and results The goal of the experimental work is to evaluate the accuracy of intra-subject classification of epilepsy signals. In other words, to test the performance of classifiers that are trained and tested with signal samples from a single test subject that are recorded in different days.

6.1. Data sets The signals from epileptic episodes were recorded from four test subjects (four rodents) on whom the seizures are elicited and the signals recorded; call them subjects S1, S2, S3 and S4. For each subject, a level-5 seizure is elicited and recorded on five consecutive days; call them Day1, Day2, Day3, Day4 and Day5. Afterwards, the signal is classified manually by a human expert, who specifies where each epilepsy stage begins and ends. This manual classification establishes the ground-truth for the supervised-learning problem. The signal is divided into N segments, each constituting a sample from the corresponding stage.

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

We test three different segment lengths: 1, 2 and 3 s. When the signal is divided, we allow for different amounts of overlap between consecutive segments given by a percentage of signal duration. Five different overlaps are tested: 0%, 20%, 40%, 60% and 80%. Signal segments that lie on two adjacent stages are removed from the data-sets. In total, this gives 600 different experimental configurations, and by performing 30 runs for each configuration, that is a total of 18,000 independent runs used to test the proposed classifiers. The number of data samples for each class depends on how the signal is analyzed. For instance, if we consider signal segments of 1 s, on average that would produce 30 samples for the Pre-Ictal stage, and 19 for the Ictal stage if overlap is not allowed, and for 2 s segments the number of samples would be half those amounts. Similarly, by increasing the amount of overlap considered, the number of data samples is roughly increased by the same percentage. Notice then that the problem is slightly imbalanced, since the Ictal stage is shorter than the other two. 6.2. GP implementation Both GP classifiers use a standard Koza style tree-based representation, with subtree-crossover and sub-tree mutation. The basic parameter values of both systems are presented in Table 1. For both GP classifiers, the terminal elements are basic statistical features computed for each signal segment x that is classified. Specifically, the terminal set T contains the following features: mean value xμ , median xm, standard deviation xs , maximum xmax, minimum xmin, skewness xs and kurtosis xk. For each experimental configuration, 30 independent runs are of each GP classifier are executed. This is necessary, since GP is a nondeterministic search process. 6.3. Intra-subject classification As stated above, the task of automatic epileptic stage identification is posed a supervised learning problem. The training data consists of all of the signal segments taken from a single recording day. Then, for intra-subject classification the evolved classifiers are tested on the signal segments recorded on each additional day. This is done for all subjects (S1, S2, S3 and S4) and all possible combinations of training and testing days, as well as for all segment durations (1, 2 or 3 s) and overlaps (0%, 20%, 40%, 60% and 80%). Fig. 4 shows the median classification error achieved by the SRS-GPC over all 30 runs, by considering the best solution found in each run. The error is plotted with respect to the day used to train the classifier and the amount of overlap. Each row corresponds to each subject (top row is S1 and so on) and each column corresponds to a different segment duration (leftmost row is 1 s and so on). Fig. 5 represents the same for PGPC. The plots shown in Figs. 4 and 5 represent an extensive empirical evaluation of the proposed approach for stage identification. The best way to read these plots is as follows. In all cases, the flatness of the plotted surfaces is related to the robustness of the proposed approach to the different experimental configurations. These figures illustrate strong statistical trends regarding the behavior of the system, from which several conclusions can be derived. First, overall the best performance is achieved by the simpler SRS-GPC classifier, compared with PGPC. This was slightly unexpected, given the strong assumptions made by the SRS-GPC method. Nonetheless, across all configurations and all test subjects, SRS-GPC exhibits a lower median error on most tests. Second, it appears that the segment length is an important determining factor in classification performance. In particular, smaller segments (1 or 2 s) are easier to classify than longer ones (3 s). Also, it appears that a higher amount of overlap between

PRE−ICTAL ICTAL POS−ICTAL

Probability

1720

−1

−0.8

−0.6

−0.4

−0.2

0

0.2

0.4

0.6

0.8

1

L(x) Fig. 6. The probability density function associated to each class (Pre-Ictal, Ictal and Post-Ictal) by classifier L; the parameters of each Normal distribution are derived from the training data recorded on the Day 1.

segments induces better results. A reasonable explanation for both observations can be given based on the machine learning approach followed by the proposal. In particular, shorter segments and a larger overlap produce more training data with which to train the classifier and produce a better supervised learning process. Third, while the results described above hold for all test subjects, there are still some subtle and important differences among them. The GP classifiers achieve good results for three of the test subjects (S1, S2 and S3). However, the performance is considerably worse for test subject S4. This result is coherent with the general observation that different subjects can produce quite different signal patterns, even if all signals appear similar at a coarser scale. Moreover, the plots show that classifier performance can sometimes vary depending upon which training day is used. In effect, this shows that, of course, the recorded signals should not be expected to be homogeneous given the different (even if slight) conditions under which each recording session is carried out. Finally, considering results obtained in previous work using GADs [36,37], where stage identification accuracy was as low as 15% and never above 25% for any particular stage, here we find that the ability to identify epilepsy stages are substantially better with the GP classifiers. The results suggest that the time domain properties used in GP were more efficient for stage classification than the signal complexity computed from the time–frequency domain. Moreover, such features, and the evolved classifiers, are of low computational complexity and can be easily implemented on dedicated hardware, a prototype of which is presented in the following section. 6.4. Prototype of a hardware implementation for stage identification In the previous section we analyzed the overall statistical performance and tendencies of the GP-based classifiers, describing the behavior of the system based on thousands of runs. Nonetheless, for a real-world scenario only a single classifier would be required and used. Therefore, in this section a fine grained look is taken, choosing a single evolved classifier and describing how it operates, what features are used, and analyzing its performance in detail. Moreover, we are interested in assessing how these classifiers would fare in a real embedded system, where the signal is read and processed in real time. To this end, a working prototype is developed and implemented using a field programmable gate array (FPGA). 6.4.1. The evolved classifier Given the results described in the preceding section, there are some clear tendencies regarding the best parametrization for the

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

1721

Fig. 7. Block diagram of the experimental setup and design for the proposed embedded system.

Table 2 Statistical performance of the evolved mapping function L tested on the embedded system. Measure

Sensitivity Specificity

Stage Pre-Ictal

Ictal

Post-Ictal

0.9761 0.9940

0.6561 0.9172

0.8305 0.8609

depicted in Fig. 6. The parameters of each distribution are: N 1 ðμ1 ¼ 0:23; s1 ¼ 0:13Þ; N 2 ðμ2 ¼ 0:64; s2 ¼ 0:03Þ and N 3 ðμ3 ¼ 0:48; s3 ¼ 0:07Þ.

Fig. 8. The embedded system implemented on the ELVIS II development platform from National Instruments. On the left is the FPGA board, on the right is the DAQ card and on the top is the laptop with development software.

proposed system. For instance, the best performance was achieved using the SRS classifier and segments of 2 s duration. However, we are not interested in testing the best possible configuration, the goal is to assess the performance of an “average” classifier. Therefore, an evolved classifier is chosen from the set of 30 runs that were executed with the following configuration: the PGPC classifier trained for subject 1 using the signal from Day1, segments of 1 s in duration and considering 0% overlap. The best solutions from all 30 runs were sorted based on training fitness and a classifier was chosen from the top 50 percentile. The best solution is not taken however, rather the individual expressed by the most compact (or computationally simple) expression. The reason for this choice is based on the fact that the goal is to implement a real embedded system, were computational costs must be maintained as low as possible. The choice provides a good compromise between performance and solution complexity. The chosen classifier, hereafter referred to by the mapping function L,2 is given by LðxÞ ¼ sin ð sin ð sin ð cos ðmaxðxÞðminðxÞ þ1ÞÞÞÞÞ;

ð4Þ

where x is the signal segment. Notice the simplicity, yet somewhat unorthodox operation; i.e., it is evident that the classifier is quite simple computationally, but its interpretability is in fact quite low. Given the training data, the probability distributions for each class, Pre-Ictal ðN 1 Þ, Ictal ðN 2 Þ and Post-Ictal ðN 3 Þ, are graphically 2 In fact, the syntactic expression evolved by GP was algebraically simplified for ease of presentation and implementation on the embedded system.

6.4.2. Experimental setup and embedded system The goal is to test the classifier in the following scenario: to process a complete epilepsy signal online, classifying each signal segment in real time. However, as stated before, performing experimental tests with real epilepsy signals as they occur is highly impractical. Therefore, the experimental setup used to simulate this real-world scenario consists of two main modules: (A) signal reproduction and (B) the embedded system used to classify the data; Fig. 7 provides a graphical depiction of the system. To reproduce the epilepsy signals, a digital-to-analog converter takes the previously recorded signals as input and outputs the reconstructed analogue signal; in our tests the DAQ USB-6211 card by National Instruments is used. Afterwards, the embedded system is implemented on the Digital Electronics FPGA Board and the ELVIS II development platform with a Spartan 3E FPGA from National Instruments. The embedded system is composed of three sub-modules, respectively these are for: (i) reading the signal; (ii) computing statistical features; and (iii) performing the classification. The signal is read through the analogue-todigital converter on the FPGA board, and all of the computation is done on the FPGA programmed using the LabVIEW system-design platform. It is important to point out that the classifiers were originally evolved using Matlab, so they had to be ported completely to LabVIEW, given the restrictions imposed by the FPGA environment. The complete system is shown in Fig. 8.

6.4.3. Results The experiments for this classifier replicate those described in Section 6.3 for intra-subject stage classification, where the classifier is trained on the signal from one seizure, and then tested on all the other signals from the same patient. However, the signal reconstruction process does add some noise, thus the classification process is repeated 30 times and the average is used to asses classification performance. The average classification error is 14.55%, or an 84.45% accuracy, and system response was achieved in real-time. To provide a more detailed description, we employ the sensitivity Se and specificity Sp measures, given by Se ¼

TP ; TP þ FN

ð5Þ

1722

Sp ¼

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

TN ; TN þ FP

ð6Þ

where TP, TN, FP and FN respectively represent the true positives, true negatives, false positives and false negatives of the confusion matrix for each class. Results based on these measures are summarized in Table 2, that shows the average sensitivity and specificity achieved with respect to each class. Notice that while the performance is quite good, this classifier does not represent the best possible solution found by GP. Nonetheless, the results are meant to confirm the following: (1) it is practical to port the evolved classifiers and to implement them on a real system; (2) the results are consistent with the average performance of the GP search; (3) and the classifiers are efficient and simple, able to process the signal in real-time.

7. Conclusions This paper presents an automatic method for identifying the three main stages of an epileptic seizure from ECoG signals. The proposal is based on posing the problem as a supervised learning problem and solving it with GP. The results exhibit strong statistical tendencies of the GP classifiers that suggest that the approach is able to solve the intra-patient classification problem. These results are unique and show substantial improvement when compared with previous methods [36–38]. Moreover, the classifiers are ported to an FPGA embedded system, and tested in a realistic online scenario. The results are encouraging, since the GP classifier was simple and efficient, and the performance was consistent with that achieved entirely in simulation, while operating in real-time. Future work will center on extending the test data, particularly to focus on signals from human patients and to consider other recording strategies, particularly EEG data. Moreover, since the variables used to classify the signals are simple time-domain features, and the classifiers themselves are constructed from a set of simple mathematical functions, the computational system can easily be implemented on a dedicated embedded system, opening a pathway to develop automatic systems that can detect and track the evolution of a seizure over all of its stages. Such a system can aid researchers in their understanding and analysis how seizures develop over time in response to new treatment methods.

Conflict of interest statement None declared.

Acknowledgments Funding for this work provided by CONACYT (Mexico) Basic Science Research Project no. 178323 and DGEST (Mexico) Research Project no. TIJ-ING-2012-110. Fifth author is supported by a CONACYT (Mexico) doctoral scholarship no. 226981. Thanks are extended to Francisco Sancho from Hospital Universitario de Valencia, for his collaboration and support during the signal recording process. Finally, thanks are given to Moises Zonta, Ivan Garcia and Enrique Naredo from Instituto Tecnológico de Tijuana for their collaboration and support in the development of experimental work and graphical content of this paper. References [1] O.C. Cockerell, Epilepsy, current concepts, 2003. [2] B. Marchesi, A.L. Stelle, H.S. Lopes. Detection of epileptic events using genetic programming. IEE, vol. 3, October 1997, pp. 1198–1201.

[3] J. Barcia, P. Rubiuo, Anticonvulsant and neurotoxic effects of intracerebroventricular injection of phenytoin, phenobarbital and carbamazepine in an amygdala-Kindling model of epilepsy in the rat, Epilepsy Res. 33 (1999) 159–539. [4] H. Witte, L.D. Iasemidis, B. Litt, Special issue on epileptic seizure prediction, IEEE Trans. Biomed. Eng. 50 (2003) 537–539. [5] B. Litt, J. Echauz, Prediction of epileptic seizures, Lancet Neurol. 1 (May (1)) (2002) 22–30. [6] P.J. Franaszczuk, G.K. Bergey, Time–frequency analysis using the matching pursuit algorithm applied to seizures originating from the mesial temporal lobe, Electroencephalogr. Clin. Neurophysiol. 106 (1998) 513–521. [7] C. Bigan, W.S. Woolfson, Time–frequency analysis of short segments of biomedical data, IEE Proc. Sci. Meas. Technol. 147 (6) (2000) 368–373. [8] P.J. Franaszczuk, Gregory K. Bergey, Time–frequency analysis using the matching pursuit algorithm applied to seizures originating from the mesial temporal lobe, Electroencephalogr. Clin. Neurophysiol. 106 (6) (1998) 513–521. [9] John R. Koza, Genetic Programming II: Automatic Discovery of Reusable Programs, MIT Press, Cambridge, MA, USA, 1994. [10] William B. Langdon, Riccardo Poli, Foundations of Genetic Programming, Springer Publishing Company, Incorporated, 1st ed., 2010. [11] Arturo Sotelo, Enrique Guijarro, Leonardo Trujillo, Luis Coria, Yuliana Martínez, Analysis and classification of epilepsy stages with genetic programming, in: Oliver SchÃijtze, Carlos A. Coello Coello, Alexandru-Adrian Tantar, Emilia Tantar, Pascal Bouvry, Pierre Del Moral, Pierrick Legrand (Eds.), EVOLVE, Advances in Intelligent Systems and Computing, vol. 175, Springer, 2012, pp. 57–70. [12] W. Loscher, C. Rundfeldt, Kindling as a model of drug-resistant partial epilepsy: selection of phenytoin-resistant and non-resistant rats, J. Pharmacol. 258 (1991) 438–489. [13] D.A. Coulter, D.C. McIntyre, W. Loscher, Animal models of limbic epilepsies: what can they tell us? Brain Pathol. 2 (12) (2002) 240–256. [14] D.M. Durand, M. Bikson, Suppression and control of epileptiform activity by electrical stimulation: a review, Proc. IEEE 89 (7) (2001) 1065–1082. [15] K. Morimoto, M. Fahnestock, R.J. Racine, Kindling and status epilepticus models of epilepsy: rewiring the brain, Prog. Neurobiol. 73 (May) (2004) 1–60. [16] R.J. Racine, Modification of seizure activity bye electrical stimulation. ii motor seizure, Clin. Neurophysiol. 32 (3) (1972) 281–294. [17] E. Bertram, Time-frequency analysis using the matching pursuit algorithm applied to seizures originating from the mesial temporal lobe. The relevance of Kindling for human epilepsy, Epilepsia 48 (2) (2007) 65–74. [18] J. Barcia, P. Rubiuo, Anticonvulsant and neurotoxic effects of intracerebroventricular injection of phenytoin, phenobarbital and carbamazepine in an amygdala-Kindling model of epilepsy in the rat, Epilepsy Res. 33 (1999) 159–167. [19] M.M. D'Alessandro, G.J. Vachtsevanos, R.Esteller, J. Echauz, A. Koblasz, B. Litt, Spectral entropy and neuronal involvement in patients with mesial temporal lobe epilepsy, in: International Conference on Mathematics and Engineering Techniques in Medicine and Biological Sciences, 2000. [20] H.P. Zaveri, Time frequency representation of electrocortigrams in temporal lobe epilepsy, IEEE Trans. Biomed. Eng. 39 (1992) 502–509. [21] C. Jouny, P. Franaszczuk, G. Bergey, Characterization of epileptic seizure dynamics using Gabor atom density, Clin. Neurophysiol. 114 (3) (2003) 426–437. [22] L. Vezard, M. Chavent, P. Legrand, F. Faita-Ainseba, L. Trujillo, Detecting mental states of alertness with genetic algorithm variable selection, in: IEEE Congress on Evolutionary Computation, CEC'13, 2013. [23] W. Loscher, E. Reissmüller, Anticonvulsant effect of fosphenytoin in amygdalakindled rats: comparison with phenytoin, Epilepsy Res. 30 (1998) 69–76. [24] M. Jeub, H. Beck, E. Sie, C. Ruschenschmidt, E.J. Speckmann, U. Ebert, H. Potschka, C. Freichel, E. Reissmuller, W. Loscher, Effect of phenytoin on sodium and calcium currents in hippocampal ca1 neurons of phenytoinresistant Kindled rats, Neuropharmacology 42 (1) (2002) 107–116. [25] M. Teplan, Fundamentals of EEG measurement, Meas. Sci. Rev. 2 (2) (2002) 1–11. [26] S. Mingui, M. Scheuer, Time–frequency analysis of high-frequency activity at the start of epileptic seizures, Proc. IEEE/EMBS 3 (1997) 1184–1187. [27] A. Chiu, S.S. Jahromi, H. Khosravani, P.L. Carlen, B.L. Bardakjian, The effects of high-frequency oscillations in hippocampal electrical activities on the classification of epileptiform events using artificial neural networks, J. Neural Eng. 3 (1) (2006) 9–20. [28] L. Iasemidis, D.S. Shiau, J.C. Sackellares, P.M. Pardalos, A. Prasad, Dynamical resetting of the human brain at epileptic seizures: application of nonlinear dynamics and global optimization techniques, IEEE Trans. Biomed. Eng. 51 (3) (2004) 493–506. [29] M. D'Alessandro, R. Esteller, G. Vachtsevanos, A. Hinson, J. Echauz, Brian Litt, Epileptic seizure prediction using hybrid feature selection over multiple intracranial EEG electrode contacts: a report of four patients, IEEE Trans. Biomed. Eng. 50 (5) (2003) 603–615. [30] J.C. Sackellares, Seizure prediction, Epilepsy Curr. 8 (3) (2008) 55–59. [31] Heitor S. Lopes, Genetic programming for epileptic pattern recognition in electroencephalographic signals, Appl. Soft Comput. 7 (1) (2007) 343–352. [32] Otis Smart, Hiram Firpi, George Vachtsevanos, Genetic programming of conventional features to detect seizure precursors, Eng. Appl. Artif. Intell. 20 (8) (2007) 1070–1085. [33] Ling Guo, Daniel Rivero, Julián Dorado, Cristian R. Munteanu, Alejandro Pazos, Automatic feature extraction using genetic programming: an

A. Sotelo et al. / Computers in Biology and Medicine 43 (2013) 1713–1723

[34]

[35]

[36]

[37]

[38]

[39]

[40]

[41]

[42]

application to epileptic eeg classification, Expert Syst. Appl. 38 (8) (2011) 10425–10436. B. Litt, Implantable devices for epilepsy: a clinical perspective, in: Engineering in Medicine and Biology, 2002. Proceedings of the Second Joint 24th Annual Conference and the Annual Fall Meeting of the Biomedical Engineering Society EMBS/BMES Conference, vol. 3, 2002, pp. 2035–2036. M.T. Salam, M. Sawan, Dang Khoa Nguyen, A novel low-power-implantable epileptic seizure-onset detector, IEEE Trans. Biomed. Circuits Syst. 5 (6) (2011) 568–578. A. Sotelo, E. Guijarro, M.J. García, C.E. Vázquez, Epoch parameterization by Gabor atom density in experimental epilepsy, in: Fourth International Conference on Electrical and Electronics Engineering, September 2007, pp. 61–64. A. Sotelo, E. Guijarro, M.J. García, C.E. Vázquez, Epilepsy stages diagnosis by Gabor atom density according to their aspect ratio, in: Fourth International Conference on Electrical Engineering, Computing Science and Automatic Control, November 2008, pp. 158–163. Arturo Sotelo, Enrique Guijarro, Luis N. Coria, Leonardo Trujillo, Paul A. Valle, Epilepsy ictal stage identification by 0–1 test of chaos, in: Proceedings of the Third IFAC CHAOS Conference, CHAOS 2012, IFAC, 2012, pp. 223–228. Jeroen Eggermont, Joost N. Kok, Walter A. Kosters, Genetic programming for data classification: partitioning the search space, in: Proceedings of the 2004 ACM Symposium on Applied Computing, SAC '04, New York, NY, USA, ACM, 2004, pp. 1001–1005. Mengjie Zhang, Will Smart, Using gaussian distribution to construct fitness functions in genetic programming for multiclass object classification, Pattern Recognition Lett. 27 (August) (2006) 1266–1274. Leonardo Trujillo, Yuliana Martínez, Edgar Galván-López, Pierrick Legrand, Predicting problem difficulty for genetic programming applied to data classification, in: Proceedings of the 13th Annual Conference on Genetic and Evolutionary Computation, GECCO '11, New York, NY, USA, ACM, 2011, pp. 1355–1362. G. Paxinos, C. Watson, The Rat Brain in Stereotatic Coordinates, 4th ed., Academic Press, Sydney, 1986.

Arturo Sotelo received a degree in Electronic and Communications Engineering from the Escuela Superior de Ingeniería Mecánica y Eléctrica del Instituto Politécnico Nacional in 1990, and a Master's degree in Science in the field of Digital Systems in 1997 from the Centro de Investigación en Desarrollo de Tecnología Digital del Instituto Politécnico Nacional, México. He is currently a Ph.D. candidate from the Universidad Politécnica de Valencia, Spain, and a research professor at ITT. He has worked in automation and robotics, and is currently involved in biomedical engineering and interpretation of brain signals from epileptic patients, the topic of his doctoral dissertation.

Enrique Guijarro has a doctorate in Industrial Engineering from Universidad Politécnica de Valencia, Spain. Currently, he is a titular professor at the Electronica Engineering Department in UPV, and a researcher at the Instituto Interuniversitario de Investigación en Bioingeniería y Tecnología Orientada en el SerHumano (i3BH).

1723

Leonardo Trujillo received an Electronic Engineering (2002) and a Masters in Computer Science (2004) from the Technical Institute of Tijuana in México. He then received a doctorate in Computer Science from CICESE research center, Mexico (2008), developing Genetic Programming (GP) applications for Computer Vision problems; focusing on feature extraction and image description. He is currently professor at the Technical Institute of Tijuana in México (ITT), where he is President of the Doctorate Program and head of the Cybernetics research group. He actively collaborates with other institutions, particularly the University of Bordeaux and INRIA in France where he has received the distinction of invited professor and researcher in 2010 and 2011. Moreover, he is actively involved with research work with the University of Extremadura in Spain and Trinity College Dublin in Ireland. Currently, he is involved in interdisciplinary research in the fields of Evolutionary Computation, Computer Vision, Image Analysis, Pattern Recognition and Autonomous Robotics. He is a level 1 member of the National Research System of México and head of a basic science project from CONACYT México, studying problem difficulty and prediction of expected performance for GP systems, developing new theory and applications. The work of him has led to the publishing of over 40 works in international journals and conferences, receiving several awards and distinctions.

Luis N. Coria is a researcher of the Cibernetics group at Instituto Tecnologico de Tijuana (ITT) and invited professor at Centro de Investigación y Desarrollo de Tecnología Digital (CITEDI-IPN), both in Tijuana, Mexico. He participates in the Robotics and Mechatronics network at Instituto Politecnico Nacional (IPN). He is Electronics Engineer (1999) and M.Sc. in Digital Systems (2005). His Ph.D. is in Electronics and Communications (2010) and now his interest is focused in Chaotic Systems, Nonlinear Systems Analysis and Signal Analysis with nonconventional methods. His recent contributions are on the analysis of Nonlinear Cancer and Biological models. Furthermore, currently he is focused on the analysis and characterization of brain signals from epileptic models.

Yuliana Martínez is first year student in the Ph.D. program in Engineering Sciences from the Instituto Tecnologico de Tijuana. She received a Masters degree in Computer Science from the Instituto Tecnologico de Tijuana (2009–2011), and an Engineering Degree in Computing from instituto Tecnologico de Los Mochis (2002– 2007). Her main area of research is genetic programming, in the study of one of the most important problems within the community, prediction of performance and difficulty of problems in GP.