Computational models for two neuron/astrocyte networks are developed to explore mechanisms underlying the astrocytes’ role in maintaining neuronal firing patterns. For the first network, a single neuron receives periodic excitatory inputs at varying frequencies. We consider the role played by several astrocytic dendritic processes, including the Na+-K+ ATPase pump, K+ channels and gap junctions in maintaining extracellular ion homeostasis so that the neuron can faithfully sustain spiking in response to the excitatory input. The second network includes two neurons coupled through mutual inhibitory synapses. Here we consider the role of astrocytic dendritic processes in maintaining anti-phase or synchronous oscillations. Dynamical systems methods, including bifurcation theory and fast/slow analysis, is used to systematically reduce the complex model to a simpler set of equations. In particular, the first network, consisting of differential equations for the neuron and astrocyte membrane potentials, channel state variables and intracellular and extracellular Na+ and K+ concentrations, is reduced to a one dimensional map. Fixed points of the map determine whether the astrocyte can maintain extracellular K+ homeostasis so the neuron can respond to periodic input.

**Introduction**

**Methods**

*Network**Neuron*: The neuron’s membrane potential, VN, satisfies the equation

The first two terms on the right hand side of this equation correspond Na^{+} and K^{+} currents. These are defined as

\(I_{Na}=(g_{Na} m_∞^3 (V_N )h+g_{NaL} )(V_N-E_{Na} )\)

\(I_K=(g_k n^4+g_{KL} )(V_N-E_K )\)

Note that there are both Na^{+} and K^{+} leaks. As in (Huguet *et al*., 2016), we assume that *h = 1 − n* and n satisfies a differential equation of the form

\(\frac{dn}{dt}=ϕ(n_∞ (V_N )-n)/τ_n (V_N)\) (2)

If *X *= *m *or *n*, then

\(X_∞ (V)=\frac{1}{1+e^{-(V-θ_X)/σ_x }}\) and \(τ_n (V)=τ_0+\frac{τ_1-τ_1}{1+e^{-(V-θ_n0)/σ_n0 }}\) (3)

The Nernst potentials are given by

\(E_K=\frac{RT}{F} ln \frac{K_e}{K_i}\) and \(E_{Na}=\frac{RT}{F} ln \frac{Na_e}{Na_i} \) (4)

where R,T and F are the gas constant, temperature and Faraday’s constant, respectively, and K_{i}, K_{e}, Na_{i} and Na_{e} are the K^{+} and Na^{+} concentrations in the neuron’s cytosol and extracellular space.

As in (Kager *et al*., 2002) the Na^{+}-K^{+} ATPase pump current is given by

\(I_{PN}=ρ_N (\frac{K_e}{2+K_e })^2 (\frac{Na_i}{7.7+Na_i} )^3\) (5)

where *r** _{N}* represents the maximal pump current.

The term *I _{exc}*corresponds to periodic, excitatory synaptic input and is given by

\(I_{exc}=g_{exc} s(t)(V_N-E_{syn}) .\)

Suppose that the input is at *f _{r}* hz. For each integer

*j*and

*t*

_{j}= j***

*f*we let

_{r}/1000,*s(t*and

_{j}) = 1*s’(t) = **-**b** s(t)* for *t _{j} < t < t_{j+1.}*

The neuron model parameters are: C_{m} = 1 mF/cm^{2}, g_{Na} = 20, g_{K} = 3, g_{NaL} = .03, g_{KL} = .2, j = .1, q_{m} = -37, s_{m} = 10, s_{n} = 10, q_{n} = -55, t_{0} = .1, t_{1} = 1, q_{n0} = -40, s_{n0} = -12, g_{exc} = 2, b = 1 ms^{-1}. Units for the maximal conductances (g’s) are mS/cm^{2} and half activation variables (q’s) are mV.

**Astrocyte:** The astrocyte’s membrane potential, *V _{A}*, satisfies the equation

\(C_m^A \frac{dV_A}{dt}=-I_K^A-I_{Kir}-I_{Na}^A-I_{P}A-I_{gap}.\) (6)

where \(C_m^A\)= 1u\(\)F/cm^{2} and \(I_K^A\) and \(I_Na^A\) correspond to K^{+} and Na^{+} leak currents; these are given by

\(I_K^A=g_K^A (V_A- E_K^A) \) and \(I_{Na}^A=g_{Na}^A (V_A- E_{Na}^A)\)

The term *I _{Kir}* corresponds to an inward rectifying K

^{+}current (Ransom, 1996), and is given by

\(I_{Kir}=g_{Kir} \left({\frac{K_e^{1/2}}{1+e^{(V_A-E_K^A)/19.2}}}\right)(V_A-E_K^A)\) (7)

The Nernst potentials, \(E_K^A \) and \(E_{NA}^A\), as well as the ATPase pump current *I _{PA}*

_{,}are defined very similar to (4) and (5).

The term *I _{gap}* corresponds to electrical coupling with another astrocyte that is electrically coupled with other astrocytes within a syncytium. We assume that this other astrocyte remains in steady state; that is, its membrane potential, , and intracellular K

^{+}and Na

^{+}concentrations, and , are constant. We assume that = -90 mV, = 135 mM and = 12 mM.

