Modeling implicates inhibitory network bistability as an underpinning of seizure initiation

A plethora of recent experimental literature implicates the abrupt, synchronous activation of GABAergic interneurons in driving the sudden change in brain activity that heralds seizure initiation. However, the mechanisms predisposing an inhibitory network toward sudden coherence specifically during ictogenesis remain unknown. We address this question by comparing simulated inhibitory networks containing control interneurons and networks containing hyper-excitable interneurons modeled to mimic treatment with 4-Aminopyridine (4-AP), an agent commonly used to model seizures in vivo and in vitro. Our in silico study demonstrates that model inhibitory networks with 4-AP interneurons are more prone than their control counterparts to exist in a bistable state in which asynchronously firing networks can abruptly transition into synchrony due to a brief perturbation. We further show that perturbations driving this transition could reasonably arise in vivo based on models of background excitatory synaptic activity in the cortex. Thus, these results propose a mechanism by which an inhibitory network can transition from incoherent to coherent dynamics in a fashion that may precipitate seizure as a downstream effect. Moreover, this mechanism specifically explains why inhibitory networks containing hyper-excitable interneurons are more vulnerable to this state change, and how such networks can undergo this transition without a permanent change in the drive to the system. This, in turn, potentially explains such networks’ increased vulnerability to seizure initiated by GABAergic activity. Author summary For decades, the study of epilepsy has focused on the hypothesis that over-excitation or dis-inhibition of pyramidal neurons underlies the transition from normal brain activity to seizure. However, a variety of recent experimental findings have implicated a sudden synchronous burst of activity amongst inhibitory interneurons in driving this transition. Given the counter-intuitive nature of these findings and the correspondingly novel hypothesis of seizure generation, the articulation of a feasible mechanism of action underlying this dynamic is of paramount importance for this theory’s viability. Here, we use computational techniques, particularly the concept of bistability in the context of dynamical systems, to propose a mechanism for the necessary first step in such a process: the sudden synchronization of a network of inhibitory interneurons. This is the first detailed proposal of a computational mechanism explaining any aspect of this hypothesis of which we are aware. By articulating a mechanism that not only underlies this transition, but does so in a fashion explaining why ictogenic networks might be more prone to this behavior, we provide critical support for this novel hypothesis of seizure generation and potential insight into the larger question of why individuals with epilepsy are particularly vulnerable to seizure.

explains why inhibitory networks containing hyper-excitable interneurons are more vulnerable to this state change, and how such networks can undergo this transition without a permanent change in the drive to the system. This, in turn, potentially explains such networks' increased vulnerability to seizure initiated by GABAergic activity.

