The Edge of Stability: Response Times and Delta Oscillations in Balanced Networks

The standard architecture of neocortex is a network with excitation and inhibition in closely maintained balance. These networks respond fast and with high precision to their inputs and they allow selective amplification of patterned signals. The stability of such networks is known to depend on balancing the strengths of positive and negative feedback. We here show that a second condition is required for stability which depends on the relative strengths and time courses of fast (AMPA) and slow (NMDA) currents in the excitatory projections. This condition also determines the response time of the network. We show that networks which respond quickly to an input are necessarily close to an oscillatory instability which resonates in the delta range. This instability explains the existence of neocortical delta oscillations and the emergence of absence epilepsy. Although cortical delta oscillations are a network-level phenomenon, we show that in non-pathological networks, individual neurons receive sufficient information to keep the network in the fast-response regime without sliding into the instability.

gamma oscillatory instability when ∆q is very large and positive. However, both of 11 these types of dynamics have been described previously [1,2] and are not the focus of 12 our work. 13 The single population model can be reduced to a third order system by recognizing 14 that each of the four synaptic equations has the same input, R(t). Therefore, setting 15 S ampa + = S ampa − and S nmda + = S nmda − reduces the system to third order. After canceling 16 out terms the new single population model is: where R represents the firing rate of the population, with intrinsic time constant τ e , and 18 w is the synaptic weight. S l represents the synaptic activation which decays 19 exponentially with time constant τ l . l is the synapse type, either AMPA or NMDA. ∆q Instabilities of the network can be found by computing the roots of the denominator in 28 Eq (S3). Any root with positive real part is unstable while a root with a real part of 29 zero is marginally stable. Rather than look for all of the roots of the polynomial we 30 search for a set of parameters which make the system marginally stable by assuming 31 that γ = ω 0 j where j is imaginary and ω 0 is the oscillatory frequency in radians per 32 second. Inserting ω 0 j into the polynomial gives, as well as the value of w∆q at instability, 36 w∆q = τ nmda τ ampa τ e (τ nmda − τ ampa ) (τ nmda τ ampa + τ nmda τ e + τ ampa τ e ) − τ nmda + τ ampa + τ e τ nmda − τ ampa (S6) The value of w∆q given by Eq (S6) is used for the simulation of undamped oscillations 37 in the main text, Fig 1C, left column, middle row. In the main paper we compare the dynamical response of the full balanced network to a 40 damped oscillator. This approximation is derived directly from the single population 41 model which we described in the previous section. We begin by rewriting Eq (S3) .

42
Multiplying the top and bottom of Eq (S3) by τ e τ nmda τ ampa then dividing the top and 43 bottom by (τ nmda γ + 1)(τ ampa γ + 1) yields, 44 R(γ) = 1 We then assume that the AMPA connections act instantaneously and therefore take the 45 limit as τ ampa approaches zero then rearrange into the standard form, The denominator of this transfer function describes a damped oscillator with dynamics 47 given by the differential equation where ω 0 and ζ are defined as, and I(t) is some time dependent input. Note that although we have taken the limit as  As can be seen from Eq (S12), the damping coefficient ζ is a linear function of w∆q. 55 This implies that ∆q can be used as a reasonable approximation to ζ and that its 56 impact on the network should be roughly equivalent. Given this damped oscillator 57 formulation we can also compute the value of w∆q where the network is critically 58 damped by setting ζ = 1 yielding w∆q = (2 τ nmda τ e − τ nmda − τ e )/τ nmda .