As in (Ma *et al*., 2016; Huguet *et al.*, 2016), we model *I _{gap}* using the Goldman-Hodgkin-Katz equation. That is,

*I*where

_{gap }= I_{Kgap}+ I_{Nagap}\(I_{Kgap}= P_{Kgap} F_{φA} \left(\frac{K_{iA} e^{(-φ_A )}-K_iA^0}{e^{-φ_A} -1} \right)\) (8)

\(I_{Nagap}= P_{Nagap} Fφ_A \left( \frac{Na_{iA} e^{(-φ_A )}-Na_{iA}^0}{e^{-φ_A} -1} \right)\)

Here, *\(\varphi{A} = (F/RT)(V_A-V_A^0)\).* The constants *P _{Kgap}* and

*P*are the K

_{Nagap}^{+}and Na

^{+}permeabilities, averaged over the entire cell membrane. Following (Ma

*et al*., 2016), we assume that

** P_{Kgap} = d_{gap} P_{K}** and

**P**_{Nagap }= 0.8 P_{Kgap}where *P _{K}* = 6 *10

^{−5 }cm/s. We will vary

*d*and

_{gap}, g^{A}_{K}*g*in the analysis to determine how the dynamics depends on gap junction coupling strength and K

_{Kir}^{+}currents.

**Ion concentrations: **The neuron’s intracellular K^{+} and Na^{+} concentrations satisfy the equations

\(\frac{dK_i}{dt}=-\frac{10 S_N}{F Ω_N }(I_K-2I_{PN})\) (9)

\(\frac{dNa_i}{dt}=-\frac{10 S_N}{F Ω_N }(I_{Na}-2I_{PN})\)

where *S _{N}* is the neuron’s surface area and W

_{N}is the volume of the neuron’s cytoplasm. The factor ’10’ is needed for consistency of units.

The astrocyte’s intracellular K^{+} and Na^{+} concentrations satisfy the equations

\(\frac{dK_iA}{dt}=-\frac{10 S_A}{F Ω_A }(I_K^A+I_{Kir}+I_{Kgap}-2I_{PA})\)

\(\frac{dNa_iA}{dt}=-\frac{10 S_A}{F Ω_A }(I_{Na}^A+I_{Kir}+I_{Nagap}-2I_{PA})\)

where *S _{A}* is the astrocyte’s surface area and Ω

_{A}is the volume of the astrocyte’s cytoplasm.

The extracellular K^{+} and Na^{+} ion concentrations satisfy the equations

\(\frac {dK_e}{dt}=\frac {10 S_N}{F Ω_N } (I_K-2I_PN )+\frac {10 S_A}{F Ω_A }(I_K^A+I_{Kir}+I_{Kgap}-2I_{PA})\) (10)

\(\frac{dNa_e}{dt}=\frac{10 S_N}{F Ω_N} (I_{Na}+3I_PN )+\frac{10 S_A}{F Ω_A }(I_{Na}^A+I_{Nagap}+3I_PA)\)

where Ω_{E} is the volume of the extracellular space, which is assumed to be a fixed fraction of the neuron’s volume; that is, Ω_{E} = α_{0} Ω_{N}*.* In the simulations, we let *S _{N}*= 10

^{4}μm

^{2}, Ω

_{N}= 5 10

^{3}μm

^{3},

*S*= 1.6 10

_{A}^{3}μm

^{2}, Ω

_{A}= 2 10

^{3}μm

^{3}and α

_{0}= .3.

**Network II**

Network II consists of two mutually coupled neurons. Each neuron shares extracellular space with an astrocyte, which, as before, is electrically coupled with another astrocyte. Both ’units’, consisting of a neuron, astrocyte and shared extracellular space, are modeled precisely as Network I above, except we replace the Na^{+} and K^{+} leaks in (1) by a single leak current, *I _{L} = g_{L}(V_{N }− E_{L}),* where

*g*= .02 and

_{L}*E*= −60 mV. Moreover, we replace the excitatory current,

_{L}*I*in (1) by a term corresponding to inhibitory synaptic input from the other neuron. That is, if j = 1 or 2, then we replace

_{exc},*I*in the voltage equation for neuron j with

_{exc }\(I_{syn}^j=g_{syn} s_k (V_N^j-E_{syn})\)

where \(k \neq j \) and *s _{k}* satisfies the equation

\(s_k^,=α(1-s_k ) s_∞ (V_N^k )-βs_k\)

Here, *E _{syn}* = −85 mV is the synaptic reversal potential and

*s*is defined as in (3). The parameter values for Network II are the same as Network I, except

_{∞}*θ*=−55, t

_{n}_{1}=2,

*φ*=.2,

*g*=.05 and

_{syn}*β*=.18.

**Results**

**Solutions of the two networks**

Solutions of Network I are shown in Fig. 1. This figure demonstrates that whether the model can maintain repetitive spiking depends on several factors including the input frequency (*f _{r})*, the strengths of the Na

^{+}-K

^{+}ATPase pumps (

*ρ*and

_{N}*ρ*) and the strength of gap junction coupling (

_{A}*d*). Fig. 1A shows that the neuron can maintain a 30 hz firing rate for 5 seconds for moderate levels of gap junction coupling and Na