Author summary
For decades, the study of epilepsy has focused on the hypothesis that over-excitation or 1 dis-inhibition of pyramidal neurons underlies the transition from normal brain activity 2 to seizure. However, a variety of recent experimental findings have implicated a sudden 3 synchronous burst of activity amongst inhibitory interneurons in driving this transition. 4 Given the counter-intuitive nature of these findings and the correspondingly novel 5 hypothesis of seizure generation, the articulation of a feasible mechanism of action 6 underlying this dynamic is of paramount importance for this theory's viability. Here, we Epilepsy research is divided into studies focusing on seizure initiation, propagation, or termination (illustrated by the rows in this figure and the example studies cited). This paper is interested in seizure initiation in the context of a "GABAergic initiation hypothesis", schematized on the top row of the figure. The focus of the current work is the sudden transition of interneurons into synchrony, as the articulation of a potential mechanism explaining this transition, alongside a justification as to why networks in a seizure state are more vulnerable to this transition, should be identified in order for this overall hypothesis of seizure initiation to be viable. environment to suddenly transition into synchrony, the necessary initial step in this 33 hypothesis. We thus focus on the earliest time in the transition to seizure and not 34 aspects of propagation and termination. 35 The study of inhibitory network synchrony is decades old, dating back to the work of 36 Wang and Rinzel [14]. Various mechanisms have been proposed to explain the 37 generation of oscillations in purely inhibitory networks, the most prominent of which 38 may be the Interneuron Network Gamma (ING) mechanism [15][16][17][18][19][20]. Previous work has 39 shown that inhibitory networks built to examine population activity in an in vitro 40 hippocampal preparation manifest "sharp transitions" into coherent population activity 41 caused by a small, permanent increase to the external drive to the network [21]. an analogue to this increased drive, potentially explaining why hyper-excitable systems 53 transition to synchronous states of a seizure or an inter-ictal spike (IIS). Consistent with 54 this suggestion is that interneuronal firing increases before pyramidal cell firing prior to 55 IIS generation and seizures in animal epilepsy models [1,[26][27][28][29][30]. Analogously in humans, 56 putative interneurons increase their firing before seizure onset [13]. 57 The existing computational insights into inhibitory network synchrony, combined 58 with the experimental literature implicating this dynamic in seizure initiation, motivate 59 this computational study. To explore the role of interneuronal synchrony in seizure 60 initiation, randomly connected, purely inhibitory network models are developed. These 61 networks utilized cell models mimicking properties exhibited by neurons treated with 62 4-Aminopyridine (4-AP), a commonly used experimental model to generate 63 seizures [31,32], or properties of a healthy, control interneuron. Utilizing these tools, 64 this investigation articulates a mechanism explaining how a sudden transition from 65 asynchronous to synchronous firing might arise in an inhibitory network that also offers 66 an explanation for the predisposition of hyper-excitable networks towards this transition. 67 This mechanism was uncovered by comparing the tendency of control and 4-AP 68 inhibitory networks to synchronize, both from random initial conditions and following a 69 perturbation biasing the network towards synchrony. This revealed that 4-AP networks 70 are much more likely to transition from asynchronous to synchronous dynamics 71 following a perturbation, due to a significantly larger regime of network parameters 72 supporting bistability. The existence of a "bistable transition" driving an inhibitory 73 network into synchrony expands upon existing literature probing such mechanisms, 74 especially in the context of epilepsy. As control models do not exhibit the same 75 predisposition for "bistable transitions" when compared to 4-AP networks, this finding 76 may provide important insight into mechanisms underlying seizure generation. In turn, 77 these findings provide paramount in silico support for the integral role of inhibitory 78 interneurons in seizure initiation. 79 Materials and methods Neurons were modeled via a two dimensional system of ordinary differential equations 93 first described by Izhikevich [34]. This model has two variables: V , which represents the 94 membrane potential in mV; and u, which represents the slow "recovery" current in pA. 95 The model utilized here is slightly altered in the fashion described by Ferguson et. 96 al. [21], and is given by: In the above equations, C m represents the membrane capacitance in pF, v r represents 98 the resting membrane potential in mV, v t represents the instantaneous threshold 99 potential in mV, v peak is the spike cut-off value in mV, I syn is sum of all incoming 100 synaptic current to the neuron in pA (described in detail below), I app represents the 101 external applied current in pA (described in detail below), a is the recovery time

106
The use of Izhikevich model neurons was motivated by the goals of this study: 107 namely, here we do not strictly constrain our neuron model with experimental results, 108 but rather create a model that more "generally" matches the properties of an cortical interneurons. is also implied by the literature [32] and correspondingly influenced the determination of 131 these parameters. As the model of Ferguson et. al. [21] was used as a reasonable model 132 of a fast-firing inhibitory cell, the slope of the frequency-current (FI) of that neuron was 133 used for the control case. Except for the changes caused by a shifted rheobase and the 134 presence of adaptation, this slope was kept approximately the same for the 4-AP model. 135 With the different rheobases, this means that the firing frequency is larger in the 4-AP 136 model relative to control for a given input current.

137
The parameter values for both what will hereafter be referred to as the "control" 138 model and what will hereafter be referred to as the "4-AP" model are included in Curves are shown for frequencies calculated using the initial (blue and red) and final (cyan and magenta) inter-spike intervals to illustrate the tendency for spike-frequency adaptation (SFA). These comparisons show that the neuron models utilized in this study match the decreased rheobase and increased excitability and SFA of 4-AP treated neurons in comparison to control neurons (with the rheobases determined from unpublished in-house experiment for control and 4-AP neurons highlighted on the figure by the colored dots).

145
Similar to inhibitory network models developed by Ferguson et. al. [21], the neurons in 146 the networks modeled here were randomly connected by synapses utilizing a first-order 147 kinetic model. Each synapse is modeled by where g syn is the maximum inhibitory synaptic conductance in nS, s is the gating 149 variable, V is the membrane potential of the post-synaptic cell in mV, and E syn is the 150 inhibitory reversal potential in mV. As this value of E syn is set at an inhibitory value of 151 −75 mV for every possible synapse, this study includes only inhibitory synaptic 152 connections. Furthermore, g syn is uniform for each network studied, meaning each 153 connection in a given network has the same strength.

