Interaction between synaptic inhibition and glial-potassium dynamics leads to diverse seizure transition modes in biophysical models of human focal seizures

Journal of Computational Neuroscience, Aug 2016

How focal seizures initiate and evolve in human neocortex remains a fundamental problem in neuroscience. Here, we use biophysical neuronal network models of neocortical patches to study how the interaction between inhibition and extracellular potassium ([K +] o ) dynamics may contribute to different types of focal seizures. Three main types of propagated focal seizures observed in recent intracortical microelectrode recordings in humans were modelled: seizures characterized by sustained (∼30−60 Hz) gamma local field potential (LFP) oscillations; seizures where the onset in the propagated site consisted of LFP spikes that later evolved into rhythmic (∼2−3 Hz) spike-wave complexes (SWCs); and seizures where a brief stage of low-amplitude fast-oscillation (∼10−20 Hz) LFPs preceded the SWC activity. Our findings are fourfold: (1) The interaction between elevated [K +] o (due to abnormal potassium buffering by glial cells) and the strength of synaptic inhibition plays a predominant role in shaping these three types of seizures. (2) Strengthening of inhibition leads to the onset of sustained narrowband gamma seizures. (3) Transition into SWC seizures is obtained either by the weakening of inhibitory synapses, or by a transient strengthening followed by an inhibitory breakdown (e.g. GABA depletion). This reduction or breakdown of inhibition among fast-spiking (FS) inhibitory interneurons increases their spiking activity and leads them eventually into depolarization block. Ictal spike-wave discharges in the model are then sustained solely by pyramidal neurons. (4) FS cell dynamics are also critical for seizures where the evolution into SWC activity is preceded by low-amplitude fast oscillations. Different levels of elevated [K +] o were important for transitions into and maintenance of sustained gamma oscillations and SWC discharges. Overall, our modelling study predicts that the interaction between inhibitory interneurons and [K +] o glial buffering under abnormal conditions may explain different types of ictal transitions and dynamics during propagated seizures in human focal epilepsy.

A PDF file should load here. If you do not see its contents the file may be temporarily unavailable at the journal website or you do not have a PDF plug-in installed and enabled in your browser.

Alternatively, you can download the file locally and open with any standalone PDF reader:

https://link.springer.com/content/pdf/10.1007%2Fs10827-016-0615-7.pdf

Interaction between synaptic inhibition and glial-potassium dynamics leads to diverse seizure transition modes in biophysical models of human focal seizures

