GABA Regulation of Burst Firing in Hippocampal Astrocyte Neural Circuit: A Biophysical Model

It is now widely accepted that glia cells and gamma-aminobutyric acidergic (GABA) interneurons dynamically regulate synaptic transmission and neuronal activity in time and space. This paper presents a biophysical model that captures the interaction between an astrocyte cell, a GABA interneuron and pre/postsynaptic neurons. Specifically, GABA released from a GABA interneuron triggers in astrocytes the release of calcium (Ca2+) from the endoplasmic reticulum via the inositol 1, 4, 5-trisphosphate (IP3) pathway. This results in gliotransmission which elevates the presynaptic transmission probability rate (PR) causing weight potentiation and a gradual increase in postsynaptic neuronal firing, that eventually stabilizes. However, by capturing the complex interactions between IP3, generated from both GABA and the 2-arachidonyl glycerol (2-AG) pathway, and PR, this paper shows that this interaction not only gives rise to an initial weight potentiation phase but also this phase is followed by postsynaptic bursting behavior. Moreover, the model will show that there is a presynaptic frequency range over which burst firing can occur. The proposed model offers a novel cellular level mechanism that may underpin both seizure-like activity and neuronal synchrony across different brain regions.


INTRODUCTION
Spiking neural networks (SNNs) are considered to be the most biologically plausible representation of brain function (Ghosh-dastidar and Adeli, 2009). Additionally, SNNs capture a Hebbian type learning paradigm where the timing between pre-and post-synaptic spikes dictates whether synaptic depression or potentiation occurs (Song et al., 2000). SNNs have also been shown to be effective in time series prediction (Reid et al., 2014), spatiotemporal pattern recognition (Hu et al., 2013), and system control  in various application domains. In SNNs, the neurons and synapses are fundamental components in the network, where the information is encoded in spikes or action potentials for transmission between neurons (Izhikevich, 2003). In the central nervous system neurons receive input stimuli and respond by firing spike patterns such as bursting, which has been observed in the hippocampus of rodents (Miles and Wong, 1986), electric fish (Gabbiani et al., 1996), and in the primary motor cortex, brainstem and thalamus within the somatomotor system of humans (Arichi et al., 2017). The bursts can, in some cases, represent normal brain function and in other cases abnormal brain function (e.g., epilepsy) (Araque et al., 1999;Halassa et al., 2007).
Research has shown that astrocytes, one type of glial cell, modulate neuronal activity (Halassa et al., 2007;Breslin et al., 2018;Flanagan et al., 2018) where a single astrocyte may enwrap a large number of synapses (∼ 10 5 synapses), and connect to several neighboring neurons (four-eight). The interplay between an astrocyte and the neighboring neurons is believed to occur at the tripartite synapse (Araque et al., 1999), which is bi-directional and serves, in some cases, to modulate the synaptic transmission probability rate (PR): via the direct/indirect retrograde signaling messenger endocannabinoids (Wade et al., 2012). This gives rise to re-modeling of the SNN connectivity (Wade et al., 2011;Naeem et al., 2015;Johnson et al., 2018;Liu et al., 2018).
It has also been reported that gamma-aminobutyric acidergic (GABA) interneurons participate in astrocyte-mediated control of excitatory synaptic transmission (Perea et al., 2016) and exercises control over the firing frequency of pyramidal cells. Furthermore GABA release synchronizes principal cell population discharge contributing to the generation of rhythmic activity in neuronal networks, such as theta and gamma frequency oscillations (Kullmann, 2011). A recent paper reported that GABA released in proximity to a tripartite synapse can activate GABA-B receptors on the astrocyte leading to gliotransmission, which is known to regulate synaptic transmission probability . The research reported in Kurosinski and Götz (2002), Kullmann (2011), Liu et al. (2018) provides the underpinning for the work presented here.
In this paper, we investigate the coupling between a GABA interneuron, an astrocyte terminal and the pre and postsynaptic terminals. The main contributions of this paper include (i) a novel biophysical model that describes the signaling pathways at the tripartite synapse and (ii) a novel mechanism that can potentially explain postsynaptic neuron burst firing. The rest of paper is organized as follows. Section 2 presents the biophysical model while section 3 provides simulation results that demonstrate the bursting. Section 4 concludes the paper and discusses future work.