154
The gating variable models the proportion of open synaptic channels, with its 155 dynamics given by where α represents the inverse of the rise time constant and β represents the inverse of 157 the decay time constant [35].
[T ] models the concentration of neurotransmitter released 158 April 11, 2019 7/36 following a pre-synaptic action potential.
[T ] is represented as a unitary pulse lasting 1 159 ms, from the time of the pre-synaptic spike (t 0 ) to the end of the pulse (t 1 ). With this, 160 the dynamics of s can be simplified to the following two equations, where s ∞ = α α + β and τ s = 1 α + β .

162
Network model parameters 163 A fast rise time rate constant of α = 3.7037 ms −1 is used here as in [21]. Values for the 164 inhibitory reversal potential (-75 mV) and the synaptic decay rate constant (β = 0.3333 165 ms −1 ) were taken from Traub et. al. [36], and the range of inhibitory synaptic 166 conductances explored (0 to 10 nS) encompasses cortical estimates [36]. al. [37]. Given that a single large basket cell synapses onto approximately 23 other large 173 basket cells in this brain region [37], and assuming random connectivity, the probability 174 of connection between such cells would be at least 0.04. Considering networks with 175 connectivity densities lower than 0.04 would be unlikely to exhibit coherent dynamics in 176 a randomly connected inhibitory network [23,38], this density is used as a lower bound 177 for this investigation. Based on such estimations, this study utilized networks of 500 178 neurons with the neurons randomly connected with connection probabilities of 0.04, 179 0.08, 0.12 and 0.16. Connection probabilities much larger than this were not needed as 180 they would be clearly unrealistic relative to the biological estimations.

181
As done in Ferguson et. al. [21], cell heterogeneity in the networks was implemented 182 by varying the amplitude of the tonic external input current, I app , to each neuron. The 183 input currents were selected from a normal distribution with a mean value of I µ , with 184 the degree of heterogeneity in the input currents determined by the standard deviation, 185 σ. σ = 3, 6 and 12 pA were studied.

186
Simulations 187 The code underlying these simulations was written in the C programming language and 188 run on a Linux-based high-performance computing cluster utilizing Compute Canada did not trigger synaptic current until 100 milliseconds into the simulation (via a simple 193 manipulation in the code) to allow initial transients to decay.

194
In order to uncover other potential dynamical states of the network, a brief, large 195 amplitude current pulse was delivered uniformly to each cell in the network to perturb 196 the system and potentially bias it towards the synchronous dynamical state. This 2 ms 197 pulse had an amplitude of 1000 pA and was delivered at 1000 ms. This is analogous to 198 imposing homogeneous initial conditions causing instantaneous spiking of all neurons in 199 the network, in contrast to the randomized initial conditions that begin the simulations. 200 To identify networks that exhibited bistability between asynchronous and clustered The measure used to quantify coherent activity in the simulated networks, here termed 216 a Synchrony Measure, is a slight adaptation of a commonly used measure created by 217 Golomb and Rinzel [40,41] that quantifies the degree of spiking coincidence in the 218 network. This particular implementation of this measure has been utilized in previous 219 studies [23,42,43] 220 Briefly, the measure involved convolving a Gaussian function with the time of each 221 action potential for every cell to generate functions V i (t). The population averaged 222 voltage V (t) was then defined as and April 11, 2019 9/36 where < · > indicates time averaging over the interval for which the measure is taken. 227 The Synchrony Measure S was then defined as  A: An example network that exhibits "messy" synchronous dynamics both before and after the perturbation is delivered at 1000 ms, resulting in a moderate value of the Synchrony Measure in each case. Dynamics before the perturbation are shown in the left panel, while dynamics following the perturbation are shown in the right. Although the Synchrony Measure following the perturbation is larger than that before the perturbation, this increase does not indicate a bistable transition from asynchronous to synchronous dynamics, but rather qualitatively "tighter" synchrony. The choice of 0.3 as the "threshold value" in the articulation of the Bistability Measure prevents cases such as this from contributing positively to the measure. B: An example network exhibiting asynchrony before the perturbation (left panel) and very clear synchrony afterwards (right panel), along with the corresponding Synchrony Measures. Very low synchrony measures (typically less than 0.25) indicate asynchrony, while higher Synchrony Measures illustrate synchrony, with higher values indicating more ordered and less "messy" synchrony.
This research was interested not merely in the degree of synchronous firing in the  This measure was calculated in three steps. First, the difference between the 243 Synchrony Measure following the modeled perturbation and the Synchrony Measure 244 from randomized initial conditions was taken for each network in some context (i.e. for 245 a particular parameter range encapsulated by a given heatmap). Second, only cases 246 when this difference was > 0.3 were included in the summation (this choice is justified 247 in detail below). Finally, the sum of the Synchrony Measure differences exceeding this 248 threshold value was taken to yield the final Bistability Measure.