_{gap}^{+}-K

^{+}ATPase pump strengths. However, as shown in Fig. 1B, if the input rate is increased to 40 hz, then the neuron stops firing at around 3 seconds, at which time it goes into so-called depolarization block. If we remove gap junction coupling by setting

*d*, then, as shown in Fig. 1C, the neuron cannot maintain periodic firing, even at a lower input rate of 10 hz. If we increase the Na

_{gap}= 0^{+}-K

^{+}ATPase pump strength, then, as shown in Fig. 1D, the neuron is able to maintain periodic firing.

**Figure 1:** Solutions of Network I for different values of input frequency (fr), the strengths of the neuron’s Na+ -K+ ATPase pump (rN) and the strength of gap junction coupling (dgap). For each plot, rA = .5, g^{A}_{K}, b = 3 and gKir = 0. A) With gap junctions, the neuron can maintain firing for moderate pump strengths and input frequencies. B) Even with gap junctions, sufficiently high input rates lead to depolarization block. C) Without gap junctions, the neuron cannot maintain firing at 10 hz, unless, as shown in D), the pump strength is sufficiently high.

In Fig. 2 we illustrate the range of values of the Na^{+}-K^{+} ATPase pump strengths (with *ρ _{N} = ρ_{A}*) and input frequencies for which the neuron is able to maintain spiking for 10 seconds. We consider the two cases:

*(*

*, g*and

_{Kir}^{A}) = (3, 0)*(*

*, g*We choose a smaller value for

_{Kir}) = (0, 1).*g*to account for the K

_{Kir}^{+}dependence of

*I*open probability. The neuron exhibits depolarization block for values of pump strengths and input frequencies above each curve.

_{Kir }Without gap junctions (*d _{gap} = 0*), there is almost no difference between these two cases in the frequencies at which the neuron can maintain steady spiking. In particular, without gap junctions, the inward rectifying K

^{+}current does not seem to enhance K

^{+}buffering, which is needed to prevent

*K*from rising above the threshold for depolarization block.

_{e}**Figure 2: **The range of values of the Na^{+}-K^{+} ATPase pump strengths (with r_{N} = r_{A}) and input frequencies for which the neuron is able to maintain spiking for 10 seconds. The neuron exhibits depolarization block for values of pump strengths and input frequencies above each curve. The d_{gap} = 0 (blue) curve corresponds to both cases: ( , g_{Kir}) = (0,1) and ( , g_{Kir})= (3,0).

With gap junctions (*d _{gap} = 1*), the neuron can maintain higher firing rates with just

*I*than with just if the pump strengths are sufficiently strong. This is mainly because the threshold for

_{Kir }*K*when the neuron exhibits depolarization block increases at higher pump strengths. (This will be demonstrated later.) Large

_{e }*K*values strengthen

_{e}*I*and enhance K

_{Kir}^{+}buffering.

Fig. 3 shows solutions of Network II. With gap junction coupling (*d _{gap} = 1*), the network maintains anti-phase spiking. However, without gap junction coupling (

*d*) the network switches from anti-phase to synchronous spiking at around 4.5 seconds. For this simulation,

_{gap}= 0*g*and

_{Kir }= 1*= 0*. The result is almost identical if instead,

*g*and g

^{A}_{Kir}= 0^{A}

_{K}

*= 3.*

**Figure 3:** Solutions of Network II. A) With gap junctions, the two neurons exhibit anti-phase oscillations. B) Without gap junction coupling, the network switches from anti-phase to synchronous spiking at around 4.5 seconds. C) With gap junctions, the model maintains a nearly constant K_{e} level, but not without gap junctions. D) Blow up of solution shown in B).

**Isopotentiality and K ^{+} buffering**

Our results demonstrate that astrocytic gap junctional coupling plays a critical role in maintaining neuronal firing patterns. The basic mechanism underlying this behavior is so-called *K ^{+} spatial buffering:* gap junction coupling allows astrocytes to maintain a nearly constant extracellular K

^{+}concentration in the face of neuronal activity that would tend to increase it (Kofuji & Newman, 2004; Orkand

*et al*., 2006). If the gap junction coupling strength is weak, then the astrocytes are not able to clear elevated extracellular K

^{+}levels. This leads to increased neuronal excitability, which may change their firing properties. When extracellular K

^{+}concentrations rise above some threshold, the neurons exhibit depolarization block.

How does strong gap junction coupling lead to K^{+} spatial buffering? This depends on so-called isopotentiality; that is, an astrocyte’s associated syncytium provides powerful electrical coupling that equalizes the astrocyte’s membrane potential with its neighbors (Ma *et al*., 2016). To understand why isopotentiality leads to K^{+} spatial buffering, we first note that K^{+} currents are normally outward; that is, they allow for the flow of K^{+} ions from the cell’s inside to the extracellular space. In general, we can express a K^{+} current as *I _{K} = g_{K}(V_{A} − *

*).*If

*V*

_{A}>*,*then

*I*and the current is outward. In order for the K

_{K}> 0^{+}currents to be inward, the cell’s membrane potential must lie below the K

^{+}reversal potential. When neurons spike, they release K

^{+}into the extracellular space. Hence,