59
Additionally, we can determine the point at which the system becomes marginally 60 stable by setting ζ = 0. We also call this the undamped network in the main text. The 61 spring network is marginally stable when w∆q = −(1 + τ e /τ nmda ). These equations for 62 the critically damped and undamped networks were used to compute the values of w∆q 63 in Fig 1C, top row, used in the simulations. These equations also give an intuition for 64 why the temporal balance condition as described in the main text does not exactly 65 determine the stability of the network. For Eq (S12) to become less than zero, the 66 w∆qτ nmda term in the numerator of that equation must be negative and have its 67 absolute value greater than the τ e + τ nmda portion of the numerator. Both τ e and 68 τ nmda represent damping caused by both the cell bodies of the neurons and the passing 69 of the recurrent activity through the synapses. For the network to become unstable the 70 transient imbalances caused by the w∆qτ nmda term must drive the system sufficiently 71 to counteract the inherent damping in the cell body and synapses.  The full linear network is defined in Eqs (2)-(4) in the main text. The steady states of 74 the network can be found by setting the derivatives to zero and then solving for R e , R i 75 and S l mn . The steady state rates are, , (S13) As long as the balance conditions are met, the steady state equations for this network 77 have only one solution and if the input I is zero the steady state is also zero. These The roots of the characteristic polynomial of the balanced network 87 s n + a n−1 s n−1 + . . . + a 1 s + a 0 = 0. (S14) yield the eigenvalues. Using the coefficients of Eq (S14) it can be shown that the first 88 eigenvalue is small and negative when the ratio a 1 /a 0 is large and positive, yielding a 89 large network time constant τ n [3]. Conversely, when a 1 /a 0 is small and positive then 90 the first eigenvalue will be large and negative, resulting in a short τ n . We do not provide 91 sufficient conditions for small τ n , only a necessary condition. However, simulations show 92 that when a 1 /a 0 is small and positive, τ n is also small, Fig 2D. We assume that J is The ratio a 1 /a 0 is large, thereby producing a large τ n , if the following three conditions 98 are met, These three conditions are taken directly from terms in Eqs (S15) and (S16) using 100 the requirement that a 1 /a 0 is large. The first constraint, Eq (S17), states that positive 101 feedback (second term) and negative feedback (first term) should be balanced. This 102 "balance condition" is also a requirement for the system to be stable as we will see in a 103 later section. The other two constraints describe the "temporal balance condition," first 104 introduced in this contribution. Formula (S18) describes a condition which depends on 105 the overall time constants of the EE and IE connections. It produces changes in the inequality (S18) reduces to the balance condition, inequality (S17), with a constant 110 coefficient. This term has been described previously for networks with only one time 111 constant on each projection [3]. We assume that synapses of the same neurotransmitter 112 type have the same time constant on all projections, so τ ampa is the same for all 113 projections with AMPA currents, and the same applies to τ nmda and τ gaba . Therefore, 114 Eq (S18) reduces to Eq (S17) and has little impact on τ n , leaving only Eq (S19).

115
Eq (S19) describes the impact of changing the ratio of two time constants on the EE 116 and IE projections which is how our network produces large τ n as ∆q increases.

117
Fundamentally, this constraint states that changing the strength of AMPA and NMDA 118 currents on the EE or IE projections will change the average time constant of that 119 projection and thereby alter τ n . In a mathematical sense, the cross multiplication of 120 synaptic time constants for one type of synapse with the synaptic strength of the other 121 type of synapse means that the first term in Eq (S19) can be made larger than the (S20) This change increases the ratio a 1 /a 0 yielding a longer time constant of decay. In 128 addition, as can be seen from Eq (S13), this type of balanced change in synaptic 129 strength will not change the equilibrium state of the network. In addition to finding the eigenvalues of the system, the characteristic polynomial can 132 also be used to determine the stability of the system. The Routh-Hurwitz stability 133 criterion provides necessary and sufficient conditions for the stability of a time invariant 134 linear system with constant coefficients. However, due to the complexity of the 135 inequalities this yields we do not require that the polynomials be Hurwitz, but only that 136 the coefficients of the characteristic polynomial are all positive which is a necessary but 137 not sufficient condition for stability. The stability conditions below work well for the 138 delta oscillatory instability when ∆q is small but miss the gamma oscillatory instability 139 for large and positive ∆q.

140
The stability criteria are, The first condition, eq. (S21), requires that negative feedback be greater than positive 141 feedback; this is the balance condition. We use it and the assumption that the time The strict inequality in Eq (S25) implies that such a system should be unstable.

