Balanced networks of spiking neurons with spatially dependent recurrent connections

Networks of model neurons with balanced recurrent excitation and inhibition produce irregular and asynchronous spiking activity. We extend the analysis of balanced networks to include the known dependence of connection probability on the spatial separation between neurons. In the continuum limit we derive that stable, balanced firing rate solutions require that the spatial spread of external inputs be broader than that of recurrent excitation, which in turn must be broader than or equal to that of recurrent inhibition. For finite size networks we investigate the pattern forming dynamics arising when balanced conditions are not satisfied. The spatiotemporal dynamics of balanced networks offer new challenges in the statistical mechanics of complex systems.

The study of spatiotemporal dynamics and variability in complex systems is at the interface of the physical, chemical, biological, and social sciences [1,2].In the neurosciences, a longstanding topic of interest is the significant variability in cortical neuron spike train responses [3,4].Models of cortical networks capture this high variability when recurrent excitatory and inhibitory inputs are balanced.Such "balanced networks" show irregular and asynchronous spiking dynamics through a complex, sometimes chaotic, network state [5].Nevertheless, the statistics of balanced networks are amenable to mean field analysis [6][7][8][9][10], using techniques developed for spin-glass systems [11].Subsequent experiments in cortex lend support to balanced network states with measurements of large and opposing excitatory and inhibitory synaptic currents [12,13], asynchronous cortical activity [14], as well as the sensitivity of network dynamics to small perturbations [15].
The probability that two cortical neurons are connected depends on their separation in physical space or, for some sensory systems, feature space [16][17][18][19].There has been substantial theoretical work on the the spatiotemporal dynamics of phenomenological macroscale models of cortex [20,21].In contrast, theoretical work in balanced networks assumes a spatially homogeneous or discretely clustered topology [5,8,22].The capacity for pattern formation and spatial filtering in balanced networks with spatially dependent connection probabilities has not been addressed.
In this letter, we derive experimentally testable conditions on the strength and spatial profile of connection probabilities that must be satisfied for a recurrent network of excitatory and inhibitory neuron models to maintain a stable balanced state in the continuum limit.Specifically, we find that external inputs must be broader than recurrent excitation, which in turn must be broader than or equal in broadness to recurrent inhibition.Further, we investigate spatiotemporal spiking dynamics when stable balanced solutions do not exist.
Network model.We consider a network of N integrate-and-fire neurons, half of which are excitatory and half inhibitory, spaced evenly on the state space Γ = (0, 1], so that the kth excitatory or inhibitory neuron is at location x = 2k/N.The input current to the kth excitatory (α = e) and inhibitory (α = i) neuron is given by respectively, where x = 2k/N and s e, j (t) = ∑ i δ(t − t i e, j ) is the spike train of the jth excitatory neuron and similarly for s i, j (t).Static external input is provided by the terms J α (x).The synaptic weight, so that Γ has periodic boundaries and k αβ is the spatial profile of β to α connectivity.As in [5,6], we fix k αβ (x) 1 to assure asynchrony and we then consider the behavior of the network as N → ∞.
Cortical neurons receive a large number of high amplitude excitatory inputs, implying that a post-synaptic cell only requires only a fraction of excitatory presynaptic cells to drive spike responses [23].Following past studies in balanced networks [5,6,10] we model this with an O(1) distance between rest and spike threshold and consider To simplify calculations we define Under these scaling assumptions, a neuron receives recurrent input from O(N) excitatory neurons but only requires O( √ N) excitatory inputs to be active in an integration window to produce a spike.Finite firing rates are therefore only maintained in the continuum limit through a dynamically stable balance between excitation and inhibition [5,6,10].We next derive conditions under which such a stable balanced network state exists.
Conditions on the existence of a balanced state in the continuum limit.The mean firing rates of neurons in the network are denoted by ν α (x) = E[s α,k (t)] , where E[•] represents expectation over network connectivity and • the average over time.In the continuum limit, the mean input currents are related to the firing rates by for α = e, i where w αβ (x) = j αβ k αβ (x) and * denotes circular convolution on Γ .Similarly, the infinitesimal temporal variances of the input currents are given by We aim to derive conditions under which ν α (x), µ α (x) and V α (x) each converge to a finite limit as N → ∞ and ν α (x) does not become identically zero.For these conditions to be realized, we must have that Taking N → ∞ gives a Fredholm equation of the first kind whose solution, when it exists, is given in the Fourier domain by where f (n) = Γ e −2π xni f Γ (x)dx.This equality must hold at every Fourier mode, n, for which 0 at some Fourier mode, then for a solution to exist, it must also be true that j e (n) w ii (n) Requiring firing rates to be non-negative and not identically zero implies that where f = f (0) = Γ f Γ (x)dx.Note that Eq. ( 6) is equivalent to a balance condition derived in [6] for spatially homogeneous networks.We show below that Eq. ( 6) leads to a stable balanced state for large N but Eq. ( 7) does not.The solution in Eq. ( 5) is only viable if ν α has a well-defined inverse Fourier transform, which requires at least that for α = e, i.We investigate this condition for specific examples below.Example with Gaussian connectivity -The analysis of the balanced state above is valid in the N → ∞ limit for a large class of neuron models [8].To find solutions at large but finite system size, we use a leaky integrateand-fire (LIF) model [40].Steady-state firing rates can be found numerically using Monte Carlo simulations of the full network or by searching for a fixed point [ν 0 e (x), ν 0 i (x)] that satisfies where φ(µ, V) relates input mean and variance to firing rate of the LIF model in the diffusion limit [41] and where µ 0 α (x) and V 0 α (x) are given in terms of ν 0 e (x) and ν 0 i (x) by Eqs. ( 2)-(3).Numerical solutions to Eq. ( 9) were used for the curves labeled "FP" in Figs.1-2.
For ease of exposition, we consider Gaussian shaped connectivity kernels and we assume that probability (but not strength) of a connection depends only on presynaptic cell type.In particular, we set w αβ (x) = w αβ w β (x) and j α (x) = j α j(x) where , satisfy j = w β = 1 for α, β ∈ {e, i}.In this case, the balance condition in Eq. ( 8) is satisfied only if σ o > σ e , σ i .Hence, external inputs must be spatially broader than recurrent connections for a balanced solution to exist.Under this condition, taking the inverse transform in Eq. ( 5) gives the balanced solutions where ν α = ν α (0) from Eq. ( 5).Note that the peaked shape of the firing rate profile from Eq. ( 10), though spatially filtered by recurrent activity, is inherited by the peaked shape of the inputs.Flat inputs (p = 0) lead to a flat firing rate profile (ν α (x) = ν α ).
When the balanced state exists [42], simulations show asynchronous and irregular spiking dynamics (Fig. 1a).The microscopic state of the network is highly sensitive to the deletion of a single spike (Fig. 1a,c), but sufficiently small perturbations of the membrane potentials do not cause a divergence of trajectories (not pictured).These findings are consistent with previous studies showing that balanced networks can exhibit "stable chaos" characterized by exponentially long transients and insensitivity to sufficiently small perturbations [9,[24][25][26][27].
The macroscopic dynamics, measured by the network firing rates, is stable to the deletion of spikes.The firing rates are given by fixed point of Eq. ( 9), which converges to the balanced fixed point given by Eq. ( 10) as the network size increases (Fig. 1b,d).The distribution of Pearson correlation coefficients between the spike counts of neighboring neurons is approximately Gaussian-shaped with a mean near zero despite the fact that neighboring neurons share more than 5% of their inputs (Fig. 1e), consistent with the network having reached a stable asynchronous state [10].
Spatially imbalanced networks -An O(c) deviation of the firing rates away from balance yields an O(c √ N) deviation of the mean input currents, but only an O(c) perturbation of the input variance, c.f. Eqs. ( 2)-(3).When mean input is large in magnitude, the firing rate transfer of an LIF neuron can be approximated as threshold-linear, motivating the following mean-field approximation to firing rate dynamics, Here Θ(•) is the Heaviside function, τ m is the characteristic timescale of the neurons, γ > 0 is the gain of the  9).Solid blue curves computed from full network simulations.Dashed black line computed from the linear approximation given in Eq. ( 12).
neuron [43] and µ α is related to ν β through Eq. ( 2) for α, β ∈ {e, i}.Eq. ( 11) can be solved for arbitrary N and will provide intuition for network solutions when condition Eq. ( 8) is violated.If Eq. ( 11) admits a fixed point with strictly positive firing rates, it is given in the Fourier domain by ν 0 e = j e + j e w ii − j i w ei 2 − w ee + w ii + w ei w ie − w ee w ii ν 0 i = j i + j e w ie − j i w ee 2 − w ee + w ii + w ei w ie − w ee w ii (12) where = 1/(γ √ N).If Eq. ( 8) is satisfied then the fixed point in Eq. ( 12) converges to the balanced solution in Eq. ( 5) as N → ∞.If Eq. ( 8) is violated (σ o < σ e , σ i ) then the higher spatial Fourier modes, and therefore peak firing rates, from Eq. ( 12) diverge as N → ∞ (Fig. 2).Eventually this growth of higher Fourier modes causes ν α (x) < 0 for some x (Fig. 2a-c), at which point Eq. ( 12) no longer reflects a fixed point solution to Eq. (11).
Stability of the balanced state -The balanced fixed point from Eq. ( 12) is stable for the mean-field model in Eq. (11) whenever has eigenvalues with negative real part or, equivalently, when  13) as a function of n with σ e = 0.02 in (a) and σ e = 0.05 in (b) (other parameters as in Fig. 1a).c,d) Spike rasters from simulations of the LIF network.e) Relative L 2 deviation of simulated firing rates from the fixed point determined by Eq. ( 9) for various values of σ e and N (see inset).Arrows along horizontal axis mark the smallest value of σ e /σ i at which some eigenvalue of A(n) has positive real part.
at each Fourier mode, n.For the Gaussian-shaped kernels described above, stability of the balanced state as N → ∞ under this approximation requires that w ee < w ii and σ e ≥ σ i are satisfied in addition to Eq. ( 6), but networks satisfying Eq. ( 7) do not satisfy Eqs. ( 14) for large N.The mean-field model predicts instabilities of the balanced state for full network simulations reasonably well (Fig. 3).In particular, when σ e is sufficiently smaller than σ i , A(n) has eigenvalues with positive real part and the balanced fixed point loses stability and different spatial pattern is produced (Fig. 3a,c,e).
Further, for σ e /σ i near the stability transition, network exhibits waves of activity, but the time-averaged firing rates remain close to the balanced fixed point (Fig. 3b,d,e).The direction that these waves travel depends on initial conditions even when the network remains fixed (data not shown), suggesting a symmetrybreaking multistability.This spatially coherent activity is not captured by the mean field model in Eq. ( 11) and its theoretical description is outside the scope of this study.Regardless, our analysis of the mean field approximation provides a useful explanation for why the balanced state becomes unstable when excitatory projections are too much narrower than inhibitory projections.
Discussion -By taking into account the spatial dependence of connection probabilities, we have derived new conditions for the existence and stability of balanced solutions.With Gaussian connectivity, the conditions are simply σ o > σ e ≥ σ i .Consistent with this conclusion, several studies have found that thalamocortical projections are generally broader than intracortical projections [28][29][30] and circuit measurements in cortical layer 4 show that excitation projects more broadly than inhibition [19].In contrast, many previous models rely on broad inhibition to sharpen tuning curves [31,32] and promote pattern formation [20,21].Our results refute the notion that dynamical mechanisms relying on such broad inhibition can coexist with a balanced state in the continuum limit.Nevertheless, Eq. ( 10) reveals that recurrent connections in our model sharpen tuning curves even when σ e ≥ σ i .
For simplicity, we used a one-dimensional singlelayer model with periodic boundary conditions.Our methods can easily be adapted to different spatial topologies.The analysis of a balanced network on the entire real line is identical to that given in Eqs. ( 2)-( 10) except that a continuous Fourier transform takes the place of the discrete transform.Similarly, if a twodimensional network is considered, an identical analysis with a two-dimensional Fourier transform yields analogous results.Our model should be interpreted as a model of a single cortical layer where input from other layers is accounted for by the external inputs, J α (x).Recurrent connections between layers can be represented explicitly by adding additional excitatory and inhibitory populations [33], suggesting a possible direction for future work.
Spatially extended stochastic neural field models are typically constructed by appending additive noise to a deterministic model [21,34], similar to the practice of augmenting reaction diffusion systems with additive or multiplicative noise [1].Analysis of neural field models driven by external stochastic forcing shows that the spatiotemporal structure of noise is a critical determinant of the ensuing stochastic dynamics [35][36][37].In balanced networks, variability arises naturally through internal mechanisms [5,6,8,22], so that assumptions about the structure of external stochastic forcing are not required.Thus, whereas the study of spatially distributed systems with external stochastic forcing show how pattern forming systems filter noise, balanced networks with spatial interactions offer an alternative framework where complex internal dynamics is the source, as opposed to filter, of spatiotemporal variability.Our work lays a theoretical foundation for studying such networks and shows that they can exhibit rich dynamics, suggesting several directions for future study.