J Comput Neurosci Interaction between synaptic inhibition and glial-potassium dynamics leads to diverse seizure transition modes in biophysical models of human focal seizures E. C. Y. Ho 0 1 2 3 Wilson Truccolo 0 1 2 3 Action Editor: Steven J. Schiff 0 1 2 3 0 Wilson Truccolo wilson 1 E. C. Y. Ho ernest 2 This research was supported by: the National Institute of Neurological Disorders and Stroke (NINDS), grant R01NS079533; the U.S. Department of Veterans Affairs, Merit Review Award I01RX000668; and the Pablo J. Salame '88 Goldman Sachs endowed Assistant Professorship of Computational Neuroscience. The contents do not represent the views of the U.S. Department of Veterans Affairs or the United States Government 3 U.S. Department of Veterans Affairs, Center for Neurorestoration and Neurotechnology , Providence, RI , USA How focal seizures initiate and evolve in human neocortex remains a fundamental problem in neuroscience. Here, we use biophysical neuronal network models of neocortical patches to study how the interaction between inhibition and extracellular potassium ([K +]o) dynamics may contribute to different types of focal seizures. Three main types of propagated focal seizures observed in recent intracortical microelectrode recordings in humans were modelled: seizures characterized by sustained (∼ 30 − 60 Hz) gamma local field potential (LFP) oscillations; seizures where the onset in the propagated site consisted of LFP spikes that later evolved into rhythmic (∼ 2 − 3 Hz) spikewave complexes (SWCs); and seizures where a brief stage of low-amplitude fast-oscillation (∼ 10 − 20 Hz) LFPs Focal epilepsy; Seizure dynamics; Spike-wave discharges; Gamma oscillations - Department of Neuroscience & Institute for Brain Science, Brown University, Providence, RI, USA preceded the SWC activity. Our findings are fourfold: (1) The interaction between elevated [K +]o (due to abnormal potassium buffering by glial cells) and the strength of synaptic inhibition plays a predominant role in shaping these three types of seizures. (2) Strengthening of inhibition leads to the onset of sustained narrowband gamma seizures. (3) Transition into SWC seizures is obtained either by the weakening of inhibitory synapses, or by a transient strengthening followed by an inhibitory breakdown (e.g. GABA depletion). This reduction or breakdown of inhibition among fastspiking (FS) inhibitory interneurons increases their spiking activity and leads them eventually into depolarization block. Ictal spike-wave discharges in the model are then sustained solely by pyramidal neurons. (4) FS cell dynamics are also critical for seizures where the evolution into SWC activity is preceded by low-amplitude fast oscillations. Different levels of elevated [K +]o were important for transitions into and maintenance of sustained gamma oscillations and SWC discharges. Overall, our modelling study predicts that the interaction between inhibitory interneurons and [K +]o glial buffering under abnormal conditions may explain different types of ictal transitions and dynamics during propagated seizures in human focal epilepsy. 1 Introduction Epilepsy is a neurological disorder characterized by recurrent seizures affecting an estimated 65 million people worldwide (Thurman et al. 2011) . Its apparent simple characterization belies the often complex and varied cellular and synaptic basis underlying the abnormal neuronal activity during epileptic seizures. As a result, despite significant advances in understanding epilepsy over the past several decades, key knowledge gaps remain with respect to the disorder. One example is the mechanism(s) concerning the transition of neural dynamics from normal to ictal as focal seizures propagate to distal regions in the brain, especially in neocortex. Seizures have been hypothesized to originate from a wide variety of cellular and network mechanisms (Richardson et al. 2008) . It is possible that, even in the same patient, different seizures might be triggered by different factors. Facing the diversity of possible triggers of seizure transition, recent theoretical studies have proposed canonical models that focus on general unifying dynamical principles (Kramer et al. 2012; Jirsa et al. 2014; Sritharan and Sarma 2014; Wang et al. 2014; Wei et al. 2014b; Naze et al. 2015) . Recent animal experiments (Trevelyan et al. 2006; Zhang et al. 2012; Cammarota et al. 2013; Zˇ iburkus et al. 2013) and human studies (Truccolo et al. 2011; Schevon et al. 2012; Ahmed et al. 2014) have also suggested that complex dynamics between different types of neurons can be in play as cortical networks transition into seizures. In particular, these studies have indicated the potential critical role of fast-spiking (FS) inhibitory interneurons (see Paz and Huguenard 2015 for a review) and their interactions with changes in extracellular potassium concentration [K+]o (e.g. Zˇ iburkus et al. 2006) . Elucidating such interplay between interneuron activity and ionic concentrations during seizure initiation and propagation can provide the much-needed clues for the development of seizure prediction, early detection and novel therapeutic approaches. Furthermore, studies by Truccolo et al. (2011, 2014) and Wagner et al. (2015) have characterized spiking in ensembles of single neurons and local field potential (LFP) activity recorded via 96-microelectrode arrays in humans during focal seizures. Three main types of seizure activity were recorded in neocortical patches distal (∼ 2 − 4 cm) to the identified seizure onset zones: seizures consisting of sustained (∼ 30 − 60 Hz) gamma-band oscillations, seizures where the onset consisted of a direct transition into LFP spike discharges which evolved later into (∼ 2−3 Hz) spikewave complexes (SWCs), and seizures where the onset included a brief transient of low voltage fast oscillations (10−20 Hz) before the evolution into SWC discharges. (For a more comprehensive classification of seizures in terms of onset patterns and relationships to pathology in larger datasets, see Perucca et al. 2014) . In the gamma seizures, despite the narrowband nature of the oscillations, neuronal spiking activity remained highly irregular and asynchronous (Truccolo et al. 2014) . Classification of recorded single neurons into putative principal and interneuron cells suggested that inhibition was preserved throughout the various ictal stages in these gamma seizures (Truccolo et al. 2011) . On the other hand, a recent investigation of SWC seizures in these datasets indicated that inhibitory activity seems to shut off just before the emergence of large amplitude spikewave discharges (Ahmed et al. 2014) . Ahmed et al. (2014) hypothesized the role of depolarization block in interneurons and examined their sensitivity to [K+]o in a single cell model of a FS inhibitory interneuron. Here, we examine how the interaction between FS inhibitory interneuron activity and extracellular potassium concentration in biophysical neuronal network models can capture the different types of propagated seizure dynamics observed in these human neocortical recordings. We used neuronal network models of neocortical patches that include pyramidal and FS interneurons (Wang and Buzsa´ki 1996) , and [K+]o diffusion coupled with glial cell activity. The neuronal network models spanned 9 × 9 minicolumns, each 25 μm wide and containing 12 pyramidal neurons and 4 FS interneurons (for a total of 9 × 9 × 16 = 1296 model neurons). Coupled [K+]o and glial buffering dynamics were adapted from previous models in Fro¨hlich et al. (2010) . Several previous computational models (Kager at al. 2000; Bazhenov et al. 2004; Traub et al. 2005; Anderson et al. 2007; Cressman Jr. et al. 2009; Ullah et al. 2009; Fro¨hlich et al. 2010; Krishnan and Bazhenov, 2011; Krishnan et al. 2013, Wei et al. 2014a, b) have explored the role of intracellular and extracellular concentrations of various types of ions (e.g. potassium, sodium, chloride, calcium, etc.), as well as inhibitory/excitatory conductances and metabolic factors, on seizure activity. Our model specification was inspired by several of these previous studies, including, for example, the formulation of glial-buffering processes (see Section 2). Nevertheless, we emphasize that, to our knowledge, the work presented here is the first to examine how the interaction between inhibitory activity and abnormal levels of extracellular potassium concentration can lead to three different types of seizures (gamma, SWC, low-voltage fast-oscillations followed by SWCs) observed in intracortical recordings from patients with focal epilepsy. Overall, we demonstrate that transition into and maintenance of both gamma-band and SWC seizures in the model can be explained by variations in inhibition strength (and therefore changes in the excitation/inhibition balance) and [K+]o levels. We consider this to be the main novel contribution reported in our study. 2 Methods 2.1 Model description and parameter choices Single neuron model We used a single-compartment conductance based model to represent both the fast-spiking inhibitory interneurons and pyramidal cells (Wang and Buzsa´ki 1996; Skinner 2006; Anderson et al. 2007) . The differential equation for the membrane potential of a given neuron corresponds to: where X can be h or n, with αh(V ) = 0.07 × exp − C × dV l dt = Iilnt (t ) + Islyn(t ) + Ielxt (t ), (1) for l = 1, 2, . . . , 1296 neurons. Iilnt denotes the intrinl sic current of the neuron in question, Isyn denotes synaptic currents and Ielxt denotes external input. The specific capacitance C is assumed to be 1 μF/cm2. Iilnt has the following general form: Iilnt = −INa − IK − IL − IKCa = −gNa × m3∞ × h × (V − ENa) −gK × n4 × (V − EK ([K+]o)) −gL × (V − EL([K+]o)) − IKCa, (2) which represent sodium, delayed-rectifier potassium, leak and calcium-mediated potassium currents, respectively. The leak reversal potential is dependent on various ionic concentrations via the Goldman-Hodgkin-Katz voltage equation: EL([K+]o) = 26.64mV 1 × [K+]o + 0.085×[N a+]o + 0.1× [Cl−]i , × loge 1 × [K+]i + 0.085×[N a+]i + 0.1× [Cl−]o (3) where [K+]o, [K+]i , [N a+]o, [N a+]i , [Cl−]o, [Cl−]i denote ionic concentrations, with subscript i and o representing intracellular and extracellular, respectively. The coefficient before each type of ionic concentration is the relative membrane permeability of that ion. The “fast” activation variable m∞ of the sodium current is assumed to have reached a steady state and has the form: αm m∞ = αm + βm , −0.1 mV−1 × (V + 35 mV) αm(V ) = exp(−0.1 mV−1 × (V + 35 mV)) − 1 V + 60 mV βm(V ) = 4 × exp − 18 mV , , and the inactivation variable h of the sodium and the activation variable n of delayed-rectifier potassium current obey the following differential equation: dX dt = 5 ms−1 × (αX × (1 − X) − βX × X), (4) (5) βh(V ) = exp(−0.1 mV−1 × (V + 28 mV)) + 1 , 0.01 mV−1 × (V + 34 mV) αn(V ) = − exp(−0.1 mV−1 × (V + 34 mV)) − 1 , βn(V ) = 0.125 × exp − V 8+04m4VmV . (6) The purpose of the calcium-mediated potassium current (IKCa) is to provide the spike-frequency adaptation of the model pyramidal cells. It has the following form: [Ca2+]i IKCa = gKCa × (V − EK ) × 1mM + [Ca2+]i , where gKCa and [Ca2+]i denote the conductance for calcium-mediated potassium currents and the intracellular calcium concentration, respectively. The [Ca2+]i dynamics are given according to Cressman et al. (2009) and Ullah et al. (2009) : (7) d[Ca2+]i dt = − 0.002 cmm·omlC × gCa × (V − ECa) 1 + exp − (V2+.525mmVV) − [C80a2m+s]i , (8) where gCa and ECa denote the conductance and reversal potential related to Ca2+, respectively. Network structure and connectivity The structure of our network of model neurons is based on the arrangement of neurons in minicolumns in the neocortical area (Buxhoeveden and Casanova 2002) . In this network model, each minicolumn consists of 16 single-compartment model neurons, of which 12 are excitatory pyramidal cells and 4 are inhibitory interneurons. The neurons are arranged into four layers within each minicolumn (Fig. 1). Each layer consists of four neurons in a square arrangement, with each vertex of the square being the location of one neuron. Pyramidal cells occupy the top three layer and the bottom layer is for the four inhibitory interneurons. (We do not attempt here to replicate the actual 6-layer neocortical structure). The entire cortical network consists of 9×9 such minicolumns for a total of 1296 model neurons. Details of the structure configuration of the cortical network model and various distance parameters between model neurons are given in Fig. 1 and Table 1. Unless otherwise stated, we chose a global random connectivity to avoid boundary effects associated with the small size of the simulated cortical network. Glial and potassium ionic dynamics Each neuron in the network is endowed with a potassium pump. An associated (9) Values where K + and B denote the extracellular potassium ions and the free buffer in the glia, respectively, and kb and kf correspond to the backward and forward glial binding rate, respectively (Tables 3 and 4). K + and B combine reversibly to form the bound form of the buffer K B. We fixed the maximal amount of free glial buffer [B]max to 500mM (Table 3, see also Kager et al. 2000) . The above process (10) can be written as a set of differential equations: ∂ [K +]o ∂ t = kb × ([B]max − [B]) 1.1 − kf × [K +]o × [B], |glia d [B] d t kf = = kb × ([B]max − [B]) − kf × [K +]o × [B], kb/1mM 1 + exp [K+]o−[K+]o(th) θ , where [K +]o(th) and θ are specified in Table 3. We vary [K +]o(th)and θ to simulate a range of physiological and pathological states for the potassium glial-buffering system. [K +]o(th) corresponds to the equilibrium level (set point) for the extracellular potassium concentration being targeted by the glial buffering system. How tightly the glial buffering system tracks this equilibrium set-point depends on the parameter θ . For example, the lower the absolute values of θ , the closer the glial buffer system attempts to keep the [K +]o concentration to the target [K +]o(th) level. Under this condition, the system may fail if there is a big enough difference between [K +]o and [K +]o(th), leading to bi-stability and abnormal functioning. In more physiological conditions, the absolute value of θ is kept at higher values so that the glial buffering system tracks [K +]o(th) in a softer fashion, allowing a broader range of [K +]o levels and more effectively avoiding bi-stability. We note that while we vary [K +]o(th) and θ in order to change the glial buffering dynamics and to induce abnormal states, we do not yet have a hypothesis about specific causes that would lead to these changes in first place. The diffusion of extracellular potassium ions across the network was modelled according to the diffusion equation: ∂ [K +]o ∂ t diff = D × ∇2[K +]o, where ∇2 is the Laplace operator ∂∂x22 + ∂∂y22 + ∂∂z22 and D is the diffusion constant. In the simulations, ∇2 is approximated by a set of finite difference equations, since the model potassium uptake dynamics related to nearby glial cells are included in the model to regulate the extracellular potassium concentration [K +]o. The potassium pump is represented by: IKpump = IKmax 1 + [K[K+]+o](eoq) 2 , where IKmax and [K +]o(eq) are defined in Table 3. The glial potassium uptake mechanism is modelled in terms of a reversible binding process between extracellular potassium ions and the glial buffer (Kager et al. 2000; (11) (12) neurons are discretely located (Fig. 1). We impose the “noflux” or Neumann boundary condition (i.e. ∂[K+]o = 0 on ∂x the y − z planes at both the x boundaries and so on) in the network simulations. Therefore, the total time derivative of [K+]o of a particular model neuron, d[K+]o (unit in mM·ms−1), is the sum of dt all the potassium, pump, glia and diffusion currents d[K+]o dt 50 mol = 96489 cm · mC × (IKCa + IK − IKpump +gL|[K+] × (V − EK ([K+]o))) + diff , where gL|[K+] (≥ 0) is the value of the leak conductance attributable to potassium. The value can be estimated by fitting the equation, gL × (V − EL([K+]o)) = gL|[K+] × (V − EK ([K+]o)) +gL|[Cl−] × (V − ECl ) +gL|[Na+] × (V − ENa ), with the constraint that gL = gL|[K+] + gL|[Cl−] + gL|[Na+], where gL|[Cl−] (≥ 0) and gL|[Na+] (≥ 0) are the values of leak conductance due to Cl− and N a+ respectively, while ECl and ENa are their respective reversal potentials as shown in equations (16) and (17). gL|[K+] values at ∼ 6– 13 % of total leak give reasonable fits over a wide range of physiological [K+]o and membrane potential values. (Various fitting scenarios indicate that >70 % of the total leak conductance value is due to chloride). Based on the above procedure, we assume a gL|[K+] value at 6 % of the gL value. The potassium reversal potential EK is affected by [K+]o via the Nernst equation: ∼ −74mV, ∼ 54mV. (13) (14) (15) (16) (17) Synaptic model Unless otherwise stated, the synaptic current Islyn corresponds to l Isyn = Γm→l × gsmy→nl × sm(t −τdelay ) × (E{e,i} −V l ), m=l with s·(t ) representing the synaptic gating variable of each model neuron, such that (18) τ{e,i} sl → sl + smax (when V l crosses 0 from below), ddstl = − sl (otherwise). In the above, Γm→l denotes the binary connectivity matrix element in which a value of 1 represents the existence of a synaptic connection in neuron l from neuron m, and a value 0 represents the absence of such connection. Table 5 lists the probabilities of connection between pyramidal neurons (P (e → e)), from pyramidal neurons to interneurons (P (e → i)), from interneurons to pyramidal neurons (P (i → e)) and between interneurons (P (i → i)). The term E{e,i} denotes the excitatory or inhibitory reversal potential value, respectively. (Whether the reversal potential is excitatory or inhibitory depends on the particular neuron–i.e. neuron m–from which neuron l receives synaptic current). As done in Ho et al. (2012) for computational efficiency, we use a discontinuous model to represent the opening of the synapses (equation set (18)) whenever there is a spike. Parameters τ{e,i} and smax (for sl ) denote the excitatory (inhibitory) synaptic decay time constant and the value of the maximal opening of the synaptic gates per spike, respectively (whether the decay time constant should be excitatory or inhibitory depends on whether neuron l itself is a pyramidal cell or an inhibitory interneuron). Table 5 provides more details on the terms of the above equations. Background synaptic activity Unless otherwise stated, each model neuron is driven by an external input Iext with the form (Rudolph et al. 2004) to mimic the activities of the fluctuating synaptic background: Ielxt = gel(t ) × (Ee − V l ) + gil (t ) × (Ei − V l ), dg{el,i}(t ) dt = g{e,i}0 −g{el,i}(t ) τ{e,i} + 2 × σ{2e,i} τ{e,i} ×χ{el,i}(t ), (19) where E{e,i} denotes the excitatory/inhibitory reversal potential and τ{e,i} is the excitatory/inhibitory time constant. The excitatory/inhibitory conductance values g{e,i} are stochastic variables following an Ornstein-Uhlenbeck process (equation set (19)), with mean g{e,i}0 and SD σ{e,i} (excitatory/inhibitory fluctuations). The stochastic term χ{e,i} is such that 0 χ{e,i}(s) × ds is a Gaussian distributed random variable with mean zero and variance > 0. EK = 26.64mV × loge ECl = 26.64mV × loge and ENa by ENa = 26.64mV × loge [K+]o [K+]i , [Cl−]i [Cl−]o [N a+]o [N a+]i where the intracellular potassium concentration [K+]i is set at 133 mM (Table 3). Similarly, ECl is given by See also Table 3 for values of [N a+]o, [N a+]i , [Cl−]o and [Cl−]i . Further details on the numerical implementation of the stochastic elements of equation set (19) can be found in Ho et al. (2012) . Table 5 provides more details on the terms in the above equation. 2.2 Implementing the simulations The code for the model cortical network simulation was written in C++ (gcc, version 4.4.7, http://gcc.gnu.org). We used open MPI (version 1.8, http://www.open-mpi.org) for parallel processing of each simulation. Each typical 1296-neuron simulation used eight CPUs. A 10-hour simulation on the Brown University computer cluster (https://www.ccv.brown.edu/) usually yielded ∼ 2-3 minutes of data. Parameter scanning was automated via the use of custom Perl (version 5.18.4) and Linux shell scripts. Post simulation analysis of data was carried out using Matlab (version R2013b, http://www. mathworks.com) , Octave (version 3.8.1, http://www.gnu. org/software/octave/) and customized Perl scripts. Single neuron simulations (Fig. 2) were performed using xppaut (Ermentrout 2002) . Graphical rendering of data was performed by gnuplot (version 4.6.3, http://www. gnuplot.info), xfig (version 3.2) and gimp (version 2.8.10, http://www.gimp.org). The parula palette of the timefrequency plots is from https://github.com/Gnuplotting/ gnuplot-palettes. The C++ code is available in the NEURON ModelDB database (https://senselab.med.yale.edu/ ModelDB/showModel.cshtml?model=190306). The ModelDB accession number for the network model reported in this paper is 190306. 3 Results We build model neocortical networks (Fig. 1) using conductance-based model neurons endowed with glial dynamics for potassium uptake. Based on the observations of human focal seizure data (LFPs and spiking activity in ensembles of single neurons), our aim is to use these model networks to reproduce the neural dynamics during the transition into ictal states of propagated focal seizures. Our presentation is structured as follows. First, based on both single neuron and network simulations (Section 2), we demonstrate that pathological glial uptake of potassium can lead to bi-stable dynamics in the model cortical network. The bi-stable dynamics consist of a “high-activity” state characterized by an elevated [K +]o level, and a state of “low-activity” with a low [K +]o level. We relate the “highactivity” state to “ictal” events and the “low-activity” state to non-epileptic, normal activity. Although network bistability has been implicated in computational studies as an avenue for seizure transition (Fro¨ hlich et al. 2010) , here we demonstrate that abnormal glial buffering can increase the propensity of single neurons (and thus the cortical network) to exhibit bi-stability (Fig. 2). The underlying mechanism for the bistabilty is the imbalance in the reversible reaction of the potassium buffer between the bound and unbound forms (10). Second, with the seizure transition mechanism in place for the model cortical network, we then examine how certain types of synaptic changes involving mostly inhibitory interneurons and their interaction with potassium extracellular concentration may lead to various types of seizure transitions as observed in the human data. Our strategy is to explore inhibitory synaptic and intrinsic parameters relevant to the occurrence of depolarization block in model inhibitory interneurons. We also note that, in our model, pyramidal neurons are not as prone to depolarization block (in comparison to inhibitory interneurons) because of their Ca2+-mediated potassium conductances (gKCa ) and the corresponding after-spike hyperpolarization effects (7). Variations of parameters related to synaptic inhibition and glial potassium buffering allow us to replicate in detail three types of seizure patterns observed in the human data (gamma seizures and two types of spike-and-wave seizures; Truccolo et al. 2014) . Thus, all three types of seizures can each be linked, via the same computational network model, to specific temporal patterns of synaptic activities and glial buffering dynamics during the evolution of propagated ictal states (Figs. 3, 4 and 5). Table 2 summarizes the results of our network simulations in terms of the paths leading to the three types of seizure transitions. All of the primary variables and parameters used in the simulations and analyses below are introduced and defined in Tables 1, 3, 4 and 5. 3.1 Firing properties of the single-neuron models under physiological and abnormal glia uptake conditions We used four single pyramidal cell model simulations to show the effects of the potassium pump and glial dynamics on neuronal excitability (Fig. 2). In each of these four simulations, Iext = 0.36, [K +]o(eq) = 3, [K +]o(th) = 15, IKmax = 7, gKCa = 5. In this section, we define the “physiological” value θ = −1.15, and the “pathological” value as θ = −0.05 (equation set (11)). Table 3 specifies the units of these parameters. Our aim is to show that the neural system tends to settle into a single steady state for a wide range of initial conditions when the increase in the forward binding rate kf is smooth with increasing [K +]o (as in the physiological case, effected by a larger absolute value of θ – see equation set (11)). However, when there is a sharper increase in kf with increasing [K +]o (as in the pathological case, effected by a decreased absolute value of θ ), some of the initial conditions can lead the neural system into a state where virtually all of the glial buffers are free (i.e. [K B] → 0). In this case, the system’s ability to return potassium to the extracellular space and maintain a target value of extracellular potassium level is impaired. The impairment occurs because Eq. (10) is no longer reversible when [K B] → 0. In Fig. 2a and b, we show two single pyramidal neuron model simulations with the same set of “physiological” glial and potassium pump parameters but with different [K +]o and [B] initial conditions. In Fig. 2a, we set [K +]o(t = 0) = 1 mM and [KB](t = 0) = 3 mM, and in Fig. 2b, [K +]o(t = 0) = 10 mM and [KB](t = 0) = 10 mM. It is clear that, despite the difference in initial conditions between the two simulations, the model neuron and the glial buffering system in both simulations converge to the same steady state. The steady state is characterized by a dynamic equilibrium established between the glial buffer’s binding (kf ) and unbinding (kb) of the extracellular potassium. This dynamic equilibrium is evidenced by the small but non-vanishing [K B] in the steady state. Moreover, all the steady state variables (including firing rates, [K +]o and [K B]) are identical between the two simulations, indicating that the steady states of the two simulations are the same. The glial buffering system enforces an upper ceiling for the extracellular potassium level in both simulations. In Fig. 2c and d, we show the same two simulations with “pathological” glial parameters. In contrast to Fig. 2a and b, we see that the same set of different initial conditions results in distinct steady states for the two simulations. In Fig. 2c, the model neuron and the glial buffer approach a state where there is a lower [K +]o and a vanishing value of [K B], while in Fig. 2d, the final [K +]o and [K B] values are higher and the neuron is excitable with a higher firing rate. Thus, bi- or multi-stable dynamics at the single neuron level are more easily obtained with a smaller absolute value of θ , resulting into multiple fixed points (steady states) having either a lower [K +]o or higher [K +]o value. The main factor behind this bi-stable dynamics is the parameter θ , which determines the level of bias for the forward binding rate kf against the backward rate kb (see equation sets (10) and (11)). With a small absolute value of θ , at a critical [K +]o level lower than [K +]o(th), kb can be overwhelmingly larger than kf . The large kb value rapidly lowers [K B]. Although this backward unbinding action releases K + back to the extracellular space, the IKpump (9) counteracts this effect and therefore [K +]o remains low. Moreover, the small kf value (due to the bias factor θ ) is not sufficient to replenish [K B] via the forward binding action. As a result, the system ends up in a state where virtually all the buffer is in the free state B with a low level of [K +]o (as in Fig. 2c). For the neuron and glial system to exit this state, as we have shown in Fig. 2d, one can either artificially increase the values of [K +]o or [K B] (for example, by setting up the initial conditions). Increasing [K +]o to a value closer to [K +]o(th) forces kb to be less biased against kf (equation set (11)), so that [K B] can be more efficiently replenished and thus allows a dynamic equilibrium between [B] and [K B]. Increasing [K B] artificially serves similar purposes by increasing the backward unbinding release of K + to the extracellular space (thus increasing [K +]o). In essence, decreasing the absolute value of θ reduces the range of [K +]o over which the glial buffering system can efficiently regulate, and thus makes the single neuron more susceptible to multi-stable behaviour. We show in the next section that this apparent multistability also carries over to the network simulations with abnormal glial parameters. As mentioned above, we will associate the state in which the [K +]o is low with the normal physiological state, and the elevated [K +]o level state with the “seizure” state. We also note that, although the above single-neuron analysis was based on pyramidal neurons, the same low and elevated [K +]o states occur in the case of FS inhibitory neurons. The effect of elevated [K +]o on pyramidal and FS interneuron spiking will, however, differ in the network simulations presented below. 3.2 Sustained gamma epileptiform activity results from an abnormal but balanced glial dynamics and inhibitory synaptic strength In this section we show how most features of gamma seizures observed in the human focal epilepsy data (Truccolo et al. 2011, 2014) can be replicated by our network model. These gamma seizures in the human data are characterized by sustained narrow band gamma LFP Seizure type Synaptic inhibition Gamma Spike-wave complex (SWC) Low-voltage fast oscillations followed by SWCs Strong synaptic inhibition Weak synaptic inhibition Strong synaptic inhibition at seizure onset, followed by a weakening of inhibition (e.g. GABA depletion). Interneuron depolarization block No Yes Yes Representative simulation tration due to the glial action ∂[K+]o |glia. Bottom panel shows the ∂t potassium pump current IKpump. (a)“Physiological” parameters. Initial conditions [K+]o(t = 0) = 1mM, [KB](t = 0) = 3mM. (b)“Physiological” parameters. Initial conditions [K+]o(t = 0) = 10mM, [KB](t = 0) = 10mM. (c)“Pathological” parameters. Initial conditions [K+]o(t = 0) = 1mM, [KB](t = 0) = 3mM. (d)“Pathological” parameters. Initial conditions [K+]o(t = 0) = 10mM, [KB](t = 0) = 10mM oscillations (∼ 30 − 60 Hz). However, few fine temporal synchrony transients exist in the spiking activities of the measured neurons. Neuronal spiking tends to be asynchronous and heterogeneous (Truccolo et al. 2011, 2014) . Individual neurons also fire at a lower rate than the gamma frequency of observed LFP oscillations. Fig. 3 An example gamma seizure simulation with pathological glial and pump parameters. (a) Depiction of neural activity and time course of various biophysical parameters during a 90-second simulation. The network starts with a “low-activity” non-seizure state (first 40 seconds of simulation). A DC stimulation is applied to every neuron in the network between 40 and 42.5 seconds, after which the system is “kick-started” into a high activity state with an elevated level of extracellular potassium (after 42.5 seconds). Topmost panel: average spike rate of pyramidal cells (red) and interneurons (blue) during the simulation. Second panel: raster plot of action potentials of neurons in 2 of the 81 simulated minicolumns. These 2 minicolumns are located in the centre of the network. Third panel: average membrane potential values of interneurons (blue) and pyramidal cells (red) during simulation. Fourth panel: average extracellular potassium level values for interneurons (blue) and pyramidal cells (red) during simulation. Fifth panel: Time-frequency plot of the average membrane potential values of the pyramidal cells (third panel-red). Clear sustained gamma band oscillations (∼40 Hz) at the population level emerge with the DC stimulation and persist even after the stimulation is terminated. (b) A 2-second segment from the 90-second simulation in a. Same conventions as in a, except that the bottom panel is a power spectral density (PSD) plot of the 2-second average membrane potential values of the pyramidal cells. The maximum of the PSD is located in the gamma range. (c) Histograms of interspike intervals (ISIs) of all the pyramidal cells (red) and interneurons (blue) during gamma population activity (t 45 seconds). Insets: Histograms of the ISI coefficients of variation (CV) for spike trains of individual pyramidal cells (red) and interneurons (blue), also during gamma population activity (t 45 seconds) Similar neural dynamics is observed in our model simulations (Fig. 3). Critical ingredients for the appearance of sustained gamma band LFP oscillations with asynchronous neuronal firing at a lower rate include high inhibitory conductance values (for both gsiy→ni and gsiy→ne) and an axonal conduction time delay (Brunel 2000; Brunel and Wang 2003) . We also require a value of [K +]o(th) that is not too high (see equation set (11)). A moderate [K +]o(th) value ensures that the model FS inhibitory interneurons have a high excitability in the “epileptic” state, while preventing the interneurons from being overly excited to enter depolarization block. We used gsey→n e = 0.0007, gsey→n i = 0.0007, gsiy→ni = 0.025, gsiy→ne = 0.025, ge0 = 0.01026, σe = 0.0025, gi0 = 0.084, σi = 0.02 for this gamma seizure simulation (Table 5 introduces and defines these synaptic parameters). As for intrinsic and glial properties of model neurons, we chose θ = −0.15, [K +]o(eq) = 3.6, [K +]o(th) = 7.5, IKmax = 1.45, gKCa = 5 for pyramidal cells and θ = −0.15, [K +]o(eq) = 3, [K +]o(th) = 7.5, IKmax = 1.9, gK = 6.8 for interneurons. The value of gK for interneurons was chosen so that they only enter depolarization block at a relatively high level of [K +]o ( 8 mM). Figure 3a shows the first 90 seconds of the temporal progression of several biophysical variables in the model and pyramidal cells (red). Fifth panel: Time-frequency plot of the average membrane potential values for the pyramidal neurons (third panel-red). The power peak at ∼5 Hz (after the DC stimulation) corresponds to the emergence and maintenance of rhythmic SWC discharges, even after the DC stimulation is terminated. (b) A 10-second selection from the 90-second simulation in (a) showing the transition (after the DC stimulation) into SWCs. Same convention as in (a), except the bottom panel is a power spectral density (PSD) plot of the 10-second average membrane potential values of the pyramidal cells (third panel-red). The cessation of spiking of interneurons at around t = 46 seconds (first and second panels) as they enter depolarization block (third panel-blue) is concomitant with the emergence of SWCs of increasing amplitude (third panel-red) and the time-frequency spectrum of the average membrane potentials across all the pyramidal neurons in the network (Fig. 3a, lowest panel). We use this average value as a proxy for the LFP activity. The model cortical network during the first 40 seconds of the simulation is in the “low activity” state (in which [K +]o converges to ∼ 3mM). This “low-activity” state at the network level is similar to the single neuron level in panel: average extracellular potassium level values for interneurons (blue) and pyramidal cells (red). Fifth panel: time-frequency plot of the average membrane potential values of the pyramidal cells (third panel). (b) A 30-second selection from the 90-second simulation in a showing both the transient low voltage fast (∼ 20 Hz) oscillations and the subsequent emergence of rhythmic SWC discharges. Same convention as in a, except the bottom panel is a power spectral density (PSD) plot of the 30-second average membrane potential values of the pyramidal cells. Transient (∼ 20 Hz) oscillations are observed after the DC stimulation but before inhibition wears off starting at t = 49.5 seconds (red triangle). As inhibition begins to wear off, the firing rate of interneurons shows an initial transient increase, after which it decreases to zero when the interneurons enter depolarization block (third panel-blue, after around t = 50 seconds). The reduced inhibition (due to both active wearing off of inhibition and the cessation of interneuron spiking) and increased [K+]o (fourth panel) also lead to the increased firing of pyramidal cells and the emergence of rhythmic SWC discharges (third panel-red, towards the end) Fig. 2c where the replenishment of [K B] is insufficient for the glial buffer to maintain a targeted level of [K +]o. At time t = 40s, a DC stimulation (2.5 μA/cm2 added to Iext ) was applied to every neuron in the model network for a short period of 2.5 seconds. This stimulation kick-started the extracellular potassium accumulation process and forced the cortical network to enter into a “high-activity” or “epileptic” state in which the glial buffer attempts to lock into a higher [K +]o level (Fig. 3a, 4th panel), as determined by the parameter [K +]o(th). Clear gamma band (∼ 40 Hz) oscillations are readily observed in the LFP power spectrum after the DC stimulation (Fig. 3a, lowest panel). Importantly, the gamma activity is sustained, even after DC stimulation is terminated, persisting until the end of the simulation. Furthermore, at the individual neuron level, both the fastspiking interneurons and pyramidal cells spike at lower rates of ∼ 3 and 1 spike per second, respectively, than the frequency of the gamma LFP activity (see Fig. 3a, topmost panel). The second panels from the top of Fig. 3a and b show raster plots of 2 of the 81 minicolumns in the simulation. The neuron indices in the raster plots are grouped by neuron types (pyramidal or interneuron). It is also clear that while the LFP proxy (average membrane potential of pyramidal cells) shows gamma oscillations (Fig. 3a and b, 3rd panels from top), the firing pattern of each individual neuron is not indicative of the global oscillation as observed in the LFP proxy. The density plots of the interspike intervals (ISIs) after the DC stimulation across each of the pyramidal and interneuronal populations (Fig. 3c, main plots) reveal a large variance of ISI values. Spiking in FS interneurons showed high irregularity, with the coefficient of variation (CV) of ISIs concentrating around 1 and higher values, while the CV for pyramidal neurons tended to be lower and more broadly distributed (Fig. 3c, inset plots). Irregularity of pyramidal cell spiking can be increased, nevertheless, by moderately increasing the unitary inhibition from the interneurons to the pyramidal cells. However, we observed a competition between the glial parameters (i.e. θ and [K +]o(th)) and the unitary inhibitory conductance values in terms of sustaining high enough [K +]o values to support gamma oscillations. Increasing the unitary inhibitory conductance values affects the sustainability of gamma oscillations. Too strong gsiy→ne or gsiy→ni values make the gamma oscillations less sustainable, with the system immediately settling back into the original “lowactivity” state. A rough mean-field estimate of the synaptic currents entering the pyramidal cells after the DC stimulation suggests that the magnitude of the inhibitory current is a few times higher than that of the excitatory current. This average inhibitory synaptic current entering a pyramidal cell can, for example, be estimated by the formula: Isiy→n e = P (i → e)×gsiy→ne ×Ni × (V−Ei) ×τi × νi , (20) where νi is the spike rate of the inhibitory population and Ni is the number of inhibitory interneurons in the network. Units μF/cm2 mS/cm2 mS/cm2 mS/cm2 mS/cm2 mS/cm2 mV mV μA/cm2 mM mM mM mM mM cm2/ms μA/cm2 mM mM ms−1 mM mM An “S” on the value column denotes that the parameter is dependent on specific simulations. Please refer to the main text for specific values An “S” on the value column denotes that the parameter is dependent on a specific simulation. Please refer to the main text for specific values 3.3 Transition into spike-and-wave complex (SWC) seizures: role of synaptic inhibition and depolarization block In this section we show how variations in synaptic inhibition and depolarization block in the examined neuronal network model can reproduce the dynamics of the observed SWC seizures. In this type of seizures, neural activity transitions into high amplitude rhythmic 2 − 3 Hz LFP discharges. Each SWC event in the LFP lasts for about 300 ms and is characterized by an initial fast “spike” followed by a slow wave potential. Neuronal spiking tends to occur during the LFP “spike” phase and is highly suppressed during the LFP wave phase (e.g. Truccolo et al. 2014) . Preliminary examination of neuronal spiking characteristics has indicated that although both putative principal and FS inhibitory interneurons tend to increase their firing rates in the initial stages of spike-wave seizures (Ahmed et al. 2014) , putative FS interneurons tend to shut down later, just before the emergence of full spike-wave discharges. Ahmed et al. (2014) observed that this termination of activity is preceded by a progressive decreasing of action potential amplitudes in FS interneurons – an indication that FS interneurons enter depolarization block, resulting eventually in the cessation of their activity. Figure 4 shows simulation results for the case where the overall inhibitory synaptic strength is kept at a low level during the entire simulation, leading to SWC seizures. Because the aim is to have FS interneurons to eventually enter depolarization block, we set a higher [K+]o(th) and a lower gK value for the interneurons than in the gamma seizure simulation (Fig. 3). Higher [K+]o(th) values increase FS interneuron excitability, while lower gK values lead FS interneurons to enter depolarization block at a lower [K+]o level. Parameters for model neurons and glial properties of model neurons corresponded to: θ = −0.15, [K+]o(eq) = 3.6, [K+]o(th) = 8.5, IKmax = 1.8, gKCa = 12 for pyramidal cells and θ = −0.15, [K+]o(eq) = 3, [K+]o(th) = 7.75, IKmax = 1.9, gK = 3.5 for interneurons. Unitary inhibitory conductance values gsiy→ni and gsiy→ne were set here to only one-tenth of that in the gamma seizure simulation (Fig. 3). Synaptic parameters corresponded to: gsey→ne = 0.0007, gsey→ni = 0.0007, gsiy→ni = 0.0025, gsiy→ne = 0.0025, ge0 = 0.01016, σe = 0.0025, gi0 = 0.084, σi = 0.02. As before, the simulation starts at a “low-activity” state (t = 0 − 40 seconds) and a DC stimulation (1.6 μA/cm2) is delivered between t = 40 and 42.5 seconds, “kicking” the system into a high-activity “epileptic” state. Since the overall inhibitory synaptic strength is low, but [K+]o(th) is higher than that in the gamma seizure case, FS inhibitory interneurons fire at a higher frequency after the DC stimulation than in Fig. 3, creating a positive feedback which further increases [K+]o. When [K+]o reaches a critical value, inhibitory interneurons enter depolarization block, which eventually leads to cessation of inhibition (at around t = 46 seconds in Fig. 4). We note that, in contrast to the gamma seizures, the attainment of the critical value for [K+]o is possible here because [K+]o(th) is set to a higher level, such that the glial potassium ceiling is higher. After cessation of synaptic inhibition, the network dynamics is dominated by the interaction amongst pyramidal cells and evolve into rhythmic SWC discharges. Throughout, abnormal glial buffering maintains the high [K+]o required for keeping the interneurons in the depolarization block regime, thus creating a long-lasting “seizure” state. In the model simulations, IKCa adaptation currents in the pyramidal cells and sufficiently strong excitatoryexcitatory (gsey→ne) couplings (Van Vreeswijk and Hansel 2001; Dur-e-Ahmad et al. 2012; Nicola and Campbell 2013b; Ferguson et al. 2015) are critical for both the “spike” and “wave” firing suppression phases of the the SWC discharges. 3.4 Transition into seizures consisting of initial low-voltage fast LFP oscillations followed by spike-and-wave discharges Next, we demonstrate the generation of the second type of SWC seizures. In this case, the emergence of rhythmic SWC discharges is preceded by low voltage higher frequency (∼ 10 − 20 Hz) LFP oscillations. Neuronal spiking activity during these fast oscillations tends to be asynchronous as in the previously examined gamma seizures. In other words, while there is a clear fast frequency oscillation at the LFP level, individual neurons spike at a lower rate without obvious synchrony with the global (i.e. LFP) oscillation patterns. The similarity of the firing dynamics of individual neurons between these transient oscillations and gamma seizures leads us to speculate that both phenomena might share a common synaptic and glial mechanism. As shown in Fig. 5, the seizure begins with low voltage fast (∼20 Hz) LFP oscillations which eventually evolve into low-frequency rhythmic SWC discharges. Similarly to the gamma seizures (Fig. 3), the initial fast oscillations in the model require a high overall inhibitory synaptic strength, while the evolution and maintenance of spike-wave discharges require the opposite. To reproduce this type of seizure dynamics, the following sequence of events was implemented in the model: a substantial increase in inhibitory synaptic strength at seizure transition (leading to low voltage fast oscillations) is followed by a rapid decrease, which is hypothesized here to occur via GABA depletion (Zhang et al. 2012) . Once synaptic inhibition is below normal levels, the same path as in the first type of SWC seizures follows: inhibitory interneurons enter depolarization block and the network evolves into rhythmic SWC discharges supported mostly by the population of pyramidal neurons. Specifically, the synaptic parameters gsiy→ni and gsiy→ne were set to constant values during the entire simulation. To increase the effective inhibitory synaptic strength during seizure transition (i.e. during the DC stimulation phase, t = 40 − 42.5s), we let the maximal value of the gating variable per spike, smax , for interneurons (equation set (18)), vary temporally instead of being a constant as in Figs. 3 and 4. The value of smax for interneurons was set to 1 at the beginning of the simulation (Table 5) and kept at this value during the “low activity” phase (t = 0 − 40s). After that, smax for interneurons asymptotically approached a terminal value smtearxm during the DC stimulation phase according to the following formula: smax (t + δt) = smax (t) × where τrise = 62.5ms is the time constant of the increase in smax and smtearxm = 70 is the asymptotic terminal value of smax (i.e. 70 times the value of original smax before the DC stimulation phase). The increased smax value for interneurons was then held constant during the fast oscillations phase (from 42.5 to 49.5 seconds). After that it was decreased to simulate the loss of inhibition. We set smax → 0.1 × smax for individual interneuron after each spike of the same interneuron during this “depletion” phase. (The start of this depletion phase is marked by a red triangle in Fig. 5.) Eventually the inhibition approaches zero, the “depleted” state. In summary, Fig. 5 shows the emergence of clear low voltage fast (∼ 20 Hz) LFP oscillations right after the onset of the DC stimulation. These fast oscillations disappear shortly after the start of the inhibition “depletion” phase (marked by a red triangle in Fig. 5). The strong build up of inhibitory synaptic strength (up to 70 times the initial strength) during the DC stimulation is necessary because it slows the rise of [K+]o during DC stimulation and keeps the interneurons from entering depolarization block. Moreover, this very high level of synaptic inhibition keeps the population frequency in the ∼ 10 − 20 Hz range. A lower, but still large inhibitory strength would have led to the emergence of gamma oscillations, as in Fig. 3. The initial decrease in [K+]o during the transient fast oscillations (Fig. 5a, fourth plot from top) is due to the strong inhibition which limits the firing of inhibitory interneurons. As inhibition begins to wear off at t = 49.5 seconds, the trajectory of [K+]o reverses its course because both pyramidal cells and interneurons fire more as a result of decreased inhibitory synaptic strength. Finally, when [K+]o reaches a sufficiently high level (∼ 7 mM for this simulation), the interneurons enter depolarization block and inhibition ceases, as clearly shown by the complete termination of interneuron firing (Fig. 5b, raster plot at around t = 63 seconds). After cessation of inhibition, SWC discharges emerge supported only by the activity of interacting pyramidal neurons. With the exception of the time varying nature of smax for interneurons and a DC stimulation here set to 2.5 instead of 1.6 μA/cm2, the numerical values of the synaptic (including gsiy→ni and gsiy→ne), intrinsic and glial parameters were the same as those for the first type of SWC seizures (Fig. 4). 4 Discussion Our study is motivated by recent microelectrode array recordings of ensembles of single neurons and high-density LFP activity in human neocortex during propagated focal seizures (Truccolo et al. 2011; Schevon et al. 2012; Ahmed et al. 2014; Truccolo et al. 2014; Wagner et al. 2015) . These recordings suggest different ways via which seizures starting in a focal site can recruit more distal neocortical areas. We used a biophysical cortical network model consisting of conductance-based neurons, coupled with glial buffer for [K+]o, to examine potential mechanisms underlying transitions into ictal states in recruited neocortical areas during seizure spread. Through the analyses and simulations of single neuron and network models, we have demonstrated that neocortical networks can be made susceptible to seizures as a result of abnormal glial potassium buffering and changes in synaptic inhibition. First, at the single neuron level, the imbalance in the conversion between potassium buffer in the bound and unbound forms leads to bi-stability such that one state has a lower [K+]o level (thus lower excitability) than the other. This bi-stability also carries over to the network level. Perturbation either by noise or by DC stimulation can “kick-start” the high [K+]o hyper-excitable state from the low [K+]o state. Second, we then showed that the three major types of seizure transitions as observed in the human microelectrode array recordings can be reproduced when these pathological glial potassium dynamics interact with various synaptic and intrinsic parameter settings for inhibitory interneurons. In our simulations, transitions into gamma (∼ 30 − 60 Hz) ictal activity were obtained with a high inhibitory conductance between cortical neurons. Such high inhibitory synaptic strength was required to prevent the interneurons from entering depolarization block during the seizure state. Moreover, since high inhibitory strength tends to decrease [K+]o, the long-lasting gamma seizure state was the result of a subtle balance between high inhibitory conductance values and the accumulation of [K+]o due to abnormal glial potassium buffering dynamics. For transitions into SWC ictal activity, the cortical network can either experience a generally low value of inhibitory conductance between cortical neurons, or an initial high inhibitory synaptic strength followed by inhibition breakdown. In the first scenario, the transition into SWC was direct and involved interneurons entering depolarization block because of the low inhibitory synaptic strength. As soon as the interneurons entered depolarization block, inhibition ceased. At that point, neuronal spiking activity supported only by interacting pyramidal cells led to recurring SWC discharges. In the second SWC scenario, there was a buildup of inhibitory synaptic strength prior to the transition into spike-wave discharges. The result of this inhibition buildup is the occurrence of transient lowamplitude fast oscillations around 10 − 20 Hz. This buildup of inhibitory synaptic strength was followed by a wearing off of inhibition (e.g. through GABA depletion), thus providing a low inhibition environment for interneurons to enter depolarization block and the population of pyramidal cells to eventually sustain recurring SWC discharges. 4.1 Abnormal glial buffering and its relationship to seizure In this work we have established a computational model of seizure through the occurrence of bi-stability at the single neuron level. This bi-stability is the result of a hypothesized pathological glial potassium dynamics in which there is an imbalance between potassium buffers in their bound and unbound forms. Bi- (or multi)stability, either at the single neuron or at the network level, has been implicated in numerous physiological and pathological conditions in the brain (Van Ooyen et al. 1992; Sasaki et al. 2007; Freyer et al. 2009; Fro¨ hlich et al. 2010; Anderson et al. 2012; Ho et al. 2012; He 2014; Ho et al. 2014; Hu¨ bel et al. 2014) . In the particular case of epileptic seizures, Fro¨ hlich et al. (2010) have demonstrated in computational models that it is possible for network bistability to exist even when the glial potassium dynamics are within a physiological range (through a higher absolute value of θ , the parameter determining the bias favouring forward versus backward binding rate). In this scenario, Fro¨ hlich et al. (2010) associate higher [K +]o states with epileptic seizures. In the models used here, we did not observe bi-stability either at the single neuron or network level when θ was set to physiological values. We only observed bi-stability with pathological (i.e. smaller in absolute value) values of θ . We have singled out (via single neuron simulations) the differential in transition rates between [B] and [K B] as the important factor underlying such bi-stability. The fact that Fro¨ hlich et al. (2010) used two-compartment neuron models, while we used single-compartment, might explain this difference between the two studies. Moreover, different synaptic connection parameters were used by Fro¨ hlich et al. (2010). Different synaptic connectivity might affect the length and strength of perturbation required to elicit the bi-stable behaviour, if at all possible. A more thorough examination of how glial potassium buffering may affect neuronal dynamics should involve detailed bifurcation analyses of the combined neuron-glial system (Touboul 2008; Nicola and Campbell 2013a; Kim and Nykamp 2014; Nicola and Campbell 2014) . Such analyses could determine the intrinsic and glial parameters that allow multistability. From the existence and stability of these states, one would expect to gain some insight into how resilient the network is to seizure initiation and propagation. Although we have focused on their role on seizure propagation, potassium and glial buffering may also play an important role in seizure termination (Kramer et al. 2012; Gonza´lez-Ram´ırez et al. 2015) . 4.2 Changes in inhibitory and excitatory conductances Our focus on changes in inhibition was motivated by several previous (in vivo/vitro) animal models (Trevelyan et al. 2006; Zhang et al. 2012; Grasse et al. 2013; Zˇ iburkus et al. 2013; Uva et al. 2015) and in vitro studies of human epileptic cortical tissue (e.g. D’Antuono et al. 2004), where inhibitory interneuron activity appears to change significantly during preictal and initial periods of the seizure, before major changes in principal cells are detected. Recent work by Zˇ iburkus et al. (2013) has also examined in detail how changes in both inhibitory and excitatory conductances may lead to imbalances in excitation/inhibition preceding and during seizure like events in hippocampal slices under the potassium channel blocker 4-aminopyridine and reduced extracellular magnesium. In that respect, we emphasize that the changes in inhibition as implemented in our simulations did not only affect inhibitory conductances. Changes in inhibition in our model were implemented in part by systematically varying the levels of unitary synaptic inhibition, which in turn affected several other network properties such as the number of neurons generating spikes at any given time, and thus firing rates in both populations of pyramidal neurons and FS interneurons. In this way, by varying unitary synaptic inhibition, the E-I balance in the simulated networks also varied dynamically, with both excitatory and inhibitory conductances being altered. In addition, after transitions into spike-wave discharges in the simulated SWC seizures (Figs. 4 and 5), these discharges were supported only by the synaptic interactions among pyramidal neurons. Furthermore, we note that there are a couple of important differences between the animal model in Zˇ iburkus et al. (2006 , 2013) and our model simulations. First, apart from temporary depolarization block periods, the activity of inhibitory neurons remained throughout the seizures (i.e. until seizure termination) in Zˇ iburkus et al.’s model. In our computational model, it is only during gamma seizures that the spiking activity of FS inhibitory interneurons is preserved throughout the seizure and depolarization block is absent. By contrast, FS interneurons enter depolarization block and remain in that state during the SWC seizures. Second, while in Zˇ iburkus et al.’s study transient depolarization block events in pyramidal neurons are also present and play an important role, in our computational model pyramidal neurons tend to be robust to depolarization block both in SWC and gamma seizures. This robustness resulted primarily from the slow calcium-mediated potassium afterhyperpolarization (AHP) currents (IKCa, Eq. (7)) in the model pyramidal neurons, which prevented them from evolving into regimes of high enough firing rates required for depolarization block. These features of the biophysical model examined here seem consistent with the dynamics of propagated seizures (in neocortical patches distal from putative onset areas) observed in our human data. First, both inhibitory and excitatory spiking activity appears to be preserved throughout the gamma seizures (Truccolo et al. 2011) . Second, in SWC seizures, spiking in putative FS interneurons appears to terminate before the transition into large amplitude spike wave discharges. This cessation of FS interneuron spiking activity is accompanied by signatures of depolarization block, while no similar signatures were found in putative principal cells (Ahmed et al. 2014) . 4.3 Network inhibition and its relationship to seizure propagation Using biophysical neuronal network models of neocortical patches, we have demonstrated how FS inhibitory interneurons may play an important role in three types of observed human seizures. In particular, we have quantified the role of FS interneurons by varying the inhibitory synaptic strength and controlling how easily interneurons enter depolarization block. For the two types of SWC seizures modeled here (Figs. 4 and 5), ultimately it is the failure of inhibition (when the FS inhibitory interneurons enter depolarization block) that leads the network to transition into SWC activity. We examined two main scenarios that led to the failure of inhibition, resulting in SWC seizures. In the first scenario, FS inhibitory interneurons enter depolarization block when there is an overall low inhibitory synaptic strength. When interneurons are excited during seizure initiation, the lower inhibitory synaptic strength allows interneurons to fire at a higher frequency, thus hastening the [K+]o accumulation process and leading eventually to depolarization block (Fig. 4). In the second scenario, an initially large inhibitory synaptic strength is followed by GABA depletion (Zhang et al. 2012) . Initially, interneurons are mutually inhibited by the large inhibitory synaptic strength and their firing frequencies are too low to trigger [K+]o accumulation. If such level of inhibition was not perturbed, seizure propagation would fail. The decreased effective inhibitory strength resulting from subsequent GABA depletion allows for higher firing rates in interneurons (and principal cells) to develop, again leading eventually to depolarization block in inhibitory interneurons (Fig. 5). Some previous studies have emphasized the role of feedforward inhibition in seizure propagation in neocortex (Trevelyan et al. 2006; Schevon et al. 2012; for a review see Paz and Huguenard 2015) . The main idea is that the propagation of initially localized ictal activity to more distal cortical sites involves predominantly feedforward inhibition. This feedforward inhibitory drive could be a mechanism for the initial strengthening of inhibition during seizure propagation assumed in our model in two of the three seizure propagation scenarios (gamma and SWC with preceding fast oscillations). Furthermore, Trevelyan and colleagues (Trevelyan et al. 2006; Schevon et al. 2012) have emphasized the role of an “inhibitory restraint” or “veto” mechanism during neocortical propagation of focal seizures. Inhibitory veto of seizure propagation has been shown in low Mg2+ mouse models of epilepsy (Trevelyan et al. 2006) and has been argued to be present also in human focal seizures (Schevon et al. 2012) . This inhibitory veto mechanism would lead to “penumbra areas” ahead of the ictal wavefront. The failure of such inhibitory veto would allow the successful propagation of the ictal wavefront into new recruited ictal areas. Otherwise, the inhibition would be strong enough to contain the seizure spread. In this way, one could also argue that the gamma seizures examined here would represent a failure in seizure propagation. This possibility remains an open and important experimental question. Given that these sustained gamma oscillations reflect abnormal dynamics during secondarily generalized focal seizures, we currently consider this activity as propagated seizures. The gamma seizures in our model (Fig. 3) represent a situation in which there is a balance between enhanced inhibition and [K+]o accumulation. In other words, the seizure propagates while inhibitory interneuronal activity is preserved. In a related experimental study, our group (Lu et al. 2015) has shown that constant, strong enough depolarization driven by optogenetic stimulation can induce sustained (∼ 50 Hz) gamma oscillations in healthy primate motor cortex. In view of the strong inhibition in simulations of gamma seizures (Fig. 3), we do not preclude the possibility that some of the “high-activity” states in the simulated gamma seizures may in fact be long-lasting transients as opposed to asymptotic stable states. In other words, the elevated [K+]o in the “high-activity” states may nevertheless settle back to a much lower value after a long enough time. However, this issue has little effect on our main claims in this paper, given the long time scale (at least in the order of minutes) over which [K+]o remains elevated and increasing. 4.4 Model limitations and future work Our analyses and simulations show how different types of seizure transitions observed in intracortical recordings of propagated human focal seizures can be accounted for by variations in synaptic inhibition and extracellular potassium concentration. Several limitations of the model, however, warrant future improvement. First, although the interaction of [K +]o and inhibitory synaptic dynamics can account for the phenomena observed in our human epilepsy data (Ahmed et al. 2014) , several other ionic or synaptic mechanisms are possible and have been studied in different animal models of epilepsy (e.g. Alfonsa et al. 2015) . Second, our model exploration was restricted to the transition into ictal states. In particular, we did not examine how seizures terminate, which is also a complex topic and may involve multiple ionic and synaptic components (such as glutamate depletion–Lado and Moshe´ 2008; GABA upregulation–Wen et al. 2015; fluctuations in ionic concentrations–Krishnan and Bazhenov 2011; Kramer et al. 2012) . A more detailed model, including more complex ionic dynamics, will be required to address this issue. Furthermore, gap junction effects were not included in the model. Although electrical synapses have been suggested as an underlying mechanism for interneuronal bursting (Skinner et al. 1999) and synchronous population gamma oscillations (Traub et al. 2001) , neither of these two features appeared to be prominent in our human data where neuronal spiking during gamma seizures is largely asynchronous. Thus, it remains an open question whether gap junction effects are critical to account for transitions to and maintenance of gamma seizures. Nevertheless, gap junctions may still play an important role during SWC discharges as previously hypothesized by Traub and colleagues (e.g. Traub et al. 2005) . We hope to address these open issues in the future. Third, in the hypothesized contributing factors to seizure transitions examined here, we have remained agnostic about what ultimately leads to transient dysfunctions in inhibitory activity and in potassium glial-buffering processes. We note that several other factors involving changes in oxygen, Na+/K+ ATP pumps and cell volume (e.g. Ingram et al. 2014; Wei et al. 2014a) , for example, can also contribute to elevations in extracellular potassium. We also emphasize that the goal in this study is to replicate the dynamics observed in propagated seizures, in other words, ictal dynamics in neocortical patches distal to the putative seizure onset zones. As described in Truccolo et al. (2011, 2014) and Wagner et al. (2015) , our recorded neocortical patches were close, but distal to the putative seizure onset zones. The dynamics involving seizure initiation in seizure onset/focus areas might involve different mechanisms. As stated above, we have focused on fast time-scale changes in synaptic inhibition which have been reported in several previous animal studies and in vitro studies of human epileptic cortical tissue. Although these studies have provided some initial evidence for transient changes in inhibitory interneuron activity preceding seizures, it remains a major open question how these dysfunctions in interneuron network activity (and potassium glial buffering) would arise in the longer time scale of epileptogenesis. Acknowledgments The authors would like to thank the editors, reviewer, Drs. F. Skinner, Liang Zhang, Felipe Gerhard and Timothe´e Proix for constructive comments on the manuscript. We also thank Drs. Sydney Cash, Omar Ahmed and Mark Kramer for valuable discussions. This research was conducted using computational resources and services at the Center for Computation and Visualization, Brown University. Compliance with Ethical Standards Conflict of interest The authors declare that they have no conflict of interest. Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http:// creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. Ahmed , O.J. , Kramer , M.A. , Truccolo , W. , Naftulin , J.S. , Donoghue , J.A. , Eskandar , E.N. , Cosgrove , G.R. , Blum , A.S. , Potter , N.S. , Hochberg , L.R. , & Cash , S.S. ( 2014 ). Inhibitory single neuron control of seizures and epileptic traveling waves in humans . BMC Neuroscience , 15 ( Suppl 1 ), F3 . Alfonsa , H. , Merricks , E.M. , Codadu , N.K. , Cunningham , M.O. , Deisseroth , K. , Racca , C. , & Trevelyan , A.J. ( 2015 ). The contribution of raised intraneuronal chloride to epileptic network activity . The Journal of Neuroscience , 35 ( 20 ), 7715 - 7726 . Anderson , W.S. , Kudela , P. , Cho , J. , Bergey , G.K. , & Franaszczuk , P.J. ( 2007 ). Studies of stimulus parameters for seizure disruption using neural network simulations . Biological Cybernetics , 97 , 173 - 194 . Anderson , W. , Azhar , F. , Kudela , P. , Bergey , G. , & Franaszczuk , P. ( 2012 ). Epileptic seizures from abnormal networks: why some seizures defy predictability . Epilepsy Research , 99 , 202 - 213 . Bazhenov , M. , Timofeev , I. , Steriade , M. , & Sejnowski, T.J. ( 2004 ). Potassium model for slow (2-3 Hz) in vivo neocortical paroxysmal oscillations . Journal of Neurophysiology , 92 , 1116 - 1132 . Brunel , N. ( 2000 ). Dynamics of sparsely connected networks of excitatory and inhibitory spiking neurons . Journal of Computational Neuroscience , 8 , 183 - 208 . Brunel , N. , & Wang , X.J. ( 2003 ). What determines the frequency of fast network oscillations with irregular neural discharges? I. Synaptic dynamics and excitation-inhibition balance . Journal of Neurophysiology , 90 , 415 - 430 . Buxhoeveden , D.P. , & Casanova , M.F. ( 2002 ). The minicolumn hypothesis in neuroscience . Brain: A Journal of Neurology , 125 , 935 - 951 . Cammarota , M. , Losi , G. , Chiavegato , A. , Zonta , M. , & Carmignoto , G. ( 2013 ). Fast spiking interneuron control of seizure propagation in a cortical slice model of focal epilepsy . The Journal of Physiology , 591 ( 4 ), 807 - 822 . Cressman , J.R. Jr ., Ullah , G. , Ziburkus , J. , Schiff , S. , & Barreto , E. ( 2009 ). The influence of sodium and potassium dynamics on excitability, seizures, and the stability of persistent states: I. Single neuron dyamics . Journal of Computational Neuroscience , 26 , 159 - 170 . D'Antuono , M. , Louvel , J. , Ko¨hling, R., Mattia , D. , Bernasconi , A. , Olivier , A. , Turak , B. , Devaux , A. , Pumain , R. , & Avoli , M. ( 2004 ). Gabaa receptor-dependent synchronization leads to ictogenesis in the human dysplastic cortex . Brain: A Journal of Neurology , 127 , 1626 - 1640 . Dur- e-Ahmad, M. , Nicola , W., Campbell , S. , & Skinner , F. ( 2012 ). Network bursting using experimentally constrained single compartment CA3 hippocampal neuron models with adaptation . Journal of Computational Neuroscience , 33 ( 1 ), 21 - 40 . Ermentrout , B. ( 2002 ). Simulating, analyzing, and animating dynamical systems: A guide to XPPAUT for researchers and students (SIAM). Ferguson , K. , Njap , F. , Nicola , W. , Skinner , F. , & Campbell , S. ( 2015 ). Examining the limits of cellular adaptation bursting mechanisms in biologically-based excitatory networks of the hippocampus . Journal of Computational Neuroscience , 39 , 289 - 309 . Freyer , F. , Aquino , K. , Robinson , P.A. , Ritter , P. , & Breakspear , M. ( 2009 ). Bistability and non-Gaussian fluctuations in spontaneous cortical activity . The Journal of Neuroscience , 29 ( 26 ), 8512 - 8524 . Fro ¨hlich, F. , Bazhenov , M. , Iragui-Madoz , V. , & Sejnowski, T.J. ( 2008 ). Potassium dynamics in the epileptic cortex: New insights on an old topic . The Neuroscientist , 14 , 422 . Fro ¨hlich, F., Sejnowski, T.J., & Bazhenov , M. ( 2010 ). Network bistability mediates spontaneous transitions between normal and pathological brain states . The Journal of Neuroscience , 30 ( 32 ), 10734 - 10743 . Gonza´ lez-Ram´ırez, L.R. , Ahmed , O.J. , Cash , S.S. , Wayne , C.E. , & Kramer , M.A. ( 2015 ). A biologically constrained, mathematical model of cortical wave propagation preceding seizure termination . PLoS Computational Biology , 11 ( 2 ), e1004065 . Grasse , D. , Karunakaran , S. , & Moxon , K. ( 2013 ). Neuronal synchrony and the transition to spontaneous seizures . Experimental Neurology , 248 , 72 - 84 . He , B.J. ( 2014 ). Scale-free brain activity: past, present, and future . Trends in Cognitive Sciences , 18 , 480 - 487 . Ho , E.C.Y. , Stru¨ber, M. , Bartos , M. , Zhang , L. , & Skinner , F.K. ( 2012 ). Inhibitory networks of fast-spiking interneurons generate slow population activities due to excitatory fluctuations and network multistability . The Journal of Neuroscience , 32 ( 29 ), 9931 - 9946 . Ho , E.C.Y. , Eubanks , J.H. , Zhang, L. , & Skinner , F.K. ( 2014 ). Network models predict that reduced excitatory fluctuations can give rise to hippocampal network hyper-excitability in meCP2-null mice . PloS One , 9 , e91148 . Hu¨bel , N., Scholl , E. , & Dahlem , M.A. ( 2014 ). Bistable dynamics underlying excitability of ion homeostasis in neuron models . PLoS Computational Biology , 10 , e1003551 . Ingram , J. , Zhang, C. , Cressman , J. , Hazra , A. , Wei , Y. , Koo , Y. , Zˇ Iburkus , J. , Kopelman , R. , Xu , J. , & Schiff , S. ( 2014 ). Oxygen and seizure dynamics: I. experiments . Journal of Neurophysiology , 112 , 205 - 212 . Jirsa , V.K. , Stacey , W.C. , Quilichini , P.P. , Ivanov , A.I. , & Bernard , C. ( 2014 ). On the nature of seizure dynamics . Brain: A Journal of Neurology , 133 , 2210 - 2230 . Kager , H. , Wadman , W.J. , & Somjen , G.G. ( 2000 ). Simulated seizures and spreading depression in a neuron model incorporating interstitial space and ion concentrations . Journal of Neurophysiology , 84 , 495 - 512 . Kim , C. , & Nykamp , D. ( 2014 ). Dynamics of a network of excitatory and inhibitory neurons induced by depolarization block . BMC Neuroscience , 15 ( Suppl 1 ), P76 . Kramer , M.A. , Truccolo , W. , Eden , U.T., Lepage , K.Q. , Hochberg , L.R. , Eskandar , E.N. , Madsen , J.R. , Lee , J.W. , Maheshwari , A. , Halgren , E. , Chu , C.J. , & Cash , S.S. ( 2012 ). Human seizures selfterminate across spatial scales via a critical transition . Proc Natl Acad Sci U S A , 109 , 21116 - 21121 . Krishnan , G.P. , & Bazhenov , M. ( 2011 ). Ionic dynamics mediate spontaneous termination of seizures and postictal depression state . The Journal of Neuroscience , 31 ( 24 ), 8870 - 8882 . Krishnan , G. , Filatov , G. , & Bazhenov , M. ( 2013 ). Dynamics of high-frequency synchronization during seizures . Journal of Neurophysiology , 109 , 2423 - 2437 . Lado , F.A. , & Moshe´, S.L. ( 2008 ). How do seizures stop? Epilepsia , 49 , 1651 - 1664 . Lu , Y. , Truccolo , W. , Wagner , F. , Vargas-Irwin , C. , Ozden , I. , Zimmermann , J., May, T. , Agha , N. , Wang , J. , & Nurmikko , A. ( 2015 ). Optogenetically induced spatiotemporal gamma oscillations and neuronal spiking activity in primate motor cortex . Journal of Neurophysiology , 113 , 3574 - 3587 . Naze , S. , Bernard , C. , & Jirsa , V. ( 2015 ). Computational modeling of seizure dynamics using coupled neuronal networks: Factors shaping epileptiform activity . PLoS Computational Biology , 11 , e1004209 . Nicola , W. , & Campbell , S. ( 2013a ). Bifurcations of large networks of two-dimensional integrate and fire neurons . Journal of Computational Neuroscience , 35 , 87 - 108 . Nicola , W. , & Campbell , S.A. ( 2013b ). Mean-field models for heterogeneous networks of two-dimensional integrate and fire neurons . Frontiers in Computational Neuroscience , 7 ( 184 ), 1 - 25 . Nicola , W. , & Campbell , S. ( 2014 ). Non-smooth bifurcations of mean field systems of two-dimensional integrate and fire neurons . arXiv:1408 . 4767 . Park , E.H. , & Durand , D.M. ( 2006 ). Role of potassium lateral diffusion in non-synaptic epilepsy: a computational study . Journal of Theoretical Biology , 238 , 666 - 682 . Paz , J. , & Huguenard , J.R. ( 2015 ). Microcircuits and their interactions in epilepsy: is the focus out of focus? Nature Neuroscience , 18 , 351 - 359 . Perucca , P. , Dubeau , F. , & Gotman , J. ( 2014 ). Intracranial electroencephalographic seizure-onset patterns: effect of underlying pathology . Brain: A Journal of Neurology , 137 , 183 - 196 . Richardson , K. , Fanselow , E. , & Connors , B. ( 2008 ). Neocortical anatomy and physiology , In Engel, J., & Pedley , T . (Eds.) Epilepsy: A Comprehensive Textbook (pp. 323 - 336 ): LippincottWilliams & Wilkins. Rudolph , M. , Piwkowska , Z. , Badoual , M. , Bal , T. , & Destexhe , A. ( 2004 ). A method to estimate synaptic conductances from membrane potential fluctuations . Journal of Neurophysiology , 91 ( 6 ), 2884 - 2896 . Sasaki , T. , Matsuki , N. , & Ikegaya , Y. ( 2007 ). Metastability of active CA3 networks . The Journal of Neuroscience , 27 ( 3 ), 517 - 528 . Schevon , C.A. , Weiss , S.A. , McKhann , G. Jr. , Goodman , R.R. , Yuste , R. , Emerson , R.G. , & Trevelyan , A.J. ( 2012 ). Evidence Of an inhibitory restraint of seizure activity in humans . Nature Communications , 3 ( 1060 ), 1 - 11 . Skinner , F.K. ( 2006 ). Conductance-based models . Scholarpedia 1 :1408 revision # 125663 . Skinner , F.K. , Zhang , L. , Perez-Velazquez , J.L. , & Carlen , P.L. ( 1999 ). Bursting in inhibitory interneuronal networks: a role for gap-junctional coupling . Journal of Neurophysiology , 81 , 1274 - 1283 . Sritharan , D. , & Sarma , S.V. ( 2014 ). Fragility in dynamic networks: application to neural networks in the epileptic cortex . Neural Computation , 26 ( 10 ), 2294 - 2327 . Thurman , D.J. , Beghi , E. , Begley , C.E. , Berg , A.T. , Buchhalter , J.R. , Ding , D. , Hesdorffer , D.C. , Hauser , W.A. , Kazis , L. , Kobau , R. , Kroner , B. , Labiner , D. , Liow , K. , Logroscino , G. , Medina , M.T. , Newton , C.R. , Parko , K. , Paschal , A. , Preux , P.M. , Sander , J.W. , Selassie , A. , Theodore , W. , Tomson , T. , & Wiebe , S. ( 2011 ). Standards for epidemiologic studies and surveillance of epilepsy . Epilepsia , 52 , 2 - 26 . Touboul , J. ( 2008 ). Bifurcation analysis of a general class of nonlinear integrate-and-fire neurons . SIAM Journal of Applied Mathematics , 68 , 1045 - 1079 . Traub , R.D. , Kopell , N. , Bibbig , A. , Buhl , E.H. , LeBeau , F.E.N. , & Whittington , M.A. ( 2001 ). Gap junctions between interneuron dendrites can enhance synchrony of gamma oscillations in distributed networks . The Journal of Neuroscience , 21 ( 23 ), 9478 - 9486 . Traub , R.D. , Contreras , D. , Cunningham , M.O. , Murray , H. , LeBeau , F.E.N. , Roopun , A. , Bibbig , A. , Wilent , W.B. , Higley , M.J. , & Whittington , M.A. ( 2005 ). Single-column thalamocortical network model exhibiting gamma oscillations, sleep spindles, and epileptogenic bursts . Journal of Neurophysiology , 93 ( 4 ), 2194 - 2232 . Trevelyan , A.J. , Sussillo , D. , Watson , B.O. , & Yuste , R. ( 2006 ). Modular propagation of epileptiform activity: evidence for an inhibitory veto in neocortex . The Journal of Neuroscience , 26 ( 48 ), 12447 - 12455 . Truccolo , W. , Donoghue , J.A. , Hochberg , L.R. , Eskandar , E.N. , Madsen , J.R. , Anderson , W.S. , Brown , E.N. , Halgren , E. , & Cash , S.S. ( 2011 ). Single-neuron dynamics in human focal epilepsy . Nature Neuroscience , 14 ( 5 ), 635 - 643 . Truccolo , W. , Ahmed , O.J. , Harrison , M.T. , Eskandar , E.N. , Cosgrove , G.R. , Madsen , J.R. , Blum , A.S. , Potter , N.S. , Hochberg , L.R. , & Cash , S.S. ( 2014 ). Neuronal ensemble synchrony during human focal seizures . The Journal of Neuroscience , 34 ( 30 ), 9927 - 9944 . Ullah , G. , Cressman , J.R. Jr ., Barreto , E. , & Schiff , S.J. ( 2009 ). The influence of sodium and potassium dynamics on excitability, seizures, and the stability of persistent states: II. Network and glial dynamics . Journal of Computational Neuroscience , 29 , 171 - 183 . Uva , L. , Breschi , G.L. , Gnatkovsky , V. , Taverna , S. , & de Curtis, M. ( 2015 ). Synchronous inhibitory potentials precede seizure-like events in acute models of focal limbic seizures . The Journal of Neuroscience , 35 , 3048 - 3055 . Van Ooyen , A. , Van Pelt , J. , Corner , M.A. , & Lopes da Silva , F.H. ( 1992 ). The emergence of long-lasting transients of activity in simple neural networks . Biological Cybernetics , 67 , 269 - 277 . Van Vreeswijk , C. , & Hansel , D. ( 2001 ). Patterns of synchrony in neural networks with spike adaptation . Neural Computation , 13 , 959 - 992 . Wagner , F. , Eskandar , E. , Cosgrove , G. , Madsen , J. , Blum , A. , Potter , N. , Hochberg , L. , Cash , S. , & Truccolo , W. ( 2015 ). Microscale spatiotemporal dynamics during neocortical propagation of human focal seizures . Neuroimage , 122 , 114 - 130 . Wang , Y. , Goodfellow , M. , Taylor , P.N. , & Baier , G. ( 2014 ). Dynamic mechanisms of neocortical focal seizure onset . PLoS Computational Biology , 10 , e1003787 . Wang , X.J. , & Buzsa´ki, G. ( 1996 ). Gamma oscillation by synaptic inhibition in a hippocampal interneuronal network model . The Journal of Neuroscience , 16 ( 20 ), 6402 - 6413 . Wei , Y. , Ullah , G. , & Schiff , S. ( 2014a ). Oxygen and seizure dynamics: II. computational modeling . Journal of Neurophysiology , 112 , 213 - 223 . Wei , Y. , Ullah , G. , & Schiff , S. ( 2014b ). Unification of neuronal spikes, seizures, and spreading depression . The Journal of Neuroscience , 34 , 11733 - 11743 . Wen , B. , Qian , H. , Feng , J. , Ge , R.J. , Xu , X. , Cui , Z.Q. , Zhu , R.Y., Pan , L.S. , Lin , Z.P. , & Wang , J.H. ( 2015 ). A portion of inhibitory neurons in human temporal lobe epilepsy are functionally upregulated: an endogenous mechanism for seizure termination . CNS Neuroscience & Therapeutics , 21 , 204 - 214 . Zhang , Z.J. , Koifman , J. , Shin , D.S. , Ye , H. , Florez , C.M. , Zhang , L. , Valiante , T.A. , & Carlen , P.L. ( 2012 ). Transition to seizure: Ictal discharge is preceded by exhausted presynaptic GABA release in the hippocampal CA3 region . The Journal of Neuroscience , 32 ( 7 ), 2499 - 2512 . Zˇiburkus , J., Cressman , J. , Barreto , E. , & Schiff , S. ( 2006 ). Interneuron and pyramidal cell interplay during in vitro seizure-like events . Journal of Neurophysiology , 95 , 3948 - 3954 . Zˇiburkus , J., Cressman , J.R. , & Schiff , S.J. ( 2013 ). Seizures as imbalanced up states: excitatory and inhibitory conductances during seizure-like events . Journal of Neurophysiology , 109 , 1296 - 1306 .


This is a preview of a remote PDF: https://link.springer.com/content/pdf/10.1007%2Fs10827-016-0615-7.pdf

E. C. Y. Ho, Wilson Truccolo. Interaction between synaptic inhibition and glial-potassium dynamics leads to diverse seizure transition modes in biophysical models of human focal seizures, Journal of Computational Neuroscience, 2016, 225-244, DOI: 10.1007/s10827-016-0615-7