BIOPHYSICAL MODEL OF A NETWORK BURSTING
In this section, a detailed discussion of the signaling pathways at the tripartite synapse is presented with a specific focus on GABA signaling between the presynaptic terminal and the nearby astrocyte. It will be shown that this interplay acts as a frequency dependent switch, which modulates the probability of release (PR) at the presynaptic terminal. Our Ca 2+ dynamics model shows that calcium (Ca 2+ ) oscillations only occur over a range of inositol 1, 4, 5-trisphosphate (IP 3 ) concentrations and furthermore this paper will show that Ca 2+ oscillations are periodic and this behavior is key to the bursting behavior.

Signaling Pathways and Activity Regulations
The conventional tripartite synapse has three terminals: the presynaptic axon, postsynaptic dendrite and the astrocyte cell (Wade et al., 2012;Liu et al., 2018). In this paper we consider earlier work where, in a hippocampal astrocyte neural network, GABA interneurons interact with excitatory tripartite synapses to dynamically change the synaptic transmission behavior from inhibitory to excitatory through modulation of PR (Perea et al., 2016). The signaling pathways between the GABA interneuron and tripartite synapse are shown in Figure 1. When an input stimulus of frequency (f pre ) is present at the excitatory presynaptic axon, neurotransmitter (glutamate) is released into the cleft and subsequently binds to receptors at the postsynaptic dendrite causing the depolarization of the postsynaptic neuron. While the authors accept that fast-spiking interneurons can fire at much higher frequencies than glutamatergic neurons, in this work we assume for simplicity that the firing rate of the GABA interneuron (f GABA ) follows f pre , as the most likely physiological condition would be the activation of GABA interneuron by activation of glutamatergic axons (Serrano et al., 2006;Covelo and Araque, 2018). While GABA initially binds to GABA-A receptors inhibiting synaptic transmission and post-synaptic neuronal activity, recent work (Perea et al., 2016) has shown that with repeated firing of GABA interneurons, GABA also binds to GABA-B receptors on the astrocyte membrane, resulting in a switch from inhibition to excitation at the presynaptic terminal, and an associated excitatory response at the postsynaptic terminal. In this paper we focus on astrocyte-mediated GABAinduced excitation since the postsynaptic inhibition was found negligible, and the transient acute presynaptic inhibition was overpowered by the astrocyte signaling during sustained activity (Perea et al., 2016). Hence, our model dos not incorporate these negligible or transient inhibitory effects, focusing in the sustained mechanisms and effects of inhibitory signaling through astrocyte activation.
As f GABA increases, the GABA concentration level in the extracellular space increases, and a level is reached whereby binding to GABA-B receptors on the astrocyte membrane commences, leading to the production of IP 3 : we subsequently refer to IP 3 due to GABA as IP GABA 3 , which contributes to the overall cytosolic IP 3 (Perea et al., 2016). IP GABA 3 is a secondary messenger which is degraded when released into the cytoplasm: initially cytosolic Ca 2+ and IP 3 levels are low and therefore degradation of IP 3 will also be low, as the degradation rate correlates with both Ca 2+ and IP 3 concentrations. This degradation is gradually overcome with increasing levels of GABA and the PLCδ signaling pathway, which is modulated by Ca 2+ and is accounted for in this work. Finally IP GABA 3 starts to bind to IP 3 receptors (IP 3 R s ) on the Endoplasmic Reticulum (ER). When the total cytosolic IP 3 is sufficiently high, Ca 2+ is released from the ER (De Pittà et al., 2009). At some point both IP 3 and Ca 2+ reaches a level at which an oscillating Calcium-Induced Calcium Release (CICR) occurs from the ER (Marchant et al., 1999): hereafter referred to as T CICR . Several mechanisms are believed to contribute to Ca 2+ oscillation but there is still much debate around this topic. For example, IP 3 R s have binding sites for both IP 3 and Ca 2+ , and Ca 2+ release from the ER is believed to rely on coincidence binding of these ions. The time between IP 3 and Ca 2+ binding depends on the concentration of these ions and therefore this could explain why Ca 2+ is believed to be a regulator of IP 3 R s activity: at low Ca 2+ levels IP 3 R s activity is increased, whereas the opposite is true at high Ca 2+ levels. This is in agreement with other research (Dawson, 1997). However, there is not enough experimental evidence on these receptors to formulate a sufficiently detailed model. Therefore, in this work we revert to a hitherto accepted model (Perea et al., 2016) where Ca 2+ oscillatory behavior is believed to arise from the feedback interplay between Ca 2+ , IP 3 , and IP 3 degradation. As Ca 2+ and IP 3 rapidly increase there is a complex dependency between the concentrations of both Ca 2+ and IP 3 and Ca 2+ /IP 3induced degradation of IP 3 , which is the dominant process at elevated Ca 2+ /IP 3 levels. Therefore, a transient elevation of Ca 2+ and/or IP 3 is followed by a rapid drop in IP 3 , which can reduce IP 3 to below T CICR . At this point degradation of IP 3 is weak because both the Ca 2+ and IP 3 levels have fallen and therefore IP 3 starts to increase again due to IP GABA 3 . When the T CICR level is reached again a transient elevation of Ca 2+ re-occurs. We will demonstrate that our results support this behavior. This oscillatory behavior causes the release of the glutamate from the astrocyte (gliotransmitter) into the synapse [see Glu (e-SP) pathway in Figure 1], which binds to pre-synaptic group I metabotropic Glutamate Receptors (mGluRs) at the presynaptic terminal. This signaling pathway results in an increase in PR at the presynaptic terminal (Navarrete and Araque, 2010).
As PR increases more glutamate is released into the cleft and potentiation/depression of the synaptic weight can commence with the availability of glutamate to bind to Nmethyl-D-aspartate (NMDA)-type glutamate receptors (Lüscher and Malenka, 2012). While we acknowledge that the biophysical mechanisms regulating the functional dependency between PR and plasticity are complex and not fully understood, we propose that PR acts as a "switch" which can turn on/off potentiation/depression at synaptic sites. To formulate a tractable mathematical model that captures this relationship we modulate the height of the Spike Timing Dependent Plasticity (STDP) associated plasticity window using PR: with PR ≥ PR * (PR * is defined as the plasticity activation level) the plasticity window fully opens and with PR < PR * the plasticity window closes. The decision on whether potentiation or depression occurs is governed by the STDP rule (Magee and Johnston, 1997) where potentiation occurs when the presynaptic spike precedes postsynaptic spike, otherwise depression occurs. Additionally, we consider the case where the postsynaptic neuron is sufficiently depolarized such that the retrograde messenger 2-arachidonyl glycerol (2-AG) is released from the postsynaptic neuron. Since the contribution of 2-AG signaling to the observed GABA-mediated regulatory effects of astrocytes on excitatory transmission is negligible (Perea et al., 2016), the authors take the view that 2-AG signaling onto GABAergic terminals would not be a significant factor in network bursting. However, we do consider 2-AG binding to type 1 Cannabinoid Receptors (CB1Rs) on the astrocyte membrane which then initiates the release of the IP 3 into the cytoplasm of the astrocyte: we denote this secondary messenger as IP AG 3 . During the synapse learning phase, the frequency of the postsynaptic neuron, f post , is increasing, as is the 2-AG signal and consequently IP AG 3 . As IP AG 3 contributes to the total IP 3 , IP 3 will eventually reach a level where degradation of IP 3 no longer reduces IP 3 (and therefore Ca 2+ ) to below T CICR . In this instance, both the oscillatory Ca 2+ transient and Glu (e-SP) pathways cease ) (see Figure 2). In addition, the released 2-AG also binds to CB1Rs on the presynaptic terminal triggering the Suppression of Excitation (DSE) pathway, and results in a decrease in PR (Alger, 2002). Due to the reduction in the Glu (e-SP) pathway and the increase in the DSE pathway, PR decreases at the presynaptic terminal, the level of neurotransmitter in the cleft then falls to baseline and the frequency of the postsynaptic f post diminishes which in turn causes IP AG 3 to reduce. Furthermore, the total IP 3 degrades due to cytosolic degradation pathways including IP 3 3-kinase IP 3K 3 , and dephosphorylation by inositol polyphosphate 5-phosphatase (IP 5P 3 ) (see Equation 11). Together, these processes reduce IP 3 levels below T CICR , and the rate of degradation diminishes sufficiently to allow IP 3 to increase again due to IP GABA 3 . When the T CICR level is again exceeded, Ca 2+ oscillations re-commence, the Glu (e-SP) pathway is reestablished, PR increases and the level of neurotransmitter in the cleft is raised. In this post-learning phase PR cannot be elevated to a level where the plasticity window opens (PR < PR * ), as the postsynaptic neuron is active and therefore the 2-AG pathway leads to a reduction in PR due to the DSE pathway. Consequently, the postsynaptic neuron firing rate reaches a maximum when the T CICR level is reached but it subsequently falls afterwards: a postsynaptic burst has occurred. This is followed by repeated bursts at each Ca 2+ oscillatory period. We therefore propose that neuronal burst firing directly correlates with astrocytic Ca 2+ oscillation.
Moreover, it should be noted that the duration of the burst correlated with the frequency of presynaptic terminal f pre at the excitatory presynaptic axon. As f pre increases so does the rate of increase of IP 3 and the burst period is reduced. Therefore, network bursting is f pre dependant and will only occur over a range of f pre .

