Keywords: human TLE, hippocampus, intracerebral EEG, neuronal population model, parameter identification, ictogenesis mechanisms

model identification in human TLE are in agreement with some results obtained experimentally, either on animal models of epilepsy or on the human epileptic tissue.

Abstract

In human partial epilepsies, as well as in experimental models of chronic and/or acute epilepsy, the role of inhibition and the relationship between the inhibition and excitation and epileptogenesis has long been questioned. Besides experimental methods carried out either in vitro (human or animal tissue) or in vivo (animals), pathophysiological mechanisms can be approached by direct recording of brain electrical activity in human epilepsy. Indeed, in some clinical pre-surgical investigation methods like stereoelectroencephalopraphy (SEEG), intracerebral electrodes are used in patients suffering from drug resistant epilepsy to directly record paroxysmal activities with excellent temporal resolution (in the order of one millisecond). The study of neurophysiological mechanisms underlying such depth-EEG activities is crucial to progress in the understanding of the interictal to ictal transition.

In this study, we relate electrophysiological patterns typically observed during the transition from interictal to ictal activity in human mesial temporal lobe epilepsy (MTLE) to mechanisms (at a neuronal population level) involved in seizure generation through a computational model of EEG activity. Intracerebral EEG signals recorded from hippocampus in five patients with MTLE during four periods (during interictal activity, just before seizure onset, during seizure onset and during ictal activity) were used to identify the three main parameters of a model of hippocampus EEG activity (related to excitation, slow dendritic inhibition and fast somatic inhibition). The identification procedure used optimization algorithms to minimize a spectral distance between real and simulated signals. Results demonstrated that the model generates very realistic signals for automatically identified parameters. They also showed that the transition from interictal to ictal activity can not be simply explained by an increase in excitation and a decrease in inhibition but rather by time-varying ensemble interactions between pyramidal cells and local interneurons projecting to either their dendritic or perisomatic region (with slow and fast GABA A kinetics). Particularly, during preonset activity, an increasing dendritic GABAergic inhibition compensates a gradually increasing excitation up to a brutal drop at seizure onset when faster oscillations (beta and low gamma band, 15 to 40 Hz) are observed. These faster oscillations are then explained by the model feedback loop between pyramidal cells and interneurons targeting their perisomatic region. These findings obtained from

Introduction

In human partial epilepsies, as well as in experimental models of chronic and/or acute epilepsy, the role of inhibition and the relationship between the inhibition and excitation and epileptogenesis has long been questioned ( [START_REF] Dichter | Basic mechanisms of epilepsy: targets for therapeutic intervention[END_REF], [START_REF] Dalby | The process of epileptogenesis: a pathophysiological approach[END_REF]). Recent advances have suggested that different types of inhibition (depending on the location of synapses/receptors, for example) are present in a given neuronal tissue and that they may be differently involved in epileptic processes.

Particularly, in experimental focal models (kainate or pilocarpine-treated rat), it has been shown that GABAergic inhibition is impaired, but not uniformly: dendritic inhibition is reduced whereas somatic inhibition is preserved [START_REF] Cossart | Dendritic but not somatic GABAergic inhibition is decreased in experimental epilepsy[END_REF][START_REF] Houser | Vulnerability and plasticity of the GABA system in the pilocarpine model of spontaneous recurrent seizures[END_REF]. This intact inhibition could be due to selective survival of perisomatic inhibitory interneurons, as shown in human epileptogenic tissue (sclerotic hippocampus in temporal lobe epilepsy -TLE -) [START_REF] Wittner | Surviving CA1 pyramidal cells receive intact perisomatic inhibitory input in the human epileptic hippocampus[END_REF][START_REF] Wittner | Preservation of perisomatic inhibitory input of granule cells in the epileptic human dentate gyrus[END_REF].

Besides experimental methods carried out either in vitro (human or animal tissue) or in vivo (animals), pathophysiological mechanisms can be approached by direct recording of brain electrical activity in human epilepsy. Indeed, in some clinical pre-surgical investigation methods like stereoelectroencephalopraphy (SEEG) [START_REF] Bancaud | Methodology of stereo EEG exploration and surgical intervention in epilepsy[END_REF], intracerebral electrodes are used in patients suffering from drug resistant epilepsy to directly record paroxysmal activities with excellent temporal resolution (in the order of one millisecond). The study of neurophysiological mechanisms underlying such depth-EEG activities is crucial to progress in the understanding of the interictal to ictal transition. Here, we propose to perform this study using not only descriptive signal analysis techniques but also computational modeling, a research methodology that has been largely developed for the study of various (patho)physiological mechanisms and that attracts growing interest in neuroscience. We address the question whether EEG dynamics observed during the transition from interictal to ictal activity could result from variable global interactions between populations of pyramidal cells and inhibitory interneurons and whether these interactions could be revealed by appropriate model-based processing of intracerebral EEG signals. Hence, we relate electrophysiological patterns of seizure generation in human hippocampus to excitatory and inhibitory interactions in a neuronal population model of hippocampal EEG activity.