4 FIG. 1 :
FIG.1: Balanced network dynamics.(Color online) (a) Raster plots from two Monte-Carlo simulations with identical initial states (blue and red dots respectively, N = 10 5 ).In one simulation (red dots), the first spike after 150 ms was skipped, revealing a sensitivity to perturbations.(b) Population firing rates ν e (x) and ν i (x) calculated from full network simulations (solid blue line, N = 2 × 10 6 ), the fixed point from Eq. (9) (dotted red line, N = 2 × 10 6 ), and the balanced solution from Eq. (10) (dashed black line).(c) The normalized L 2 deviation between the vector of membrane potentials of the two network simulations from (a) as a function of time elapsed since a spike was skipped.(d) Relative L 2 distance between the balanced state and the fixed point (dashed black line with squares), and between the fixed point and network simulations (solid blue line with circles).(e) Histogram of spike count correlations between neighboring neurons (N = 10 5 ; mean of 4.84 × 10 −4 and a standard deviation of 0.0711).Correlations computed by counting spikes over a 300 ms window in each neuron from 200 simulations of the same network realization with different initial conditions.

FIG. 2 :
FIG. 2: Loss of balanced state if external inputs are narrower than recurrent connections.(Color online) Firing rates of the excitatory population when external inputs are narrower than recurrent inputs for system sizes (a) N = 10 5 , (b) N = 7.5 × 10 5 , and (c) N = 5 × 10 6 .(d) Peak firing rate of the excitatory population as a function of system size.In all panels σ o = 0.1, σ e = σ i = 0.2, and other parameters are as in Fig. 1.Dotted red curves computed by numerically solving Eq. (9).Solid blue curves computed from full network simulations.Dashed black line computed from the linear approximation given in Eq. (12).

2 FIG. 3 :
FIG. 3:The balanced state is unstable if recurrent excitation is too narrow compared to inhibition.(Color online) a,b) Maximum eigenvalue of the matrix A(n) from Eq. (13) as a function of n with σ e = 0.02 in (a) and σ e = 0.05 in (b) (other parameters as in Fig.1a).c,d) Spike rasters from simulations of the LIF network.e) Relative L 2 deviation of simulated firing rates from the fixed point determined by Eq. (9) for various values of σ e and N (see inset).Arrows along horizontal axis mark the smallest value of σ e /σ i at which some eigenvalue of A(n) has positive real part.