Postsynaptic Neuron Model
In this paper, the Leaky Integrate and Fire (LIF) model (Gerstner and Kistler, 2002) is used due to the relatively low computing requirement and minimal parameters tuning. The LIF model is given by where τ m is the neuron membrane time constant, v is the neuron membrane potential, R m is the membrane resistance, I i syn is the current injected to the neuron membrane by ith synapse, and n is the total number of synapses associated with the neuron. When the neuron membrane potential v is greater than the firing threshold value, v th , the neuron fires and outputs a spike followed by a reset state or a refractory period (∼ 2ms). The release of 2-AG correlates with the postsynaptic neuron activity (Naeem et al., 2015) and this is expressed as where AG denotes the released amount of 2-AG, τ AG , and r AG are the 2-AG decay and production rates, and t sp is the postsynaptic spike time. The released 2-AG binds to the CB1Rs at the presynaptic terminal and at the astrocyte terminal, and this will be discussed in section 2.4.

GABA Interneuron
The spike train at the presynaptic axon also presents at the GABA interneuron causing the release of GABA neurotransmitter (Perea et al., 2016), which can be described by where GABA denotes the released amount of the neurotransmitter GABA, τ GABA , and r GABA are the GABA decay and production rates, and t sp is the presynaptic spike arrival time. GABA binds to the GABA-B receptors at the astrocyte cell and this is modeled in the next subsection.