*K*increases, as does = (RT/F) ln(

_{e}*K*). However, because of isopotentiality, the astrocyte’s membrane potential remains nearly constant. This allows for

_{e}/K_{iA}*V*, so that the astrocyte’s K

_{A}<^{+}currents become inward and are, therefore, able to clear K

^{+}from the extracellular space.

In Fig. 4, we consider Network I and plot the response of *V _{A}* and

*K*to different excitatory input rates with (Fig. 4A,B) and without (Fig. 4C,D) gap junction coupling. When

_{e}*d*, both the astrocyte’s membrane potential and extracellular K

_{gap}= 1^{+}concentration remain nearly constant. However, when there are no gap junctions (

*d*), there is a steady rise in both

_{gap}= 0*V*and

_{A}*K*even at low input rates.

_{e}**Figure 4:** Solutions of Network I showing V_{A} and K_{e}. Here, g_{Kir} = 3 and = 0. A,B) With gap junctions the network can maintain nearly constant V_{A} and K_{e} , even at 20 hz input. C,D) This is not the case if there are no gap junctions.

Fig. 5 shows plots of *V _{A} ,* and

*I*, with and without gap junctions, for the same solutions shown in Fig. 4 with an input rate of 10 hz. With gap junctions,

_{Kir}*V*falls below and

_{A}*I*reverses (

_{Kir}*I*< 0). Without gap junctions,

_{Kir}*V*tracks very closely with and

_{A}*I*remains negligible.

_{Kir}Fig. 3C shows that for Network II, with strong gap junction coupling, *K _{e}* remains nearly constant; however, without gap junction coupling, there is a steady rise in

*K*during neuronal firing.

_{e}**Figure 5:** Plots of A) V_{A} and , and B) I_{Kir}, with and without gap junctions, for the same solutions shown in Fig. 4 with an input rate of 10 hz. With gap junctions, V_{A} falls below and I_{Kir} reverses (I_{Kir} < 0). Without gap junctions, V_{A} tracks very closely with and I_{Kir} remains negligible.

**Analysis**

We mathematically analyze Network I by first making some simplifying assumptions and then reducing the full model to a simpler set of equations. The analysis leads to a single equation for just *K _{e}*. This will be used to determine the

*K*threshold for when the neuron exhibits depolarization block and help explain the response of the model to excitatory input.

_{e}We begin by noting that the total amounts of K^{+} and Na^{+} ions are conserved. That is

\(Ω_e K_e+Ω_N K_i+Ω_A K_iA=K_{tot}\) (11)

and

\(Ω_e {Na}_e+Ω_N Na_i+Ω_A Na_iA=Na_{tot}\)

are constant. We assume that

\(K_{tot}=Ω_e K_e^0+Ω_N K_i^0+Ω_A K_iA^0\) and \(Na_{tot}=Ω_e Na_e^0+Ω_N Na_i^0+Ω_A Na_iA^0\)

where

\(K_e^0=4 \) mM, \(Na_e^0\) =135 mM, \(K_i^0=K_iA^0=\)135 mM and \(Na_i^0=Na_iA^0=\)12 mM.

We now make several assumptions in order to derive a reduced model. We first assume intracellular neuron electroneutrality; that is, *K _{i} + Na_{i}* = + is constant. Our next assumption is that the astrocyte’s intracellular K

^{+}and Na

^{+}concentrations are constant; that is, and . With these assumptions we can solve for all of the ion concentrations in terms of

*K*and the Network I model reduces to equations for just

_{e}*V*and

_{N}, n, V_{A}*K*.

_{e}We can reduce the model further using fast/slow analysis. The membrane potential *V _{A}* clearly evolves on a time-scale much faster than the ion concentration

*K*. Since the astrocyte does not ’spike’, we may assume that

_{e}*V*is close to steady state; that is, the right hand side of (6) is, to leading order, zero and we can solve for

_{A}*V*in terms of the other variables. Note that if there is no

_{A}*I*current (

_{Kir}*g*0) then the right hand side of (6) is, in fact, a linear function of

_{Kir}=*V*. This is because if \(K_{iA}=K_{iA}^0\) and \(Na_{iA}=Na_{iA}^0\) , then the gap junction current can be written as

_{A}\(I_{gap}=g_{gap} (V_A-V_A^0)\)

where

\(g_{gap}=\frac{F^2}{RT} d_{gap} P_K (K_{iA}^0+.8Na_{iA}^0 ).\) (12)

Setting the right hand side of (6) equal to zero, we find that

\(V_A=\frac{g_K^A E_K^A+g_{Na}^A E_{Na}^A+I_{PA}-I_{exc}}{g_K^A+g_{Na}^A-g_{gap}}\) (13)

The full Network I model is now reduced to equations for just *V _{N}, n* and

*K*, which we write as

_{e}\(C_m \frac{dV_N}{dt}=-I_{Na}-I_K-I_{PN}-I_{exc}\)

\(\frac{dn}{dt}=\frac{ϕ(n_∞ (V_N )-n)}{τ_n (V_N )}\)

\(\frac{dK_e}{dt}=Φ(V_N,n,K_e )\) (14)