This type of computational model (also referred to as "neuronal rate model" or "neuronal population dynamics" [START_REF] Wilson | Excitatory and inhibitory interactions in localized populations of model neurons[END_REF] has been demonstrated to produce realistic epileptiform activities [START_REF] Wendling | Relevance of nonlinear lumped-parameter models in the analysis of depth-EEG epileptic signals[END_REF]. Recently, it has provided a means for linking some pathophysiological characteristics of hippocampus at a neuronal population level to the fast oscillations observed at seizure onset at the EEG level [START_REF] Wendling | Epileptic fast activity can be explained by a model of impaired GABAergic dendritic inhibition[END_REF].

Here, we report results about identification of the three main parameters of this model (namely excitation, slow dendritic inhibition and fast somatic inhibition) from recordings performed in patients with mesial TLE (inverse problem). Parameter identification (based on evolutionary algorithms) is performed in order to produce signals that match those actually recorded during interictal to ictal transition. From the analysis of parameters evolution during this transition, hypotheses are generated about excitatory and inhibitory interactions that take place between subsets of cells represented within the modeled hippocampal neuronal population. In particular, we found that during preonset activity (a few tens of seconds before seizure), an increasing dendritic inhibition (produced by interneurons mediating slow GABA A IPSCs) compensates a gradually increasing excitation up to a brutal drop at seizure onset. Faster ictal oscillations (beta and low gamma band) are then explained by the model feedback loop involving pyramidal cells and interneurons targeting their perisomatic region (mediating fast GABA A IPSCs). As described in the discussion, these findings obtained from model identification in human TLE are in agreement with some experimental studies conducted in animal models of epilepsy.

Material and methods

Intracerebral EEG signals recorded from hippocampus in patients with MTLE

Five patients undergoing presurgical evaluation of drug-resistant mesial temporal lobe epilepsy (MTLE) were selected. All patients had a comprehensive evaluation including detailed history and neurological examination, neuropsychological testing, routine magnetic resonance imaging (MRI), surface electroencephalography (EEG) and stereoelectroencephalography (SEEG). As illustrated in figure 1, SEEG recordings were performed using intracerebral electrodes placed intracranially according to Talairach's stereotactic method [START_REF] Bancaud | Methodology of stereo EEG exploration and surgical intervention in epilepsy[END_REF]. These patients were selected for the present study because they satisfied the following criteria: (i) seizures involved the mesial temporal regions at the onset; (ii) MRI was normal or suggestive of hippocampal sclerosis (HS) (hippocampal atrophy and increased T2 signal); (iii) one electrode was placed in the hippocampus; and iv) the electrophysiological pattern during the transition to seizure was characterized by the emergence of a low-frequency high-amplitude spiking activity followed by a socalled "rapid discharge" in the hippocampus. This pattern has been reported to be typical in MTLE [START_REF] Engel | Surgical treatment of epilepsy: opportunities for research into basic mechanisms of human brain function[END_REF], [START_REF] Spencer | Morphological patterns of seizures recorded intracranially[END_REF], [START_REF] Velasco | Functional and anatomic correlates of two frequently observed temporal lobe seizure-onset patterns[END_REF], [START_REF] Bartolomei | Pre-ictal synchronicity in limbic networks of mesial temporal lobe epilepsy[END_REF].

Analyzed signals were recorded on a 128 channel Deltamed TM system. They were sampled at 256 Hz and recorded on a hard disk (16 bits/sample) using no digital filter. Two hardware filters are present in the acquisition procedure. The first one is a high-pass filter (cut-off frequency equal to 0.16 Hz at -3 DB) used to remove very slow variations that sometimes contaminate the baseline. The second one is a 1 st order low-pass filter (cut-off frequency equal to 97 Hz at -3 DB) to avoid aliasing.

Definition of periods of interest

For the five patients, signals recorded from the hippocampus during the transition from interictal to ictal activity were visually analyzed. In each patient, four segments of EEG activity (10 seconds duration) were defined at different stages of this transition using visual criteria (figure 2 and figure 3), as follows. During the interictal period ("intICTAL" segment, chosen 1 minute before seizure onset at least), EEG signals do not exhibit high amplitude transient spikes. The preonset period (preceding seizure onset) is marked by the appearance of high amplitude spikes ("preONSET" segments, from 10 to 50 seconds before seizure onset). The electrical onset of the seizure ("ONSET" segment) is characterized by the appearance of the rapid discharge reflecting faster oscillations, typically in the beta and low gamma frequency bands (15-40 Hz). Finally, as seizure develops, the rapid discharge slows down and changes into a more rhythmic activity ("ICTAL" segment), typically in the theta or alpha frequency band (4 -10 Hz).

Computational model of hippocampus neuronal population and simulation of EEG signals

A macroscopic approach (neuronal population, adapted to the macroscopic size of SEEG electrodes) was used to build the EEG model which has been detailed in previous reports [START_REF] Wendling | Epileptic fast activity can be explained by a model of impaired GABAergic dendritic inhibition[END_REF].

Briefly, this model is based on the global cellular organization of the hippocampus. Firstly, recurrent excitatory connections from pyramidal cells to pyramidal cells have been demonstrated to occur in CA1 [START_REF] Thomson | Excitatory Connections Between CA1 Pyramidal Cells Revealed by Spike Triggered Averaging in Slices of Rat Hippocampus are Partially NMDA Receptor Mediated[END_REF], [START_REF] Whittington | Recurrent excitatory postsynaptic potentials induced by synchronized fast cortical oscillations[END_REF]. Secondly, a series of studies based on a variety of techniques [START_REF] Miles | Differences between somatic and dendritic inhibition in the hippocampus[END_REF] demonstrated that there are two types of GABA A synaptic responses in CA1 pyramidal neurons: a fast one near the soma and a slow one at the dendrites. The first one -GABA A,fast -is a rapidly activated and decaying IPSC mediated by somatic synapses and the second one -GABA A,slow -is a slowly rising and decaying IPSC mediated by dendritic synapses. More recent works [START_REF] Banks | The synaptic basis of GABAA,slow[END_REF][START_REF] White | Networks of interneurons with fast and slow gamma-aminobutyric acid type A (GABAA) kinetics provide substrate for mixed gammatheta rhythm[END_REF] suggested that two separate classes of interneurons (possibly basket cells and interneurons in stratum lacunosum-moleculare respectively called, for simplicity, GABA A,fast interneurons and GABA A,slow interneurons) give rise to these two IPSCs. Moreover, as suggested in [START_REF] Banks | Interactions between distinct GABA(A) circuits in hippocampus[END_REF], both classes interact: GABA A,slow cells inhibit not only pyramidal cells but also GABA A,fast interneurons.

The model was designed to represent this functional organization of interacting subsets of principal cells and interneurons, as summarized in figure 4-a. It consists in three subsets of neurons, namely the main cells (i.e. pyramidal cells), the slow dendritic-projecting inhibitory interneurons (GABA A,slow receptors) and the fast somatic-projecting inhibitory interneurons (GABA A,fast receptors). Interneurons receive an excitatory input (AMPA receptor-mediated) from pyramidal cells. The influence from neighboring or more distant populations is represented by an excitatory input (modeled by a positive mean gaussian white noise) that globally describes the average density of afferent action potentials. The model output corresponds to the post-synaptic activity of the subset of pyramidal cells that mainly contributes to the EEG signal (field potential).

) t ( n

In each subset, a linear transfer function is used to transform the average pre-synaptic pulse density of afferent action potentials (the input) into an average post-synaptic membrane potential (the output). This transfer function is of order 2. The associated impulse response , and , respectively determines the excitatory, slow 4-c). In this work, these key parameters are automatically identified from real EEG signals recorded in the hippocampus.