Astrocyte Cell
When GABA binds to GABA-B receptors on the astrocyte membrane, the amount of IP 3 released is given by where IP GABA 3 is the quantity of IP 3 generated by GABA within the cytoplasm, IP GABA * 3 is the baseline GABA level, τ GABA ip3 is the decay rate of IP GABA 3 and r GABA ip3 is the production rate of IP GABA 3 .
When the postsynaptic neuron fires, the released 2-AG can also trigger IP 3 generation (Wade et al., 2012) and this is modeled by where IP AG 3 is the quantity of IP 3 generated by 2-AG within the cytoplasm, IP AG * 3 is the baseline level, τ AG ip3 is decay rate of IP AG 3 , and r AG ip3 is production rate of IP AG 3 . In addition, the IP 3 production is also increased by the hydrolysis of the highly phosphorylated membrane lipid phosphatidylinositol 4, 5-bisphosphate (PIP 2 ), such as the phosphoinositide-specific phospholipase C (PLC) isoenzyme of PLCδ (De Pittà et al., 2009). The PLCδ signaling is agonist independent and modulated by Ca 2+ (De Pittà et al., 2009), and its activation rate can be modeled by where the maximum PLCδ-dependant IP 3 production rate (De Pittà et al., 2009) can be modeled by and K δ is the inhibition constant of PLCδ activity. The Hill function (De Pittà et al., 2009) is described by where n is the Hill coefficient and K is the midpoint of the Hill function, namely the value of x at which Hill(x, K, n)| x=K = 1/2. The degradation of IP 3 mainly occurs through phosphorylation into inositol 1, 3, 4, 5-tetrakisphosphate (IP 4 ), catalyzed by IP 3 3-kinase (3K), and dephosphorylation by inositol polyphosphate 5-phosphatase (5P). The rate of IP 3 degradation by IP 5P 3 (De Pittà et al., 2009) can be modeled by where r 5P is the IP 3 degradation rate by IP-5P. The activity of IP 3K 3 is regulated by Ca 2+ in a complex fashion (De Pittà et al., 2009). The rate of IP 3 degradation by IP 3K 3 can be modeled by where v 3K is the maximum degradation rate by IP 3K 3 , K D is the Ca 2+ affinity of IP 3K 3 , and K 3 is the IP 3 affinity of IP 3K 3 . Based on the previous contributions of IP 3 , the total IP 3 is given by The Li-Rinzel model (Li and Rinzel, 1994) is used to model the Ca 2+ dynamics within the astrocyte cell. The model consists of three channels, J chan , J leak , and J pump , where J chan models the Ca 2+ channel opening based on the mutual gating of the Ca 2+ and IP 3 , J leak models the Ca 2+ leakage from the ER into the cytoplasm and J pump models how Ca 2+ is pumped out from the cytoplasm into the ER via Sarco-Endoplasmic-Reticulum Ca 2+ -ATPase (SERCA) pumps. The Ca 2+ model in the approach of De Pittà et al. (2009) is used in this work, and it is described by where J chan is Ca 2+ release depending on the Ca 2+ and IP 3 concentrations, J pump is the amount of stored Ca 2+ within the ER via the SERCA pumps, J leak is the Ca 2+ leaking out of the ER and h is the fraction of activated IP 3 R s . The parameters h ∞ and τ h are given by where J chan is given by where r C is the maximal Calcium-Induced Calcium Release (CICR) rate, C 0 is the total free Ca 2+ cytosolic concentration, C 1 is the ER/cytoplasm volume ratio, and m ∞ and n ∞ are the IP 3 Induced Calcium Release (IICR) and CICR channels respectively, which are given by J leak and J pump are described by and where r L is the Ca 2+ leakage rate, v ER is the maximum SERCA pump uptake rate and k ER is the SERCA pump activation constant. The intracellular astrocytic calcium dynamics are used to regulate the release of glutamate from the astrocyte: the Glu pathway. To model this release, it is assumed that when Ca 2+ crosses the CICR threshold, a quantity of glutamate is released (Wade et al., 2012). It is described by where Glu is the quantity of released glutamate, τ Glu is decay rate of glutamate, r Glu is production rate of glutamate, and t Ca is the time at which Ca 2+ crosses the threshold. The released glutamate drives the generation of e-SP (Wade et al., 2012). The level of e-SP is modeled by where τ eSP is the decay rate of Glu, and m eSP is a constant weight used to control the height of e-SP. It shows that the e-SP level depends on the glutamate released from the astrocyte cell. The model of DSE in the approach of Wade et al. (2012) is used to describe the relationship between the DSE and the released 2-AG from postsynaptic neuron. The DSE is assumed to change linearly with the cytosolic concentration of 2-AG, which is described by where AG is the concentration of 2-AG and K AG is the scaling factor for the DSE.