Assume, for now, that there is no excitatory input. Then fast/slow analysis is used to analyze the reduced system. If we consider the slow variable *K _{e}* to be a bifurcation parameter in the fast subsystem for

*(V*, then the resulting bifurcation diagram is shown in Fig. 6. Note that there is a stable fixed point for \(K_e<K_HB^1≈10.35\) and \(K_e>K_{HB}^2≈25.75\). There is a subcritical Hopf bifurcation at and a supercritical Hopf bifurcation at . Moreover, there are stable periodic orbits for \(K_{HB}^1<K_e<K_{HB}^2.\)

_{N}, n)We next consider the evolution of the slow variable *K _{e}*. This is done using the method of averaging. Denote the fixed points of the fast subsystem as

*(V*and the stable periodic orbits as

_{fp}(K_{e}), n_{fp}(K_{e}))*(V*. Let

_{P}(t, K_{e}), n_{P}(t, K_{e}))*T(K*be the period of the periodic orbits. Near the branch of stable fixed points,

_{e})*K*satisfies, to leading order,

_{e}\(\frac{dK_e}{dt}= Φ(V_{fp} (K_e ),n_{fp} (K_e ),K_e ) =Φ_{fp} (K_e )\)

**Figure 6:** Bifurcation diagram for (14) with bifurcation parameter K_{e}. There are both supercritical and subcritical Hopf bifurcations and a branch of stable periodic orbits.

For \(K_{HB}^1<K_e<K_{HB}^2, K_e\) satisfies, to leading order,

\(\frac{dK_e}{dt} = \frac{1}{T(K_e)} ∫_0^{T(K_e)}Φ(V_P (t,K_e ),n_P (t,K_e ),K_e ) dt = Φ_{ave} (K_e )\)

We compute \(Φ_{fp} (K_e )\) and \(Φ_{ave} (K_e )\) numerically and the result is shown in Fig. 7A.

We first consider the *K _{e}* dynamics near the branch of stable fixed points of the fast subsystem with . Note that there exists

*K*<

_{P}*K*

^{1}_{HB}such that \(Φ_{fp} (K_e )=0\) . Moreover, \( Φ_{fp} (K_e )>0\) for

*0 < K*and \(Φ_{fp} (K_e )<0\)for

_{e}< K_{P}*K*

_{P}< K_{e}<*K*

^{1}_{HB}. This implies that

*K*as

_{e}→ K_{P}*t → ∞.*

We next consider the spiking regime when . As shown in Fig. 7A, . Hence, *K _{e}* must increase and the neuron must approach depolarization block.

This analysis demonstrates that is the threshold for when the neuron exhibits depolarization block. Fig. 7 shows that this threshold is an increasing function of the neuron’s Na^{+}-K^{+} ATPase pump strength.

We next consider the neuron’s response to periodic excitatory input. Whenever the neuron spikes, there is also a spike in *I _{K}*. This leads to a fast rise in

*K*. We now assume that at the spike times,

_{e}*K*increases by a fixed amount, which we denote as

_{e}*K*

*. Based on numerics (see Fig. 4), we let*

_{d}*K*

*= .38 mM.*

_{d}**Figure 7: **A) Plots of Ф_{fp}(*K _{e}*) and Ф

_{ave}(

*K*). The threshold for depolarization block is

_{e}*K*=K

_{e}^{1}

_{HB}. B) Plot of the depolarization block threshold versus the neuron pump strength

*r*

*.*

_{N}Between spike times, *K _{e}* satisfies (14), which depends on

*V*and

_{N }*n*. We use fast/slow analysis to express these other variables in terms of

*K*. This is done by first setting the right hand sides of (1) and (2) equal to zero and then solving for

_{e}*(V*in terms of

_{N}, n)*K*. However, the right hand side of (1) is a nonlinear function of

_{e}*V*and

_{N}*n*. In order to obtain an explicit formula for

*V*, we note that between spikes,

_{N}*V*is near a resting state and the channel activation terms, \(m_∞^3 (V_N) \) and

_{N}*n*, are very small. We, therefore, use the approximation

^{4}** I_{Na} ≈ g_{NaL}(V_{N} − E_{Na})** and

**(15)**

*I*._{K}≈ g_{KL}(V_{N}− E_{K})In this case, the right hand side of (1) is linear in *V _{N}* and does not depend on

*n*. Setting the right hand side of (1) equal to zero, we find that

\(V_N=\frac{g_{NaL} E_{Na}+g_{KL} E_K-I_PN}{g_{NaL}+g_{KL} }\) (16)

We can now express F in (14) as a function of just *K _{e}*; we denote this function as \(\psi(Ke)\)

*.*This is done by first replacing

*I*by the approximation given in (15) and then using (16).

_{K}In summary, the reduced model is the following. Suppose that the excitatory input is at *f _{r}* hz. For each integer

*j*and

*t*/1000,

_{j}= j · f_{r }\(K_e (t_j^+ )=K_e (t_j^- )+K_δ\)

and

\(\frac{dK_e}{dt}=Ψ(K_e ) \) for *t _{j} < t < t_{j+1}* (17)