159
However, in practice equality in the stability conditions also gives a stable network ))))) However, in the linear case, our balanced network can be shown to be equivalent to one 189 with different combinations of AMPA and NMDA receptors at its synapses. In addition, 190 we will show that when STD is added, the STD model we presented is a special case of 191 a more general model where each synapse has one type of STD and a ratio of AMPA to 192 NMDA receptors.

193
Without STD, we consider two synapses that are allowed to have different ratios of 194 NMDA to AMPA receptors. We use q again for the fraction of synaptic strength All variables are the same as previously described except that there are now four types 202 of synapses labeled A − D which have associated NMDA percentages q A − q D . If we 203 reorganize these equations we find that they are exactly the same as Eqs (2)- (3) where 204 the synaptic strengths are now defined such that, Therefore, the linear model without STD is equivalent to one where each synapse can 206 have a different ratio of NMDA to AMPA receptors.

207
When STD is added to Eqs (S39)-(S40) the coefficients are no longer equivalent.

208
Now the coefficients representing the synaptic strengths include STD and hence are not 209 constant, where x 1 and x 2 represent the STD equations defined for the balanced network with 211 STD. This new model will produce similar responses to the model with pure AMPA or 212 NMDA synapses. When q A → 0, q B → 1, q C → 0 and q D → 1, Eqs (S45)-(S48) are 213 equivalent to the model with STD presented in the main paper.
Each subscript corresponds to the respective power of R e in the polynomial. We will 221 refer to steady state values of R e as R ss e . Given that we only consider situations where 222 I > 0 then the coefficient a 0 must always be less than zero. Since coefficient a 3 is 223 positive there is at least one change of signs in the coefficients. Thus, Descartes' rule of 224 signs implies that there will always be at least one positive root for non-oscillatory 225 solutions. Therefore, the number of positive and negative roots will be determined by 226 coefficients a 1 and a 2 . The right hand side of each of these equations is similar to the 227 balance equation implying that both a 1 and a 2 will tend to switch signs together 228 maintaining the single sign change within the coefficients. Consequently, there will 229 generally be one positive root and two negative roots when I > 0 and the balance 230 condition is met. Another way of saying this is that the firing rate of the excitatory 231 population has only one steady state solution. This reduces finding the steady state of 232 the system to just finding the largest root of the polynomial.

233
When the stimulus is removed such that I = 0 then a 0 = 0 giving a zero root.

234
Therefore, as long as there are no positive roots, the system will always decay to zero. If 235 I = 0, the system will have one sign change and consequently one positive root when a 1 236 and a 2 are negative. This occurs when positive feedback is greater than negative 237 feedback as determined by the balance equation. Therefore, when recurrent excitation is 238 large compared to the negative feedback the system will decay to a non-negative steady 239 state. This analysis shows that the steady states of this system act exactly as we would 240 expect based on the linear portion of the network. A constant input causes the system 241 to move to one positive steady state. Removing the stimulus causes the system to decay 242 to zero unless there is unbalanced positive feedback. We will not list steady state Although changing w and k does cause some some changes in the quantitative response 264 characteristics of the network the overall qualitative response remains the same.

265
The resonant frequency at the AMPA instability is highly robust to changes in k and 266 w with a total change of approximately 1.4 Hz for a wide range of synaptic strengths, to the real axis control both the AMPA dominated instability and the increase of τ n . As 291 ∆q becomes increasingly negative, the two poles, visible in the right plot, separate from 292 the real axis then approach and cross the imaginary axis as the system becomes 293 unstable. This is what causes the delta oscillations and the subsequent instability. In 294 the other direction, as ∆q becomes large and positive, both poles move to the real line. 295 One pole then moves off to infinity and the other pole approaches the imaginary axis 296 again. As the second pole approaches the imaginary axis it increases τ n . Although this 297 pole gets very close to instability it never crosses the axis even as ∆q becomes large.

298
The instability at large ∆q is due to the two poles visible in the 60 Hz range in the left 299 plot. As ∆q becomes larger they approach and then cross the imaginary axis. This is 300 what causes the oscillatory instability in the gamma range.