Synapse Model
For the synapse, a probabilistic model is employed which is based on the failure and success mechanisms of synaptic neurotransmitter release (Navarrete and Araque, 2010;Wade et al., 2012). A uniformly distributed pseudo-random number generator is used. If the generated random number rand is less than or equal to the PR, a current I inj is injected into the neuron which is shown by where r I is the current production rate, and w i syn is the weight of the ith synapse. The associated PR of each synapse is determined by the DSE and e-SP together, which is given by where PR(t 0 ) is the initial PR for each synapse. As discussed in section 2.1, the PR can switch on/off learning at the synaptic terminal by modulating the height of the plasticity learning window. The authors are not aware of any biophysical model that relates PR to the plasticity window weighting parameter A 0 and therefore in this work it is assumed that A 0 is modulated according to where PR * is the learning activation level and r is a constant value which controls the maximum height of the learning window. The STDP rule used in this approach to update the synaptic weights according to the timing difference between the post and presynaptic spikes is described by where δw( t) is the weight update, t is the time difference between the post and presynaptic spikes, A 0 is the height of the plasticity window which limits the maximum levels of weight potentiation and depression, and τ + and τ − control the width of the plasticity window. A symmetrical plasticity window is assumed in this approach, and τ + = τ − = 40ms. From the proposed models, it can be seen that if the f pre is large enough, IP 3 is generated sufficiently to cause Ca 2+ oscillations. Then the Ca 2+ -induced glutamate binds to the mGluRs receptors at the presynaptic terminal resulting in an increase of the synaptic transmission probability PR. Note that the authors wish to point out that astrocytes are believed to gate LTP and LTD by regulating glutamate levels in the synaptic cleft (Foncelle et al., 2018). Since there are many complex biophysical mechanisms involved in the regulation of glutamate, which are still under debate, the authors take the view that modulating the STDP plasticity window using PR is an effective way to capture this gating function. Elevating PR opens the synaptic plasticity learning window and over time f post gradually increases which, via the 2-AG pathway, contributes to astrocytic IP 3 level until the Ca 2+ oscillation stops. This is accompanied by a reduction in the Glu (e-SP) pathway and PR falls causing a reduction in f post . Therefore, the bursting activity of the postsynaptic neuron is regulated by the GABA interneuron and the astrocyte cell. The results in the next section show the signaling pathways leading to a bursting postsynaptic neuron.