The inhibitory feedback loop from the subset of fast somatic-projecting inhibitory interneurons uses a faster impulse response (i.e. producing faster IPSP) than . Here, the important parameter is the average somatic time constant 1/g which was chosen lower than the average dendritic time constant 1/b, consistently with data provided in (Traub et al., 1999a) is also used to model threshold and saturation effects in the relationship between the average postsynaptic potential of a given subset and the average pulse density of potentials fired by the neurons.

] e /[ e ) v ( S ) v v ( r - + = 0 1 2 0
Finally, interactions between main cells and local neurons are summarized in the model by seven connectivity constants C 1 to C 7 which account for the average number of synaptic contacts.

Each impulse response

, and introduces a pair of first order ordinary differential equation. Consequently, the model can be represented by a set of 14 first order ordinary differential equations which can be reduced to an equivalent set of 10 equations [START_REF] Wendling | Epileptic fast activity can be explained by a model of impaired GABAergic dendritic inhibition[END_REF]. This set must be solved by numerical integration methods adapted to its stochastic nature. We used the Euler method with a time step equal to 1/256 ms (corresponding to the sampling period of real signals).

) t ( h EXC ) t ( h SDI ) t ( h FSI

Model output can be represented as

where M denotes the function operated by the model to produce Xs (the simulated EEG signal) from a realization n r of input noise n(t) and for given vector P of model parameters.

( )

P , n M X r s =

Identification of model parameters EXC, SDI and FSI using evolutionary algorithms

Problem statement

Model-based interpretation requires identification of model parameters which can be formalized as an optimization problem [START_REF] Hernandez | Model-based interpretation of cardiac beats by evolutionary algorithms: signal and model interaction[END_REF]. Since model parameters are based on physiological considerations, hypotheses for the interpretation of real EEG activity could be obtained from the optimal set of parameters P ˆ= [EXC, SDI, FSI], i.e. the one that minimizes an error function (

)

) P ( F , F ˆS o ε
where and are estimators of a feature vector F respectively computed on the observed EEG (X o F ˆ) P ( F ˆS O ) and synthesized EEG (X S ) using parameter vector P.

In the present study, F was primary based on spectral properties of the EEG signal. More precisely, the power within specific frequency bands B1 (0-4 Hz, delta), B2 (4-12 Hz, theta and alpha) and B3 (12-64 Hz, beta and gamma) was retained as three first features. As these spectral features did not distinguish well between interictal activity and interictal activity mixed with high amplitude spikes, we added a fourth feature based on a morphologic criterion indicating the presence of sporadic high extrema values. It is simply defined as the difference between the first and last α-order quantiles. The four features were computed on simulated and real signals normalized by their standard deviation (this provides invariance with respect to different amplitude scale factors in the model and in reality).

In the expression of the error function, we use an estimation instead of the exact value that can not be analytically computed as a function of P due to the nonlinear nature of the model. This estimation is computed over a finite duration ∆ and consequently is dependant on realization n

) P ( F ˆS ) P ( F S r of model input noise.
We are facing an optimization problem presenting multiple local optima in which the error function i) is not deterministic and ii) is not computable (and thus not differentiable). Consequently, deterministic optimization algorithms, such as those based on gradient descent, do not apply.

Two types of optimization methods are well-suited to this problem: combinatory and exhaustive search methods on the one hand and stochastic search methods on the other hand. As the former present prohibitive computational costs, we used a particularly interesting family of stochastic search methods, known as evolutionary algorithms (EA). EA have been the subject of intensive research during the last decade and have shown to be useful in the solution of hard identification problems, including different biomedical applications [START_REF] Hernandez | Model-based interpretation of cardiac beats by evolutionary algorithms: signal and model interaction[END_REF][START_REF] Pena-Reyes | Evolutionary computation in medicine: an overview[END_REF].

The optimization procedure using EA is detailed in appendix.

Variance of identified parameters and significance of parameter changes

In order to evaluate the dispersion of triplets of identified values, the whole identification procedure was applied thirty times to the same segment, leading to a set of thirty solutions for the optimal parameter set P ˆ= [EXC, SDI, FSI]. The distributions of parameter EXC, SDI and FSI values obtained for the thirty iterations were represented in a so-called "boxplot" form for statistical analysis and comparison. This compact graphical display provides a way to picture how parameter values are distributed relative to defined percentiles, to identify if a distribution is skewed and if outliers or unusual data values are present. They are used to graphically determine if significant statistical differences do exist between sets of parameter values. This graphical test consists in comparing boxes representing the lower and upper percentiles. If boxes do not overlap, parameter values are significantly different.

Results

Frequency content of segments of EEG activity chosen during intICTAL, preONSET, ONSET and ICTAL periods