249
The choice of the "threshold" value of 0.3 in the second step above merits further 250 explanation. The Synchrony Measure is not a binary differentiation between 251 asynchronous and synchronous dynamics, but rather a quantitative measure of the 252 degree of synchronous firing. This means that increases in the Synchrony Measure, 253 particularly subtle ones, do not necessarily indicate a differentiation of asynchrony from 254 synchrony, but instead could indicate the presence of qualitatively "tighter" synchrony. 255 An example of such a case is seen in Fig 3A. However, large increases in the Synchrony 256 Measure are almost always indicative of entirely different dynamical states, as shown by 257 the example in Fig 3B). After a thorough investigation of the correspondence between a 258 qualitative assessment of synchrony (i.e. visual inspection of raster plots) and the  Ornstein-Uhlenbeck process

282
The perturbation described above is motivated primarily by the desire to uncover a 283 mechanism for the transition from asynchrony to synchrony from the perspective of 284 dynamical systems. In order to assess whether this mechanism is biologically reasonable, 285 an analogous perturbation that might arise in more biologically grounded models was 286 sought.

287
An Ornstein-Uhlenbeck process [44] is used in the literature to model background 288 synaptic input into a network [45,46], and is used in this study to determine whether 289 "perturbation-like" activity might arise naturally from this model of external synaptic 290 input. This process, used to determine the conductance of excitatory synaptic input in 291 this context, is described mathematically by the following equations [45] with an initial 292 condition g e (0) = g e0 : where N (0, 1) is a normal random number taken from a distribution with 0 mean and a 295 standard deviation of 1.

296
The insights from Piwkowska et. al. [46] allowed for the choice of parameters 297 constrained by cortical data. The parameters used in the Ornstein-Uhlenbeck process 298 utilized in this study were g e0 = 3 nS, τ e = 2 ms, and D e = 2 (a unitless diffusion 299 coefficient), and the integration time step was h = 0.01 ms.

300
Code Accessibility

301
The code/software described in the paper is freely available online at 302 https://github.com/FKSkinnerLab/CorticalInhibitoryNetwork.

304
In order for inhibitory interneurons to serve as the impetus for seizures, as posited by a 305 "GABAergic initiation hypothesis", there should exist a mechanism that explains both 306 how a network of such interneurons can suddenly transition into synchronous dynamics 307 (which in turn provides a large burst of synaptic inhibition to pyramidal neurons in the 308 behaving animal) and why model seizure networks are more prone to exhibit this 309 transition and the resultant pathological brain activity. To examine this, an 310 experimental 4-AP hyper-excitable state, as used to model seizures by Chang et. 311 al. [12], is used as a basis for these mathematical model explorations, with results from 312 such networks compared to model control networks.

313
Bistability between coherent and incoherent states is exhibited 314 more robustly by 4-AP inhibitory networks 315 The concept of bistability arises primarily from the mathematical study of dynamical 316 systems. In this context, a "stable" state is one which will be preserved by the system 317 for all time in the absence of any perturbations to the conditions defining the system.

318
In nonlinear systems it is possible for multiple stable states to exist, and for the network 319 to naturally settle into any one of these stable states depending upon a variety of 320 factors including the initial conditions and any perturbations that might be delivered.

321
In a biological system, this could manifest from the history of inputs from different 322 brain structures along various pathways to the network in question. Such a system is 323 defined to be "bistable" or "multistable" given the existence of more than one stable 324 solution to the mathematical equations [47].

325
The results presented in Fig 4 show that many of the networks within the parameter 326 regime considered in this work exhibit bistability. In panels A and B the Synchrony