RESULTS
This section provides simulation results which highlight the dynamic behavior at the synapse terminals and how the interactions between an astrocyte and GABA interneuron can give rise to bursting behavior. The MATLAB simulation platform is used in this work together with the Euler method with the time step of 1 ms. Tables A1, A2 give all the model parameters.

Bursting Output Spike Pattern
In this simulation both the presynaptic excitatory neuron and the GABA interneuron are stimulated by the same spike train at frequency f pre = f GABA which causes the release of GABA and glutamate (Perea et al., 2016). The presynaptic excitatory neuron/GABA interneuron stimulus is 40 Hz in the following simulations as this is sufficient to produce a cytosolic [IP 3 ] > 0.5µM. With PR > PR * (see Figure 3), a significant increase occurs in the level of neurotransmitter in the cleft, the learning window opens (Figure 4) and weight potentiation starts (Figure 5) resulting in postsynaptic firing. Note that in Figure 4, the plasticity window height parameter A 0 increases periodically with a corresponding potentiation of the synaptic weight ( Figure 5). In our model the resting level for PR is 0.1 and based on the model in Equation (27), if PR > PR * (PR * = 0.45 in this work), the STDP learning window opens (A 0 > 0) at ∼80 s, as shown in Figure 4. After the synaptic weight is potentiated, the synapse generates a depolarising current only when the input stimulus is presented at the presynaptic terminal and the PR value is greater than the value of a random number (see probabilisticbased synapse model in Equation 25): this current is injected into the postsynaptic LIF neuron. This injected current increases the postsynaptic potential and the neuron fires a spike if the membrane potential is greater than the firing threshold, v th .
After a period of learning the postsynaptic neuron activity has stabilized and PR drops sufficiently, toward the end of the first set of PR "spikes" (see Figure 3), closing the plasticity window (A 0 = 0) and the weight stabilizes to ∼610 at 110 s, as shown in Figure 5: note that because the postsynaptic neuron is now active, PR < PR * for all subsequent Ca 2+ oscillations as the DSE pathway is also active. Figure 6 shows the amount of GABA released by the GABA interneuron as a function of time where, as expected, GABA increases gradually and then stabilizes at 0.027 µM under the input spike stimulus. IP GABA 3 is shown in Figure 7 (blue) as a function of time and stabilizes at ∼ 0.58µM which is consistent with the input stimuli profile. Figure 7 shows the other IP 3 sources that contribute to the total IP 3 in the cytosol.
Initially the total IP 3 increases with IP GABA 3 until the T CICR level is reached triggering the release of Ca 2+ , as shown in Figure 8. We have observed from our model that T CICR is consistent with an IP 3 level of approximately 0.5 µM and whenever IP 3 exceeds this threshold a transient elevation in Ca 2+ occurs, as can be seen in Figure 8. Note however that as the IP 3 level increases with IP AG 3 the degradation in IP 3 due to elevated Ca 2+ /IP 3 levels is insufficient to reduce IP 3 to below 0.5 µM and consequently the transient elevations of Ca 2+ stops just after 100 s followed by a relatively slow degradation of IP 3 and Ca 2+ : these periodic bursts in Ca 2+ gives rise to a Ca 2+ oscillatory wave where the initial Ca 2+ burst is longer due to synaptic   potentiation. At the onset of each subsequent Ca 2+ burst the IP 3 level drops sharply and we attribute this to strong dependence of IP 3K 3 on Ca 2+ (Equation 10). As the Ca 2+ level drops J chan (Equation 17) reverses direction perturbing the rate of change in Ca 2+ (Equation 12) and this causes a rapid increase in IP 3K 3 and a corresponding decrease in IP 3 (Equation 11).
The Ca 2+ oscillation is initiated at ∼20 s (T CICR is exceeded; Figure 8) and this triggers the release of glutamate targeting group I mGluRs on the presynaptic terminal, i.e., Glu (e-SP) pathway is activated (Figure 9). Figure 9 shows that the Glu (e-SP) signal accumulates at each CICR and rapidly decays after the Ca 2+ transients have ceased at ∼120 s. Also the DSE pathway increases as the activity of the postsynaptic neuron is increasing, and competes with the Glu (e-SP) pathway to restrict PR to a relatively stable low value for all subsequent Ca 2+ oscillations that occur post-learning, as shown in Figure 10. Again note a longer period of elevation of the DSE signal at the start due to synaptic potentiation but thereafter the DSE profile repeats in time. Figure 11 shows the firing rate of the postsynaptic neuron, which is calculated based on a sliding time window of 10 s (blue) and 40 s (red), respectively. Note that the first burst reaches a higher level of postsynaptic neuron activity when compared to subsequent bursts. Also, between bursts FIGURE 8 | Ca 2+ oscillations in the astrocyte cell as a function of time. Note a longer oscillatory period at the start due to learning but thereafter the oscillator period stabilizes with constant on/off ratio.  the activity never falls back to zero. This is because the first burst occurs during the weight potentiation phase when the synaptic weight is continually updated and eventually stabilized, whereas in all subsequent neuronal bursts no weight potentiation occurs. Clearly from Figure 11 a continual postsynaptic bursting behavior is evident.
Referring to Figure 12, we show simulations for f pre of 20, 40, and 80 Hz where clearly only f pre = 40Hz results in repeated Ca 2+ oscillations. This is because at 20Hz the T CICR level cannot be reached whereas at 80Hz the astrocyte cytosol is quickly swamped with both IP 3 and Ca 2+ and subsequent degradation in IP 3 is insufficient to allow further CICR. Consequently our model shows presynaptic frequency selectivity which is consistent with work reported elsewhere (Bienenstock et al., 1982;Dong et al., 2015).
In addition, as the morphology of GABA interneurons and receptor density at the astrocyte cell differ, the IP GABA 3 levels vary under the same input f pre . IP GABA 3 is a main contributor to the total IP 3 , thus the greater IP GABA 3 , the longer the process of IP 3 degradation. Figure 7 shows that when IP 3 degrades sufficiently to once again enable IP 3 to cross T CICR from below, a transient elevation in Ca 2+ results and PR increases with a corresponding increase in the postsynaptic neuron burst frequency. Therefore, FIGURE 11 | Firing rate of the postsynaptic neuron where a continual bursting behavior is evident. The firing activity was calculated using a sliding time window of 10 s (blue) and again for a sliding window of 40 s (red) where the latter gives a better average. FIGURE 12 | Astrocytic Ca 2+ as a function of time with f pre of 20, 40, and 80 Hz as a parameter. Note that for the extreme cases of (20 or 80 Hz) no Ca 2+ oscillations occur: for f pre = 20 Hz T CICR can never be achieved and at 80Hz degradation of IP 3 is insufficient to allow CICR to repeatedly occur. Consequently, there is a frequency window over which oscillations can occur. level, a long time period is required to degrade the overall IP 3 and to restart the Ca 2+ oscillation. Thus, the bursting frequency is low. different IP GABA 3 levels lead to different burst frequencies of the postsynaptic neuron. To determine the dependency of neuronal burst frequency on the production rate of IP GABA 3 , r GABA ip3 , a simulation was carried out (Figure 13) which shows the firing rates of the postsynaptic neuron under different production rates with f pre fixed at 40 Hz. It can be seen that when the r GABA ip3 increases, the frequency of the bursting decreases. For example, for the first 1,000 s, there are 6, 5, 4 bursts under the IP GABA 3 production rates (r GABA ip3 ) of 1.8, 2, and 2.2, respectively. This is because a high IP GABA 3 level requires a significant time period to degrade the total IP 3 , and to restart the Ca 2+ oscillation and bursting behavior, thus the bursting frequency is low. Note that a fixed frequency of the input stimulus (i.e., f pre = 40 Hz) is used in this experiment, however the same results are observed for other f pre values such as 50 Hz, and the burst frequency variation is not constrained for specific f pre values. The results in Figures 12, 13 demonstrate the functionalities of the GABA interneuron including the presynaptic frequency selectivity and postsynaptic bursting frequency regulation.