For each patient, segments of EEG activity chosen during the interictal period ("intICTAL" segment), during the period preceding the onset rapid discharge ("preONSET" segment), during the rapid discharge itself that typically occurs at seizure onset ("ONSET" segment) and during the period that follows this rapid discharge as seizure develops ("ICTAL" segment) were analyzed for their frequency content. These segments are shown in figure 2. Normalized power spectral densities (PSDs) were computed using the standard periodogram method over 10 second duration segments in both the real and simulated case. PSDs were found to be very typical with respect to the analyzed period. The frequency content of intICTAL and preONSET segments was found to be relatively similar. Indeed, high amplitude sporadic spikes that appear during the preONSET period do not affect significantly the global spectral energy distribution. For the five patients, the first major change in PSDs is observed at seizure onset, during the rapid discharge, which is marked by the re-distribution of signal energy in a higher frequency band, typically from 15 to 40 Hz (EEG beta and low gamma bands). The second major change corresponds to the slowdown of the rapid discharge during the ictal activity, characterized by a transition from broadband to narrowband activity ranging from 3 to 10 Hz (mostly within the EEG theta band). Results are exemplified in figure 3 which also provides a time-frequency representation (spectrogram, figure 3-b) of the SEEG signal recorded during the transition from interictal to ictal activity in patient 2 (figure 3-a). In this case, the rapid discharge observed at seizure onset is relatively stationary with a maximum energy located in a frequency band ranging from 20 to 30 Hz. One can also notice the periodical aspect (around 5 Hz) of the signal recorded during ictal activity (figure 3-c, bottom right).

Identification of model parameters and simulation of EEG activity

For each patient, the model parameter identification procedure was firstly applied on each segment of EEG activity, in order to find appropriate parameter setting for the evolutionary algorithm. Several tests were performed to analyse the mean error value for each generation and the convergence rate of the EA, as a function of the main algorithm parameters. The best results were obtained with a mutation probability of p m = 0.2 and a crossover probability of p c = 0.9. The error function was computed for each individual by comparing spectral and morphologic features computed on simulated and real signal segments, as described in section Methods. The size of the population has been fixed to 200 individuals. The stopping criterion of the EA is met if no solution improvement is found during 10 consecutive generations, or if the maximum number of iterations (fixed to 200) is reached. The best individual of the last population was chosen as the identified solution, leading, in each case, to a triplet of values for the three model parameters: EXC (excitation), SDI (slow dendritic inhibition,), FSI (fast somatic inhibition). Some identification results found during this step are presented in figure 5 At this step, we also observed some common tendencies (among patients) in the evolution of model parameters for transitions of dynamics. In order to check whether these tendencies were significant or only due to the variability of identified parameters related to the randomness inherent to this identification procedure and to the randomness of the observed data, we conducted an analysis of the dispersion of estimated parameters EXC, SDI and FSI for each type of activity in each patient.

Dispersion of identified parameters

The identification procedure was repeated thirty times with random initialization on each segment of activity leading, in each case, to thirty triplets of values (EXC, SDI, FSI) for each patient (P1 to P5) and for each period (intICTAL, preONSET, ONSET, ICTAL). These results are displayed in figure 6, in the compact form of so-called "boxplots" (see section "material and methods", §"Variance of identified parameters"). Taken as a whole, the evolution of parameters (with respect to considered period of time) was found to be partly reproducible from patient to patient. One can also notice that the height of the boxes that represent the dispersion of the identified EXC, SDI and FSI parameters was small in most cases (except for parameter FSI during the preONSET period, see explanation below). This means that the identification procedure generally converged to similar triplets of values. Indeed, detailed examination showed the error (quantified using the error function) between the segments of activity simulated for these triplets and the real ones was low. The variance of the identification error was also found to be low over each analyzed period. This is shown by the values of mean error and its standard deviation computed for all segments of a given activity (intICTAL: ε=0.0092±0.0052; preONSET: ε=0.013±0.0095; ONSET: ε=0.013±0.0082; ICTAL: ε=0.027±0.016).

Note also that in the examples previously provided in figure 5, the error between synthesized and real activities corresponds to the above average errors. In each patient, the distribution of parameter FSI was found to be slightly broader when identified on the preONSET segment. During this period where parameter SDI increases, higher variability of the FSI parameter estimation might be explained by the fact that fast inhibitory processes are depressed by slow ones in the model and that model output becomes less sensitive to FSI parameter variations (see §"Computational model of hippocampus neuronal population" in material and methods). Finally, one can also notice that only few outliers (with respect to the thirty trials) are present in these boxplots. This means that the identification procedure was able to converge to a reproducible and small region in the parameter space for each period of interest and for each patient.

Evolution of model parameters with respect to transitions

The evolution of model parameters with respect to transitions is given in table 2.For the five patients, an increase of parameter EXC was observed during the transition from the intICTAL to the preONSET period and from the preONSET to the ONSET period. This parameter stayed either constant (P1) or increased (P2 to P5) during the following transition (ONSET to ICTAL). Regarding parameter SDI, evolution was found to be similar in the five patients. An increase was first always observed from intICTAL to preONSET. Then, the transition from the preONSET to the ONSET period was marked by a noticeable decrease of parameter SDI. Finally, this parameter was found to significantly re-increase during the transition from ONSET to ICTAL but without reaching again its initial value observed during the intICTAL period. Finally, regarding parameter FSI, no significant change was measured from intICTAL to preONSET and from preONSET to ONSET. Then, during the transition from ONSET to ICTAL period, we observed a noticeable decrease of parameter FSI.

Invariance of parameter evolution with respect to transitions

In order to determine whether reproducible changes do exist in the evolution of EXC, SDI, FSI parameters with respect to the type of transition, we compared results obtained in the five patients.