Solutions of this reduced model are shown in Fig. 8. With gap junctions (Fig. 8A) the neuron can maintain firing at 10 and 20 hz, since, in both cases, *K _{e}*remains below the threshold for depolarization block. Without gap junctions (Fig. 8B) the neuron cannot maintain firing at 10 hz if

*ρ*= .5, since

_{N}*K*increases past this threshold. However, increasing the neuron’s Na

_{e}^{+}-K

^{+}ATPase pump strength to

*ρ*= 1 allows the neuron to maintain a 10 hz firing rate.

_{N}Finally, we construct a 1-dimensional map; fixed points of the map correspond to the asymptotic behavior of *K _{e}* for the reduced model. To define the map, fix

*K*> 0 and let

_{0}*K*be the solution of (17) with

_{e}(t; K_{0})*K*. Then the map is defined as simply

_{e}(0; K_{0}) = K_{0}\(Π(K_0 )=K_e (\frac{f_r}{1000};K_0 )+K_δ\)

Note that a fixed point of this map corresponds to a periodic solution of the reduced model (17).

In Fig. 8C,D we plot P(*K _{e}*)−

*K*for the same parameter values used in Fig. 8A,B. Fixed points of P(

_{e}*K*) correspond to zeros of the corresponding curves. Note in Fig. 8D that if

_{e}*d*=0 and

_{gap}*ρ*=.5, then P(

_{N}*K*) >

_{e}*K*for all values of

_{e}*K*below the threshold for depolarization block. Hence,

_{e}*K*must steadily increase past the threshold.

_{e}**Figure 8:** A,B) Solutions of the reduced model (17). The dashed line corresponds to the threshold for depolarization block. C,D) Plots of the map P(K_{e}) − K_{e}. Fixed points P(K_{e}) correspond to zeros of the corresponding curves.

**Discussion**

The main goals of this paper were to: 1) develop computational models to study mechanisms underlying astrocytes’ role in maintaining neuronal firing patterns; and 2) use mathematical tools to systematically reduce the complex model to a simpler system in order to characterize how solutions depend on network parameters and cellular processes. Simulations of the computational model demonstrate the importance of gap junctional coupling in K^{+} spatial buffering and, thereby preventing elevated levels of extracellular K^{+} leading to depolarization block. Using dynamical systems methods, we reduced the full Network I model to a one-dimensional map. Fixed points of the map determine whether the astrocyte can maintain extracellular K^{+} homeostasis so the neuron can faithfully respond to periodic input.

The basic mechanism for extracellular K^{+} clearance described in this study extends the concept of *K ^{+}* spatial buffering, which was introduced more than a half century ago (Kofuji & Newman, 2004; Orkand

*et al*., 1966). There have been several experimental, modeling, and analytic studies of K

^{+}spatial buffering since then (Gardner-Medwin, 1983,Gardner-Medwin & Nicholson, 1983; Chen & Nicholson, 2000). However, the classic description does not take into account isopotentiality of the astrocyte syncytium. With sufficiently strong and widespread gap junction coupling, astrocytes near the region of elevated K

^{+}concentration do not depolarize significantly and this provides a powerful driving force, , for K

^{+ }uptake through membrane K

^{+}channels. This important role of syncytial isopotentiality was speculated in previous papers (Muller, 1996); however, this was not experimentally demonstrated until (Ma

*et al*., 2016). As demonstrated in (Ma

*et al*., 2016) syncytial isopotentiality minimizes the local high K

_{e}-induced

*V*depolarization, and this maintains a sustained driving force for K

_{A}^{+}uptake. By extension, syncytial isopotentiality also increases the driving force for K

^{+}release in distant regions where

*K*remains at the physiological level. Additionally, increase in both driving forces creates a maximum driving force for intracellular K

_{e}^{+}transfer from a high K

^{+}region to remote regions with normal K

^{+}. Therefore, syncytial isopotentiality facilitates all three critical steps in K

^{+}spatial buffering: K

^{+}uptake, intercellular transfer and release (Kofuji & Newman, 2004). Furthermore, recent experiments demonstrate that syncytial isopotentiality arises in several regions throughout the central nervous system and may be a unified mechanism governing the operation of astrocyte networks (Kiyoshi

*et al.*, 2019; Kiyoshi & Zhou, 2019; Huang

*et al.*, 2018).

Some papers have proposed that the inward rectifying K^{+} current, *I _{Kir}*, is primarily responsible for K

^{+}buffering. A computational model developed in (Sibille

*et al*., 2015), for example, suggests that astrocytic Kir4.1 channels are sufficient to account for elevated extracellular K

^{+}clearance, even without gap junctional coupling. However, there are important differences between the model presented in (Sibille

*et al*., 2015) and that developed in this paper. In particular, the astrocyte membrane equation in (Sibille

*et al*., 2015) contains a nonspecific leak current with a fixed reversal potential that helps stabilize the astrocyte membrane potential at -80 mV. This keeps the astrocyte’s membrane potential sufficiently hyperpolarized during neuronal activity so that

*I*can reverse to an inward current. In our model, the astrocyte’s membrane potential remains hyperpolarized due to gap junctions and isopotentiality.