335
There appear to be a number of networks in both the control and 4-AP settings that 336 show a high Synchrony Measure, and thus coherent network states, following the 337 perturbation but not from randomized initial conditions. This is indicative of a  following traumatic brain injury [48,49] underlying seizure onset, which in this context is the transition from asynchronous to 360 synchronous firing via a "bistable transition". It is also interesting to note that the 361 bistable regime is both wider (i.e. encompassing a larger range of synaptic strengths) 362 and includes lower driving currents for 4-AP networks, although the latter is perhaps In inhibitory network models of CA1 hippocampus that were constrained in size, 406 connection probability, cellular and synaptic properties, previous work by [21] 407 demonstrated "sharp transitions" between asynchronous and synchronous firing caused 408 April 11, 2019  Fig 4, where the large, brief current pulse used as the perturbation throughout this study is replaced by a current pulse informed by the Ornstein-Uhlenbeck Process shown in panel A that represents in vivo-like activity. Despite this change, which amounts to a wider pulse with significantly lower amplitude, the control and 4-AP networks still exhibit antithetical responses to this perturbation; namely, control networks return to asynchronous firing following the perturbation while 4-AP networks transition into synchronous dynamics.
by a small, permanent increase in the external drive to the network. This additional 409 mechanism has both experimental and computational support (see the discussion in the 410 Introduction) to explain a transition into synchrony in purely inhibitory networks. We 411 thus investigated it as a potential mechanism for the initial inhibitory synchrony 412 necessary for a "GABAergic initiation hypothesis" of seizure. Indeed, potentially 413 eliminating it as a candidate mechanism in this context would provide additional 414 support for the viability of the "bistable transition" mechanism described in the 415 previous section. shown jointly in Table 2.  Cortically motivated inhibitory networks exhibit a "sharp transition" between asynchronous and synchronous dynamics driven by an increase in the external drive for various connection probabilities. A-C: Heatmaps displaying the Synchrony Measure for control networks (left) and 4-AP networks (right) with a standard deviation amongst the driving currents of 6 pA and varying connectivity densities. In these heatmaps, the inhibitory synaptic weight is varied along the x-axis and the average external applied current is varied along the y-axis, and the measure is taken from random initial conditions. D-F: Two dimensional "slices" of the heatmaps in panels A-C taken to better illustrate the sharpness of the transition between asynchronous and synchronous dynamics as well as more easily compare this sharpness both across varying connection probabilities and between control and 4-AP conditions. Panel D shows results for an inhibitory synaptic weight of 1.25 nS, panel E for an inhibitory synaptic weight of 2.0 nS, and panel F for an inhibitory synaptic weight of 2.5 nS. There is no significant difference in the tendency for 4-AP versus control networks to exhibit the "sharp transition" from asynchrony to synchrony despite differences in the parameter regime supporting synchrony. Furthermore, the differences in the "sharpness" of the transition in the two cases are not robust. Varying the heterogeneity in external driving current in modeled purely inhibitory networks largely preserves the general dynamical differences and similarities seen between the 4-AP and control cases from randomized initial conditions. A-C: Heatmaps displaying the Synchrony Measure for control networks (left) and 4-AP networks (right) with a connection probability of 0.12 and varying standard deviations amongst the driving currents. In these heatmaps, the inhibitory synaptic weight is varied along the x-axis and the average external applied current is varied along the y-axis, and the measure is taken from random initial conditions. D-F: Two dimensional "slices" of the heatmaps in panels A-C taken to better illustrate the sharpness of the transition between asynchronous and synchronous dynamics as well as more easily compare this sharpness both across varying connection probabilities and between control and 4-AP conditions. Panel D shows results for an inhibitory synaptic weight of 1.25 nS, panel E for an inhibitory synaptic weight of 2.0 nS, and panel F for an inhibitory synaptic weight of 2.5 nS. Once again, there is no significant difference in the tendency for 4-AP versus control networks to exhibit the transition from asynchrony to synchrony, nor any significant differences in the "sharpness" of this transition.  Measure going up a column in the heatmaps or looking at the two-dimensional traces), 473 with only the applied current value at which the transition occurs changing.