These are reported in figure 7-a which gives the evolution of each parameter, for each transition and for each patient (significant increase:'+', significant decrease: '-', no significant change: '='). We also averaged the identified parameter values of a given run over the 5 patients to obtain thirty mean values for each parameter and for each period. Figure 7-b displays the distribution of these mean values in the form of boxplots. Results show that significant changes (no overlap between boxes) do exist and confirm the scenario described hereafter. Regarding parameter EXC (average excitatory synaptic gain), a significant increase was observed between the intICTAL and preONSET period and between the preONSET and the ONSET period. Hence, when epileptic sporadic spikes as well as fast oscillations at seizure onset are observed, an increase of the average excitatory postsynaptic potential amplitude in the projections from pyramidal cells to interneurons as well as in the loop representing recurrent excitation is concomitantly observed in the model.

As far as model parameter SDI is concerned, three major changes were observed. First, parameter SDI increased from intICTAL to preONSET period. We noticed that this significant increase might be due to the aforementioned increase of excitation and could be related to a compensatory mechanism.

Second, parameter SDI also significantly decreased when identified over the ONSET period. Indeed, during fast oscillations activity at seizure onset, parameter SDI reached its lowest value in the five studied cases. Third, a re-increase of parameter SDI is observed during ICTAL activity. As discussed below, this result suggests that transient failure of GABA a,slow interneurons can occur during the transition to seizure.

Finally, concerning parameter FSI, two major observations were made. Over the ONSET period, in contrast to slow inhibition, FSI values were found either to be stable or to increase with respect to the preONSET period. During the transition from ONSET to ICTAL period, parameter FSI significantly decreased. In the model framework, these results suggest that fast oscillations are related to the fast feedback loop activity of gabaergic interneurons targeting pyramidal cells in their perisomatic region (GABAa,fast kinetics) and that the slowing down of the onset rapid discharge is related to an exhaustion of this fast inhibition.

Discussion

Different modeling strategies have been proposed to relate EEG signals with oscillations in neuronal populations. Besides modeling techniques lying at the cellular and network levels (Traub et al., 1999b), another approach lying at a higher level of organization, i.e. the population level, starts from the fact that neurons form populations. The resulting EEG is a reflection of ensemble dynamics rising from interconnections between a small number of subpopulations of pyramidal cells on the one hand, and interneurons on the other hand. This approach was initially proposed by Freeman and co-workers [START_REF] Freeman | Models of the dynamics of neural populations[END_REF] who made substantial progress in the understanding of perceptual processing in the olfactory system. Similar ideas developed in the same time by Lopes da Silva et al. [START_REF] Lopes Da Silva | Model of brain rhythmic activity. The alpha-rhythm of the thalamus[END_REF] led to the development of a lumped-parameter population model used to explain the alpha rhythm of the EEG and more recently to study the mechanisms underlying transitions between spike and wave (SW) discharges in non-convulsive epilepsy [START_REF] Suffczynski | Dynamics of non-convulsive epileptic phenomena modeled by a bistable neuronal network[END_REF].

On the basis of experimental results about cellular organization in the hippocampus, we followed a similar approach to build a physiological-relevant neuronal population model. Main features include an excitatory loop as well as two parallel inhibitory circuits that represent the two separate populations of GABAergic interneurons that contact anatomically segregated postsynaptic receptors (dendritic/somatic region, GABA A,slow / GABA A,fast ). As in any model, certain issues are simplified.

The major simplification of our model is related to its macroscopic level. The reason for this choice is twofold. First, the neuronal population level is adapted to the nature of our observations: simulated signals can be directly compared with real intracerebral EEG signals recorded with macroelectrodes in human temporal lobe epilepsies. Second, we wanted to determine to what extent ensemble interactions between main cells and interneurons within hippocampal neuronal populations can explain dynamics and transitions of dynamics observed in real EEG signals.

Results showed that the model can produce various types of signal dynamics that were qualitatively and quantitatively similar to those observed in real depth-EEG signals recorded in human hippocampus over interictal and ictal periods. In five patients, the identification procedure led to reproducible modifications about model parameter values of excitation, slow dendritic inhibition and fast somatic inhibition with respect of the analyzed type of EEG activity (interictal, preonset, onset rapid discharge and ictal). In the following, these results are discussed with respect to current hypotheses about possible mechanisms involved in epileptogenesis or ictogenesis.

Increased excitation from interictal to preonset and ictal activity

Causes of hyperexcitability of the epileptic tissue have been largely discussed during the past decades [START_REF] Jefferys | Models and mechanisms of experimental epilepsies[END_REF]. In the hippocampus, enhanced excitability may be caused by intrinsic and/or extrinsic factors. In the present study based one a single population, parameters changes may only be related to intrinsic mechanisms although we are aware that exogenous factors such as the enhancement of the direct entorhinal cortex input to CA1 can play a role in excitability changes in the hippocampus, as reported in [START_REF] Wu | Increased dendritic excitability in hippocampal ca1 in vivo in the kainic acid model of temporal lobe epilepsy: a study using current source density analysis[END_REF]. Results obtained from automatic identification of model parameters suggest an increase in the average excitatory synaptic gain in feedback loops that link pyramidal cells to local cells within the considered neuronal population during the transition from incterictal activity to preonset activity (interictal + sporadic spikes) and from preonset activity to ictal activity (beta and gamma oscillations at seizure onset and theta-like ictal activity). Although several studies performed in experimental models of focal epilepsy [START_REF] Bernard | Changes in neuronal excitability and synaptic function in a chronic model of temporal lobe epilepsy[END_REF], [START_REF] Gorter | Long-lasting increased excitability differs in dentate gyrus vs. CA1 in freely moving chronic epileptic rats after electrically induced status epilepticus[END_REF] or performed in slices from human brain [START_REF] Schwartzkroin | Cellular electrophysiology of human epilepsy[END_REF] have already demonstrated that epileptic discharge patterns are correlated with prolonged and/or enhanced EPSPs, very few studies have attempted to relate them to gradual changes observed at the EEG level during the transition to seizure, as done here through the model of the human hippocampus. Indeed, we infer from the model that a sustained excitation enhancement constitutes a necessary (but not sufficient) condition to explain the appearance of interictal sporadic spikes and that the excitation level must stay high to get fast oscillations and theta-like ictal activity on the EEG.