_{Kir }There have been numerous earlier papers that have addressed various issues related to modeling signals in astrocyte-neuronal networks. In particular, several papers have introduced models for spreading depolarizations, spreading depression, epilepsy, persistent activity and the propagation of Ca^{2+} waves (Cressman *et al*., 2009; Frohlich & Bazhenov, 2006; Hubel & Dahlemm 2014; Hubel *et al*., 2014; Huguet *et al*., 2016; Kager *et al*., 2002; O’Connell & Mori, 2016; Somjen *et al*., 2008; Ullah *et al*., 2009; Wei *et al*., 2014; Zandt *et al*., 2011). Moreover, several papers have used dynamical systems methods to reduce the complexity of the models and analyze how the neuronal spiking activity depends on the astrocytes’ ability to maintain ion homeostasis (Barreto & Cressman, 2011; Cressman *et al*., 2009; Frohlich & Bazhenov, 2006; Oyehaug *et al*., 2012; Zandt *et al*., 2011). Many of these previous models assumed that the role of the astrocytes is to simply buffer extracellular K^{+}; this was modeled by including a simple buffering term in the equation for extracellular K^{+}. We have built on and extended previous modeling studies by incorporating a detailed biophysical model for the astrocytes, considering the role played by gap junctions and reducing a model for the response of neurons to excitatory input to a one dimensional map.

Attachment | Size |
---|---|

omp2019.001.0063.pdf | 3.65 MB |

This work was partially supported by NIH awards RO1NS062784 and R56NS097972.

Allaman I, Belanger M, Magistretti PJ. 2011. Astrocyte-neuron metabolic relationships: for better and for worse. Trends Neurosci 34:76-87.

Amiri M, Hosseinmardi N, Bahrami F, Janahmadi M. 2013. Astrocyte- neuron interaction as a mechanism responsible for generation of neural synchrony: a study based on modeling and experiments. J Comput Neurosci 34:489-504.

Attwell D, Buchan AM, Charpak S, Lauritzen M, Macvicar BA, Newman EA. 2010. Glial and neuronal control of brain blood flow. Nature 468:232-43.

Barres BA. 2008. The mystery and magic of glia: a perspective on their roles in health and disease. Neuron 60:430-40.

Barreto E, Cressman JR. 2011. Ion concentration dynamics as a mechanism for neuronal bursting. J Biol Phys 37:361-73.

Bellot-Saez A, Kekesi O, Morley JW, Buskila Y. 2017. Astrocytic modulation of neuronal excitability through K(+) spatial buffering. Neurosci Biobehav Rev 77:87-97.

Chen KC, Nicholson C. 2000. Spatial buffering of potassium ions in brain extracellular space. Biophys J 78:2776-97.

Chen Y, Swanson RA. 2003. Astrocytes and brain injury. J Cereb Blood Flow Metab 23:137-49.

Chever O, Dossi E, Pannasch U, Derangeon M, Rouach N. 2016. Astroglial networks promote neuronal coordination. Sci Signal 9:ra6.

Cressman JR, Jr., Ullah G, Ziburkus J, Schiff SJ, Barreto E. 2009. The influence of sodium and potassium dynamics on excitability, seizures, and the stability of persistent states: I. Single neuron dynamics. J Comput Neurosci 26:159-70.

Dirnagl U, Iadecola C, Moskowitz MA. 1999. Pathobiology of ischaemic stroke: an integrated view. Trends Neurosci 22:391-7.

Ermentrout B. 2002. Simulating, analyzing, and animating dynamical systems: a guide to XPPAUT for researchers and students. Society for Industrial and Applied Mathematics.

Ermentrout GB, Terman DH. 2010. Mathematical Foundations of Neuroscience (Springer Science & Business Media).

Frohlich F, Bazhenov M. 2006. Coexistence of tonic firing and bursting in cortical neurons. Phys Rev E Stat Nonlin Soft Matter Phys 74:031922.

Gardner-Medwin AR. 1983. Analysis of potassium dynamics in mammalian brain tissue. J Physiol 335:393-426.

Gardner-Medwin AR, Nicholson C. 1983. Changes of extracellular potassium activity induced by electric current through brain tissue in the rat. J Physiol 335:375-92.

Giaume C, Koulakoff A, Roux L, Holcman D, Rouach N. 2010. Astroglial networks: a step further in neuroglial and gliovascular interactions. Nat Rev Neurosci 11:87-99.

Haydon PG, Carmignoto G. 2006. Astrocyte control of synaptic transmission and neurovascular coupling. Physiol Rev 86:1009-31.

Huang M, Du Y, Kiyoshi CM, Wu X, Askwith CC, McTigue DM, Zhou M. 2018. Syncytial isopotentiality: an electrical feature of spinal cord astrocyte networks. Neuroglia 1:271-279.

Hubel N, Dahlem MA. 2014. Dynamics from seconds to hours in Hodgkin-Huxley model with time-dependent ion concentrations and buffer reservoirs. PLoS Comput Biol 10:e1003941.

Hubel N, Scholl E, Dahlem MA. 2014. Bistable dynamics underlying excitability of ion homeostasis in neuron models. PLoS Comput Biol 10:e1003551.

Huguet G, Joglekar A, Messi LM, Buckalew R, Wong S, Terman D. 2016. Neuroprotective Role of Gap Junctions in a Neuron Astrocyte Network Model. Biophys J 111:452-462.