474
Furthermore, any comparisons of the relative "sharpness" of the transitions in 475 control and 4-AP networks are qualitative at best. In a majority of the comparisons 476 illustrated in the two-dimensional plots (Fig 7D-F), 4-AP networks displayed a higher 477 slope measure shown in Table 2 than their control counterparts. However, this feature is 478 not entirely robust (see, for example, the comparison of networks with a connection 479 probability of 0.08 in Fig 7E and Table 2), and there is no guarantee that the minor al. [21] in the hippocampus. However, there is no robust difference in the tendency for 496 4-AP versus control networks to exhibit this transition. This strongly suggests that this 497 mechanism is unlikely responsible for the initial step in a "GABAergic initiation 498 hypothesis". Instead, a mechanism driven by a "bistable transition" is more plausible in 499 the context of explaining GABAergic seizure initiation, as it is much more likely to 500 occur in 4-AP networks rather than control networks. While transitions into synchrony 501 caused by minor, permanent increases to the external drive to an inhibitory network 502 certainly could occur in the brain given the existing literature, this conclusion implies 503 that mechanism is more likely to drive non-pathological oscillations rather than the 504 pathological inhibitory synchrony potentially initiating seizure.

506
Computational models of epilepsy encompass various levels of detail and address 507 different seizure aspects [50]. To help individuals with epilepsy, mechanisms underlying 508 the initiation, propagation and termination of seizures must be discovered and fully 509 understood (see Fig 1). This paper aims to provide support for a novel hypothesis for 510 seizure initiation, a "GABAergic initiation hypothesis", by proposing a viable 511 cellular-based mechanism by which the necessary initial step of this hypothesis, the 512 sudden transition of ictogenic inhibitory networks into synchrony, might come about.