Increase and decrease of slow dendritic inhibition from interictal to seizure activity

Another striking result of this study is related to the evolution of slow inhibition mediated by GABAergic interneurons targeting the dendrites of pyramidal cells (GABA a,slow kinetics, as represented in the model). Indeed, model-based interpretation of slow dendritic inhibition evolution suggests that it significantly increases over the preonset period compared to interictal and abruptly decreases at seizure onset when faster oscillations are reflected on the EEG. At first sight, the reinforcement of inhibition may appear as paradoxical since increased inhibition could be supposed to prevent the occurrence of seizures. We interpreted this phenomenon as a compensatory mechanism that responds to increased excitation over the same preonset period. This model prediction is interesting since it relates to experimental results obtained with kainic acid (KA) known to generate seizures in the CA3 region [START_REF] Khalilov | Paradoxical anti-epileptic effects of a GluR5 agonist of kainate receptors[END_REF]: using in vitro intact hippocampus, authors demonstrated that kainate is responsible for ictal activity but at the same time augments the excitatory input to interneurons and thus lead to an enhancement of inhibition. With the model identification, we can go further and hypothesize that the increase of inhibition level is selective, i.e. it mainly relates to slow dendritic inhibition. Indeed, in the model, the two types of interneurons (GABA A slow and fast kinetics) receive positive input from pyramidal cells. However, since the subset of GABA A,fast interneurons is inhibited by the subset GABA A,slow interneurons [START_REF] Banks | Interactions between distinct GABA(A) circuits in hippocampus[END_REF], the increase of excitation level at both inputs will mainly act on this latter subset that represents interneurons projecting to the dendritic region of pyramidal cells. The second significant change regarding dendritic inhibition is a noticeable decrease at seizure onset when fast oscillations appear (ONSET period). This mechanism could be also related to the increase of excitation that takes place before seizure during the preonset period: it can be hypothesized that a point is reached where inhibitory interneurons cannot compensate anymore for the increased excitation and where provided inhibition collapses due to a "fatigue process". This hypothesis relates to the concept of "fragile inhibition" proposed by Wu and Leung [START_REF] Wu | Enhanced but fragile inhibition in the dentate gyrus in vivo in the kainic acid model of temporal lobe epilepsy: a study using current source density analysis[END_REF]) from a current source density analysis of the dentate gyrus (DG) in the kainic acid model of temporal lobe epilepsy. Authors showed that the relatively strong inhibition in the DG of KA-treated rats can breakdown readily in the presence of a low dose of GABA A receptor antagonist and concluded that the compensatory inhibition in the brain in animals and humans with temporal lobe seizures is vulnerable, although underlying mechanisms have not been clearly identified.

Role of fast somatic inhibition in the generation of fast oscillations at seizure onset

At seizure onset, faster oscillations (beta and low gamma bands) were observed on the EEG signal recorded from the hippocampus, as demonstrated by the computation of PSDs. Our identification results show that they are associated to a constant or increasing value of model parameter FSI, in contrast to the dramatic decrease of parameter SDI. These results suggest that observed fast oscillations can be mainly related to the activity of inhibitory interneurons targeting pyramidal cells in their perisomatic region (with GABAa,fast kinetics). In a previous study, we already demonstrated that a fast feedback inhibitory loop is necessary for the model to generate EEG signals with gamma oscillations [START_REF] Wendling | Epileptic fast activity can be explained by a model of impaired GABAergic dendritic inhibition[END_REF]. In the present study, results are obtained from direct automatic identification procedure applied on real intracerebral EEG signal and support our previous results. This hypothesis of EEG fast oscillations related to fast inhibitory interneurons is also supported by several other studies mainly dealing with normal neuronal tissue. Penttonen et al. [START_REF] Penttonen | Gamma frequency oscillation in the hippocampus of the rat: intracellular analysis in vivo[END_REF] suggested that rapid discharges reflect rhythmic post-synaptic potentials in hyperpolarized pyramidal cells brought about by rhythmically discharging somatic-projecting interneurons. Using computational models of interconnected neurons and interneurons as well as experimental data obtained from hippocampal slices from rats, White et al. [START_REF] White | Networks of interneurons with fast and slow gamma-aminobutyric acid type A (GABAA) kinetics provide substrate for mixed gammatheta rhythm[END_REF] also ended with the conclusion that GABA A,slow and GABA A,fast interneurons exist as two different populations and that they are useful for generating theta and gamma rhythms (CA1 region). More recently, similar findings were obtained from current source density analysis in the rat intact hippocampus during awake, attentive behavior [START_REF] Hajos | Spike timing of distinct types of GABAergic interneuron during hippocampal gamma oscillations in vitro[END_REF]. Authors conclude that their results support the hypothesis of a synaptic feedback model of gamma oscillations primarily involving pyramidal cells and perisomaticprojecting inhibitory interneurons.

Development of seizure activity

After the onset of seizure, a rhythmic activity (quasi-sinusoidal or spiking, generally in the theta band) is typically observed in recorded EEG signals. This activity corresponds to the progressive development of ictal activity that generally involves other mesial structures. In the present study dealing with a single neuronal population, hypotheses about interactions between populations distributed over distant brain structures are not taken into account. However, regarding ictal activity in the hippocampus, we ended in the five patients with a situation where i) the excitation level was high, ii) the slow inhibition level was equivalent or lower compared to interictal and iii) the fast inhibition level was very low compared to interictal activity. From this observation, two remarks can be made.

First, the transition from rapid discharge to rhythmic theta-like activity is explained, in the model, by the exhaustion of fast inhibitory processes. Second, seizure activity seems to correspond to a situation where inhibition level is globally reduced compared to excitation.

Conclusion