CONCLUSIONS
In this paper, a biophysical model is proposed where it is shown that GABA interneuron regulates the astrocytic IP 3 secondary messenger and thus the probability of release (PR) at the presynaptic terminal. In our model we propose that PR modulates the height of the plasticity window and therefore controls when synaptic potentiation/depression occurs. Specifically, the simulations show that during the weight potentiation phase, increasing IP 3 leads to a cycle of CICR events where each is followed by rapid degradation in IP 3 . Over time the firing frequency of the postsynaptic neuron continually increases and eventually the synaptic weights stabilize. Postsynaptic firing results in the release of 2-AG into the extracellular space and this messenger binds to CB1R receptors on the astrocyte membrane. The associated IP AG 3 contributes to the total cytosolic IP 3 and eventually Ca 2+ oscillations, and therefore the Glu (e-SP) pathway ceases: 2-AG also binds to CB1Rs on the presynaptic terminal causing a decrease of the synaptic transmission PR via the DSE pathway. PR therefore decreases at the presynaptic terminal which reduces the level of neurotransmitter in the cleft, and consequently the firing frequency of the postsynaptic neuron diminishes, as does IP AG 3 . Thereafter, the total IP 3 /Ca 2+ degrades significantly over time but is replenished by IP GABA 3 and a subsequent cycle of CICR events commences-the Glu (e-SP) pathway is re-established with an associated increase in PR and the level of neurotransmitter in the cleft is raised. However, in this instance weight potentiation does not occur as PR < PR * . The postsynaptic neuron firing rate increases again until the Ca 2+ transients stop and thereafter the activity of the postsynaptic neuron falls off again. A network burst has occurred and this is followed by repeated bursts where each coincides with Ca 2+ transients: the network burst frequency correlates with the Ca 2+ oscillatory wave. In addition, the GABA released by the GABA interneuron controls the frequency range within which the network bursts can occur. Future work will further explore other neurotransmitters released by astrocytes such as Dserine and ATP, and also slow inward currents at the postsynaptic terminal as a result of glutamate release by astrocytes.
The authors recognize that this study is based on biological findings of the simplest signaling mechanisms involving astrocytic GABA responses and astrocytic glutamate signaling in presynaptic terminals that regulate network function. Other factors, such as astrocytic ATP/adenosine release from astrocytes (Covelo and Araque, 2018), are not considered in the present model, but may also contribute to further shape of network activity, adding further complexity of the network effects of astrocyte signaling. Further studies incorporating these additional elements are therefore required to get a complete view of the astrocyte roles in network function. Despite this the present findings have potential implications for the generation of normal and pathologic circuit behavior in the brain, relevant to brain diseases that feature altered synaptic properties or where there is a propensity for the episodic synchronized bursting behavior of neurons. The electroencephalogram (EEG) is a composite product of population-level neuronal firing patterns of differing frequencies. Our findings suggest that GABA-B signaling via astrocytes may be relevant to the generation of certain frequencies and behaviors in the EEG. Seizures are the hallmark of the common brain disease epilepsy and are generated by hyper-synchronous discharges of populations of neurons. Notably, gene expression levels of key components modeled here, including the IP 3 receptor and GABA-B receptor, are dysregulated in human epileptic brain tissue or animal models (Matsumoto et al., 1996;Nishimura et al., 2005;Sheilabi et al., 2018) of epilepsy. Mutations in these genes have also been identified in individuals with epilepsy (Møller et al., 2017;Yoo et al., 2017). Indeed, the GABA-B receptor is a longstanding therapeutic target for the treatment of epilepsy (Bowery, 2006), and more recently the IP 3 receptor was reported to be a target of levetiracetam, one of the most effective antiepileptic drugs (Nagarkatti et al., 2008). The present model offers a novel mechanism to explain how astrocyte-neuron interactions regulate seizure-like activity (Gómez-Gonzalo et al., 2010), and how alterations in the described pathways may contribute to hyper-synchronous firing. It may also offer therapeutic insights through targeted manipulation of the astrocytic GABA-B or IP 3 systems followed by evaluation of the resting electroencephalogram (EEG) and investigating whether this alters the frequency or occurrence of pathophysiological neuronal firing and seizures.

DATA AVAILABILITY
The raw data supporting the conclusions of this manuscript will be made available by the authors, without undue reservation, to any qualified researcher.