513
Such a mechanism has not yet been presented in the literature of which the authors are 514 aware. This hypothesis proposes that synchronous activation of inhibitory interneurons 515 (the "initial step" investigated here) gives rise to a strong inhibitory signal that 516 activates the excitatory, pyramidal cell population via PIR [12] (see Fig 1). As 517 experimental support for this theory accumulates [1,[8][9][10][11][12][13], computational insights such 518 as those presented here can provide further support for its viability by proposing synchrony and predominantly GABAergic IIS [29,51], and thus 4-AP is a commonly 526 used epilepsy model [31,32]. This experimental evidence justifies not only the use of has also been shown to underlie IIS in the in vivo pilocarpine model of epilepsy [13], 531 precede seizures in both the low-Mg, high K+ model [26,28] and electrical stimulation 532 models of seizure initiation [52], and more generally precede seizures in rodents [9,29]. 533 Thus, the insights gained from this study can be applied to the general study of 534 epileptiform activity.

535
Two potential mechanisms by which inhibitory networks could suddenly transition 536 into synchronous firing were examined in this paper. In the context of the study of 537 seizure initiation, the mere existence of such a transition is not of primary concern; 538 rather, such a transition should occur appreciably more often in hyper-excitable (i.e. were compared. In contrast, transitions caused by a brief perturbation to external drive 547 to the system, termed "bistable transitions" given the correspondence of this transition 548 with the concept of bistability from dynamical systems theory, were notably more likely 549 to occur in 4-AP than control networks. This crucial difference implies that "bistable 550 transitions" are a more viable candidate mechanism that explains a sudden transition of 551 an inhibitory network into synchrony in pathological networks.

552
The general concept of bistability has been discussed previously in epilepsy 553 literature, given that epilepsy as a disease represents the sudden transition between two 554 seemingly stable brain states: the "healthy" non-seizure state characterized by largely 555 uncorrelated neural activity and the "pathological" seizure state characterized by 556 synchronous neural firing [53]. However, it is worth emphasizing that the setting in 557 which bistability is analyzed, and thus the context of the corresponding perturbation, is 558 unique in this study, driven primarily by the focus on a "GABAergic initiation 559 hypothesis". Indeed, existing studies investigate a bistability between seizure and 560 non-seizure states in settings such as intact hippocampal slices [54], a computational 561 network of both excitatory and inhibitory cells with special emphasis on the role of 562 extracellular potassium concentrations [55], or more general mathematical settings [53]. 563 In contrast, in this study bistability is analyzed solely in an inhibitory network, and the 564 bistability does not in itself represent the transition into seizure, but rather a dynamical 565 change that might precipitate seizure onset due to its downstream effects (as illustrated 566 by a "GABAergic initiation hypothesis" schematized in Fig 1). 567 We also highlight an important distinction between this work and other 568 computational work investigating the role of GABAergic signalling in epileptiform 569 activity and inter-ictal discharges (IID): while recent literature investigating this topic 570 makes use of the potential depolarizing capacity of GABA [56,57], the work presented 571 here uses purely inhibitory GABAergic synapses. Indeed, while changes in the GABA 572 reversal potential are seen during seizure propagation [3], the changes in chloride 573 concentrations necessary to elicit this feature are unlikely to exist prior to or during 574 seizure initiation, which is the focus of this research. Moreover, the mechanisms 575 proposed in the work of Chizhov et. al. [56,57] do not focus on the capacity of excitatory 576 cells for PIR, in contrast to the "GABAergic initiation hypothesis" discussed here.

577
The mechanism proposed in this paper is the first of which the authors are aware 578 that describes how the initial step of a "GABAergic initiation hypothesis" might occur 579 with both biological [12] and computational (this study) support. This, in turn, provides 580 new and convincing evidence that may help to explain how the hyper-excitability 581 induced by 4-AP causes the cortex to be more vulnerable to seizures, and more generally 582 how interneurons can be involved in the initiation of cortical seizures clinically [58,59]. 583

584
The exploration of a transition driven by a small, permanent increase to the external 585 drive was motivated by modeling studies [21] and physiological evidence [24] of is very probable that the hippocampal networks would exhibit bistability of some form. 596 However, of critical importance in the context of this study is the lack of an appreciable 597 difference in the tendency for 4-AP and control model networks to exhibit this 598 transition.

599
The mechanism articulated in this paper making use of the mathematical concept of 600 bistability addresses the shortcomings, in the context of seizure initiation, of the indicates that this transition is potentially viable in a biologically-grounded setting as 613 well.

614
The analysis of this "bistable transition" reveals that it is more likely to occur in 615 4-AP networks as opposed to their control counterparts. This result indicates that it is 616 a much more likely culprit in the initial step of a "GABAergic initiation hypothesis" of 617 seizure than a transition brought about by a small, permanent increase in external drive 618 to the network. Further support for bistability serving a role in seizure onset is found in 619 the statistics of inter-seizure intervals [60]. Moreover, the cortically-motivated differences seen in their FI curves presented in Fig 2). This finding is fairly robust over 639 all but the weakest inhibitory syanptic weights. Thus, the critical implication is that the 640 mechanism involved in the "bistable transition" involves an interplay of cellular 641 (potentially not only the hyper-excitability, but also the increased adaptation, in 4-AP 642 neurons) and network properties, and could not be replicated merely by causing the 643 neurons to fire faster in some artificial fashion.

644
Also of interest in this analysis was that both control and 4-AP networks show a 645 similar increase in average firing rate when transitioning from asynchrony to synchrony 646 (highlighted by the example raster plots and corresponding mean firing frequencies 647 presented in Fig 9). This result is analogous to a similar finding in [21] and indicates a 648 potential avenue for an experimental exploration of the results presented here: namely, 649 a substantial increase in firing rate following a perturbation to an inhibitory network is 650 likely indicative of the transition into synchronous firing. Such behavior is likely more 651 However, networks of fast-firing interneurons can also produce slow population 676 output as shown in modeling studies [63]. The ability of fast-firing inhibitory networks 677 to produce slow population activities was shown to be possible via individual cells 678 having enough of a "kink" in their frequency-current (FI) curves that allowed a bistable 679 network mechanism to be present [63]. The modeled slow population activity (< 5 Hz) 680 is seen in vitro using a hippocampal preparation [64,65], and a bistable network 681 mechanism was subsequently leveraged to explain paradoxical changes seen in Rett 682 syndrome mice from the perspective of these same slow population activities [66].

683
A critical difference between the bistable network mechanism of Ho et. al. [63] and 684 bistability related to properties of the ING mechanism (analogous to that presented 685 here) was summarized by Skinner and Chatzikalymniou [67]. In the work of Ho et.

686
al. [63], the mean excitatory drive received by inhibitory cells in the network must be 687 close to their spiking rheobase. The bistability is between states with low or high 688 numbers of fast-firing cells, and this allows slow population activities to come about due 689 to excitatory fluctuations in the system. A similar mechanism could be in play in the 690 work of Schlingloff et. al. [24] where an in vitro representation of sharp waves was 691 examined and it was suggested that sharp waves could be generated stochastically from 692 excitatory input. In contrast, for an ING-related bistability, the excitatory drive to the 693 inhibitory cells is not close to spiking rheobase, but as shown by Rich et. al. [23] and in 694 the cortically-motivated networks presented here, bistability between synchronized high 695 computational studies of inhibitory synchrony [69][70][71][72]. As such, the coherent dynamics 703 seen in our cortically-motivated inhibitory networks correspond with insights from these 704 more abstract computational studies. This literature contributed to the articulation of 705 the ING mechanism [15,16,18,19] that is most likely driving the coherent dynamics seen 706 in these cortical inhibitory networks. Another seminal study on inhibitory synchrony 707 and ING found that the synchrony promoted by the ING mechanism is most robust 708 when networks are more densely connected and cellular heterogeneity is low [73], 709 features replicated in the cortically-motivated networks presented here.

710
In this context we note that computational studies proposing mechanisms for 711 synchronous network oscillations are typically concerned either with purely inhibitory 712 networks (as presented here), purely excitatory networks [74], or networks containing 713 inter-and intra-connected subnetworks of excitatory and inhibitory cells (E-I networks). 714 Crucially, the mechanisms underlying synchrony and their dependence on features such 715 as cell excitability properties (i.e. the Type I vs. Type II distinction [68]), external 716 drive to the network, and network connectivity can vary significantly depending on the 717 type of network studied. For example, the results of Hansel et. al. [74] imply that an 718 excitatory network made up of cells of the type studied here is highly unlikely to ever 719 synchronize. Similarly, while the Pyramidal Interneuron Network Gamma (PING) 720 mechanism is commonly cited as a mechanism causing synchronous oscillations in E-I 721 networks [16,72,75,76], recent work has revealed that the predictions of this mechanism 722 are altered by varying individual cellular properties and network connectivity [42,43]. consideration. Indeed, one could consider designing a neuromodulation study using the 730 Blue Brain Project [37] to examine this given the insights gleaned from this study.

731
The network structure used in this work, a purely inhibitory network, is also 732 simplified from the biology. However, this choice was critical in allowing for the 733 articulation of a mechanism of action underlying the transition of interneurons into 734 synchrony. Such a mechanism is a paramount and necessary "first step" towards an 735 overarching mechanism of a "GABAergic initiation hypothesis" and provides initial

740
For the work here, we focused on differences between control and 4-AP neurons as 741 encapsulated in our models. It is unlikely that utilizing a more realistic noisy synaptic 742 input would affect the primary results of this work, since both noisy [77] and 743 deterministic [21] inputs were used in previous hippocampal inhibitory network models 744 without changing insights regarding the transition into synchrony.

745
While use of a simplified neuron model and network structure enables extensive 746 parameter explorations to be easily done and dynamical aspects, like bistability, to be 747 uncovered, parameter interpretation relative to details of the biological system is less 748 straightforward. However, studies such as this could help leverage understanding and 749 motivate hypothesis-driven explorations in more detailed models. We note that due to 750 the relatively sparse connectivity of the cortically-motivated inhibitory networks studied, 751 mathematical tools such as reduction to phase oscillator models as in Hansel et. al. [74] 752 that require an assumption of all-to-all connectivity and weak coupling cannot be easily 753 applied to do further theoretical analyses.

754
With the ability to obtain and model human cellular data [78,79], it will be 755 interesting to consider how one might examine and take advantage of these mechanistic 756 insights. Expanding our mechanistic insights here to include excitatory cell networks 757 and other seizure phases (see Fig 1 [5]) are exciting considerations.

759
The hypothesis that excessive inhibitory signalling serves a causal role in seizure 760 initiation has support from numerous recent experimental studies, but until now has not 761 subject to rigorous computational investigation. Here, we provide the first 762 computational support for this theory by articulating a mechanism that not only 763 potentially explains the necessary, initial step in this process (the sudden transition of 764 inhibitory interneurons from asynchronous to synchronous firing), but also why 765 networks in a model ictogenic state are more vulnerable to this transition. This latter 766 feature distinguishes this mechanism as one particularly important in the study of 767 epilepsy, as any process asserted to relate to seizure initiation should be more likely to 768 occur in an ictogenic as opposed to healthy brain. As this novel hypothesis of seizure 769 initiation is inherently counter-intuitive, especially in comparison to the more 770 historically common theory of over-excitation or dis-inhibition initiating seizure, this 771 computational work articulating potential mechanisms of action is especially important 772 to support the theory's viability. Indeed, our work presented here not only provides the 773 first evidence of which we are aware that feasible mechanisms underlying necessary 774 steps in this process exist, providing paramount support for a "GABAergic initiation