This computational modeling study offered the unique opportunity to relate electrophysiological patterns typically recorded with clinical electrodes during the transition from interictal to ictal activity in the human hippocampus to ictogenesis mechanisms and to generate hypotheses about these mechanisms. Results demonstrated that the transition from interictal to ictal activity can not be merely explained by an increase in excitation and a decrease in inhibition, but rather by a variety of complicated time-varying ensemble interactions between pyramidal cells and interneurons with slow and fast GABA A kinetics. In particular, the model predicts an increase of excitation during the preonset period (10 to 50 seconds before seizure onset), an abrupt drop of dendritic inhibition at seizure onset when a rapid discharge (fast oscillations in the gamma band) and a crucial role of perisomatically projecting interneurons in observed EEG gamma and theta ictal activity. To us, these predictions could be tested in future experiments. In slices, electrophysiological recordings combined with current source density analysis could be used to study the spatial origin of epileptic bursts (apical dendrites and soma-basal region). In epileptic patients, during long-term video-EEG recording, stimulation protocols based on a paired-pulse paradigm could be used to follow some physiological parameters such as excitability. At present, the level of detail in the proposed model (neuronal population that best corresponds to the nature of our observations, i.e. EEG signals) does not allow us to specify hypotheses about cellular and molecular mechanisms involved in the above described phenomena. However, recent works in the field of computational modeling suggest that dynamic scale changes could be performed in such situations in order to relate macroscopic phenomena to microscopic ones [START_REF] Vittorini | The OsMoSys approach to multi-formalism modeling of systems[END_REF]. Finally, as far as the clinical impact of this work is concerned, we think that the presented model-based approach is able to reveal physiological changes that take place inside recorded brain structures during the transition from interictal to seizure activity. Future works will be oriented towards the development of time-optimized identification procedures to process long duration signals (over a sliding window). Indeed, such procedures could be used to detect excitation and inhibition changes typical of those occurring before seizure onset, and therefore, might be complementary to other methods already developed in the field of seizure prediction (see [START_REF] Lehnertz | The First International Collaborative Workshop on Seizure Prediction: summary and data description[END_REF] for review).

Appendix: Evolutionary algorithms

Evolutionary algorithms (EA) are stochastic search techniques, inspired on the theories of evolution and natural selection, which can be employed to find an optimal configuration for a given system within specific constraints [START_REF] Holland | Adaptation in Natural and Artificial Systems[END_REF]. In these algorithms, each "individual" of a "population" is characterized by a set of parameters (or chromosome). An initial population is created, usually from a set of random chromosomes, and this population will "evolve", improving its global performance, by means of an iterative process. During this process, each individual is evaluated by means of a "fitness" function, and a new generation is produced by applying mutation and crossover operators on selected individuals that present high "fitness" values, with probabilities p m and p c respectively. Convergence and robustness properties of EA have been largely studied in the literature [START_REF] Beasley | An Overview of genetic algorithms: Part 2, Research Topics[END_REF][START_REF] Goldberg | Genetic Algorithms in search, optimization and machine learning[END_REF][START_REF] Michalewicz | Genetic Algotithms + Data Structures = Evolution Programs[END_REF]. These properties depend upon: i) adequate individual coding, ii)

proper definition of the fitness function and iii) selection of appropriate genetic operators for crossover and mutation.

i) Individual representation and initial population

Each individual represents an instance of the model and is characterized by the triplet P=[EXC, SDI, FSI] of the model. In order to reduce the search space, values for parameters EXC, SDI and FSI were bounded to the physiologically plausible intervals [1 10], [1 50] and [1 20], respectively. These intervals, which have been defined after an exhaustive exploration of the parameter space [START_REF] Wendling | Epileptic fast activity can be explained by a model of impaired GABAergic dendritic inhibition[END_REF], are employed by the EA during the construction of the initial population and the application of genetic operators. The initial population is constituted of randomly generated individuals. The three parameter values of a given individual are independently generated from a uniform distribution defined under the corresponding feasibility interval.

ii) Individual evaluation

As previously stated, model output not only depends on the set of parameter values P but also on the realization n r of the input random noise. For given parameter vector P, different realizations of this noise lead to different model outputs and thus to different realizations of (and consequently of ε). In our case, at least two different strategies can been defined for individual evaluation: i) one and only one realization n ) P ( F ˆS r is used for all individuals and for the whole identification process, and ii) a new realization is created for each generation. We retained this second strategy.

iii) Selection method

Once the error function has been evaluated for each individual, selection is carried out by means of the "Roulette Wheel" method, adapted for function minimization, in which the probability of selecting a given individual depends on the value of its error function, divided by the sum of all the error values of the population [START_REF] Beasley | An Overview of genetic algorithms: Part 2, Research Topics[END_REF]. Only standard genetic operators, defined for real-valued chromosomes, have been used in this work: "uniform crossover", which creates two new individuals (offspring) from two existing individuals (parents), by randomly copying each allele from one parent or the other, depending on a uniform random variable and "Gaussian mutation", which creates a new individual by randomly changing the value of one allele (selected randomly), based on a Gaussian distribution around the current value. This mutation operator respects the bounds defined for each allele by truncating the mutation value, if it is necessary. A detailed description of these genetic operators can be found elsewhere [START_REF] Michalewicz | Genetic Algotithms + Data Structures = Evolution Programs[END_REF]. Normalized power spectral densities computed on segments of EEG activity during the interictal period (intICTAL), before the beginning of seizure (preONSET), during the rapid discharge at seizure onset (ONSET) and finally after the rapid discharge (ICTAL) during seizure spread show that frequency contents of intICTAL and preONSET segments are quite similar. The ONSET period is marked by the re-distribution of signal energy in higher frequency band (15 to 30 Hz, in this case) whereas the ICTAL period is characterized by the appearance of a narrow band activity (around 5 Hz, EEG theta band).