Iadecola C, Nedergaard M. 2007. Glial regulation of the cerebral microvasculature. Nat Neurosci 10:1369-76.

Kager H, Wadman WJ, Somjen GG. 2002. Conditions for the triggering of spreading depression studied with computer simulations. J Neurophysiol 88:2700-12.

Kimelberg HK, Nedergaard M. 2010. Functions of astrocytes and their potential as therapeutic targets. Neurotherapeutics 7:338-53.

Kiyoshi CM, Du Y, Zhong S, Wang W, Taylor AT, Xiong B, Ma B, Terman D, Zhou M. 2018. Syncytial isopotentiality: a system-wide electrical feature of astrocytic networks in the brain. Glia.

Kiyoshi CM, Zhou M. 2019. Astrocyte syncytium: a functional reticular system in the brain. Neural Regen Res 14:595-596.

Kofuji P, Newman EA. 2004. Potassium buffering in the central nervous system. Neuroscience 129:1045-56.

Ma B, Buckalew R, Du Y, Kiyoshi CM, Alford CC, Wang W, McTigue DM, Enyeart JJ, Terman D, Zhou M. 2016. Gap junction coupling confers isopotentiality on astrocyte syncytium. Glia 64:214-26.

Moskowitz MA, Lo EH, Iadecola C. 2010. The science of stroke: mechanisms in search of treatments. Neuron 67:181-98.

Muller CM. 1996. Gap-junctional communication in mammalian cortical astrocytes: development, modifiability and possible functions. In: Spray DC, Austin DR Gap junctions in the nervous system Editors TX: RG Landes Company pp 203-212.

Nakase T, Fushiki S, Sohl G, Theis M, Willecke K, Naus CC. 2003. Neuroprotective role of astrocytic gap junctions in ischemic stroke. Cell Commun Adhes 10:413-7.

Nedergaard M, Dirnagl U. 2005. Role of glial cells in cerebral ischemia. Glia 50:281-6.

Nedergaard M, Ransom B, Goldman SA. 2003. New roles for astrocytes: redefining the functional architecture of the brain. Trends Neurosci 26:523-30.

Newman EA. 2003. New roles for astrocytes: regulation of synaptic transmission. Trends Neurosci 26:536-42.

O'Connell R, Mori Y. 2016. Effects of Glia in a Triphasic Continuum Model of Cortical Spreading Depression. Bull Math Biol 78:1943-1967.

Orkand RK, Nicholls JG, Kuffler SW. 1966. Effect of nerve impulses on the membrane potential of glial cells in the central nervous system of amphibia. J Neurophysiol 29:788-806.

Oyehaug L, Ostby I, Lloyd CM, Omholt SW, Einevoll GT. 2012. Dependence of spontaneous neuronal firing and depolarisation block on astroglial membrane transport mechanisms. J Comput Neurosci 32:147-65.

Pannasch U, Vargova L, Reingruber J, Ezan P, Holcman D, Giaume C, Sykova E, Rouach N. 2011. Astroglial networks scale synaptic activity and plasticity. Proc Natl Acad Sci U S A 108:8467-72.

Poskanzer KE, Yuste R. 2016. Astrocytes regulate cortical state switching in vivo. Proc Natl Acad Sci U S A 113:E2675-84.

Ransom BR. 1996. Do glial gap junctions play a role in extracellular ion homeostasis? In: Spray DC, Austin DR Gap junctions in the nervous system Editors TX: RG Landes Company pp 159-173.

Ransom CB, Sontheimer H. 1995. Biophysical and pharmacological characterization of inwardly rectifying K+ currents in rat spinal cord astrocytes. J Neurophysiol 73:333-46.

Sibille J, Dao Duc K, Holcman D, Rouach N. 2015. The neuroglial potassium cycle during neurotransmission: role of Kir4.1 channels. PLoS Comput Biol 11:e1004137.

Somjen GG. 2001. Mechanisms of spreading depression and hypoxic spreading depression-like depolarization. Physiol Rev 81:1065-96.

Somjen GG, Kager H, Wadman WJ. 2008. Computer simulations of neuron-glia interactions mediated by ion flux. J Comput Neurosci 25:349-65.

Szabo Z, Heja L, Szalay G, Kekesi O, Furedi A, Szebenyi K, Dobolyi A, Orban TI, Kolacsek O, Tompa T and others. 2017. Extensive astrocyte synchronization advances neuronal coupling in slow wave activity in vivo. Sci Rep 7:6018.

Ullah G, Cressman JR, Jr., Barreto E, Schiff SJ. 2009. The influence of sodium and potassium dynamics on excitability, seizures, and the stability of persistent states. II. Network and glial dynamics. J Comput Neurosci 26:171-83.

Wei Y, Ullah G, Schiff SJ. 2014. Unification of neuronal spikes, seizures, and spreading depression. J Neurosci 34:11733-43.

Zandt BJ, ten Haken B, van Dijk JG, van Putten MJ. 2011. Neural dynamics during anoxia and the "wave of death". PLoS One 6:e22127.

Zhao Y, Rempe DA. 2010. Targeting astrocytes for stroke therapy. Neurotherapeutics 7:439-51.