Legends

  , which displays four examples of simulated signals identified from real EEG signals along with their corresponding PSDs and error values. It can be observed that identified parameters lead the model to generate EEG signals that are close (in terms of morphological and spectral content) to real EEG signals. One can notice that signals produced by the model are very realistic: waveforms are qualitatively similar to real ones and PSDs are quantitatively comparable (as exemplified by epileptic spikes in figure 5-preONSET period).

Figure 1 :

 1 Figure 1: Depth electrodes implanted for SEEG in temporal lobe epilepsy. (a) coronal view of a preoperative plan in which electrode #1 records electrophysiological activity within the hippocampal head (internal contacts) and the middle part of the middle temporal gyrus (external contacts) and in which electrode #2 records the activity within the entorhinal cortex (internal contacts) and the anterior part of the inferior temporal gyrus (external contacts). Generally, 5 to 8 electrodes are used to spatially sample the temporal region. The implantation accuracy is per-operatively controlled by telemetric Xray imaging. A post-operative computerized tomography (CT) scan without contrast is then used to verify the absence of bleeding and the precise location of each recording lead. Intracerebral electrodes were then removed and an MRI performed, permitting visualization of the trajectory of each electrode. Finally, CT-scan/MRI data fusion can be performed to anatomically locate each contact along the electrode trajectory. (b) Electrodes include 10 to 15 contacts, length: 2 mm, diameter: 0.8 mm, 1.5 mm apart.

Figure 2 :

 2 Figure 2: intracerebral EEG signals recorded from hippocampus during the transition from interictalactivity to ictal activity in five patients suffering from mesial temporal lobe epilepsy. Solid line rectangles correspond to segments of EEG signal (10 second duration) on which model parameters were identified. intICTAL: interictal activity, preONSET: mixing of interictal activity and highamplitude spikes that arise before seizure onset, ONSET: rapid discharge occurring at the electrical onset of the seizure and in which faster oscillations are observed, ICTAL: activity observed after the rapid discharge period as seizure develops. intICTAL segments were chosen 1 minute (at least) before seizure onset (note that high-amplitude epileptic spikes are not present in these segments).

Figure 3 :

 3 Figure 3: (a) intracerebral EEG signals recorded from hippocampus during the transition from interictal activity to ictal activity in patient 2 and (b) time-frequency distribution of its energy (spectrogram, red color indicates higher energy). This sequence of electrophysiological patterns is characteristic of the hippocampus activity observed in temporal lobe epilepsy of mesial origin. (c)

Figure 4 :

 4 Figure 4: Neuronal population model based on the cellular organization of the hippocampus. a) Schematic representation. A whole population of neurons is considered inside which a subset of principal cells (pyramidal cells) project to and receive feedback from other local cells. Input to interneurons is excitatory (AMPA receptor-mediated). Feedback to pyramidal cells is either excitatory (recurrent excitation) or inhibitory (dendritic-projecting interneurons with slow synaptic kinetics -GABA A,slow -and somatic-projecting interneurons -grey rectangle -with faster synaptic kinetics -GABA A,fast -). As described in (Banks et al. 2000), dendritic interneurons project to somatic ones. b) Corresponding block diagram representation. In each subset, the average pulse density of afferent action potentials is changed into an average inhibitory or excitatory post-synaptic membrane potential using a linear dynamic transfer function of impulse response , and while this potential is converted into an average pulse density of potentials fired by the neurons using a static nonlinear function (asymmetric sigmoid curve, S(v)). The subset of somatic-projecting interneurons (grey rectangle) receives input from both subsets of pyramidal and dendritic interneurons. One of the model outputs represents the summated average post-synaptic potentials on pyramidal cells. It reflects an EEG signal. The three main parameters of the model respectively correspond to the average excitatory synaptic gain (EXC), to the average slow inhibitory synaptic gain (SDI) and to the average fast inhibitory synaptic gain (FSI). These three parameters are automatically identified. Other

Figure 5 :

 5 Figure 5: Examples of simulated segments of EEG activity after automatic identification of model parameters from displayed segments of real EEG activity. Corresponding power spectral densities (solid line: real signal, dotted line: simulated signal). EXC: excitation, SDI: slow dendritic inhibition, FSI: fast somatic inhibition (see table1for details), ε: error between simulated and real segments (see section "Methods" for detailed description of error function).

Figure 6 :

 6 Figure 6: Parameter identification results for the 20 pre-segmented EEG signals (5 patients -P1 to P5during 4 different EEG activities). Boxplots are calculated from 30 independent realizations of the identification process for each segment of real activity. Results show the stability of the identified solutions and depict a global tendency in the evolution of parameter values according to the EEG activity.

Figure 7 :

 7 Figure 7: a) Comparison of results obtained in the five patients about the evolution of EXC, SDI, FSI parameters with respect to the type of transition (significant increase:'+', significant decrease: '-', no significant change: '='). Similar changes were found in the 5 studied cases (grey boxes). b) Identified parameter values of a given run of the identification procedure were averaged over the 5 patients.Distributions of averaged values in the form of boxplots show that significant changes (no overlap between boxes) do exist in the evolution of model parameters with respect to the type of activity (interictal, preonset, onset and ictal).

Table 1 :

 1 Model parameters, interpretation and values used to produce EEG signals with the model.Standard values (except g, C5, C6 and C7) were established in(Jansen & Rit, 1995). Main parameters (EXC, SDI and FSI) are automatically identified from real intracerebral EEG signals.

Table 2 :

 2 Average values and standard deviation of identified model parameters EXC, SDI, FSI for each one of the four considered periods of time (interictal, preonset, onset and ictal) and for each patient (P1 to P5).

  parameters are detailed in table 1.

											h EXC	(	t	)	h SDI	(	t	)	h FSI	(	t	)
	c) Time-course of average post-synaptic membrane potentials: excitatory, slow inhibitory and fast
	inhibitory respectively obtained from	h EXC	(t	)	,	h SDI	(	t	)	and	h FSI	(	t	)		for standard values of EXC (3.25
	mV), SDI (22 mV) and FSI (10 mV).