- 1School of Mathematical Sciences, MOE-LSC, Institute of Natural Sciences, Shanghai Jiao Tong University, Shanghai, China
- 2Mathematical Sciences Department, Rensselaer Polytechnic Institute, Troy, NY, United States
- 3Courant Institute of Mathematical Sciences, Center for Neural Science, New York University, New York, NY, United States
- 4NYUAD Institute, New York University Abu Dhabi, Abu Dhabi, United Arab Emirates
Some previous studies have shown that chaotic dynamics in the balanced state, i.e., one with balanced excitatory and inhibitory inputs into cortical neurons, is the underlying mechanism for the irregularity of neural activity. In this work, we focus on networks of current-based integrate-and-fire neurons with delta-pulse coupling. While we show that the balanced state robustly persists in this system within a broad range of parameters, we mathematically prove that the largest Lyapunov exponent of this type of neuronal networks is negative. Therefore, the irregular firing activity can exist in the system without the chaotic dynamics. That is the irregularity of balanced neuronal networks need not arise from chaos.
1. Introduction
Neural spiking activity in the brain is highly irregular (Britten et al., 1993; Shadlen and Newsome, 1998; Compte et al., 2003; London et al., 2010). It is believed that the irregularity of the spiking activity can reflect an underlying rich coding structure for information processing (Hertz and Prügel-Bennett, 1996; Gütig and Sompolinsky, 2006; Sussillo and Abbott, 2009; Monteforte and Wolf, 2012). This viewpoint naturally leads to the investigation of the origin of the irregular neuronal activity. A number of theoretical studies have postulated that a balance between excitatory and inhibitory inputs into an individual neuron can give rise to irregular activity (van Vreeswijk and Sompolinsky, 1996; Troyer and Miller, 1997; Vreeswijk and Sompolinsky, 1998; Vogels and Abbott, 2005; Miura et al., 2007). The idea behind the theory of balanced networks is that the excitatory and inhibitory components of inputs nearly cancel each other, and the neuronal firing activity is driven by strong fluctuations that intermittently interrupt this cancellation. Consistent with the hypothesized scenario, balanced synaptic inputs have been observed in slices of the ferret prefrontal and occipital cortex (Shu et al., 2003). Moreover, it has also been found that, in vivo studies, the balanced excitation and inhibition in ferrets' prefrontal cortex can substantially influence the neuronal activity (Haider et al., 2006).
It has been shown theoretically that small perturbations of the balanced state in a network with binary neurons grow exponentially, indicating the chaotic nature of the balanced activity (Vreeswijk and Sompolinsky, 1998). Some studies suggest that highly irregular activity in the balanced state originates from chaotic network dynamics (Vogels et al., 2005; Wallace et al., 2013; Ostojic, 2014). Albeit not specifically addressing the source of irregular activity in a balanced state, studies exist that have drawn the opposite conclusions in some special situations. For example, numerical simulations of neural networks consisting of pulse-coupled spiking neurons of only inhibitory type can display irregular activity in a dynamical state with a negative Lyapunov exponent (Zillmer et al., 2006; Jahnke et al., 2008; Monteforte and Wolf, 2012). Meanwhile, in the limit of fast synaptic response, any generic trajectory was shown to be asymptotically stable in inhibition-dominated networks (only a small fraction of connections can be excitatory) with the inhomogeneous delay distribution and strong coupling (Jahnke et al., 2009). Therefore, the question of whether the irregularity in a balanced state arises from chaos in a neuronal network with both excitatory and inhibitory neurons under more realistic Poisson drive, remains an important issue to be further clarified.
In this work, we first numerically show that, over a broad range of parameters, the balanced state can exist in current-based integrate-and-fire (I&F) neuronal networks consisting of both excitatory and inhibitory neurons with delta-pulse coupling currents and pulse-like external inputs. We then mathematically prove that, driven by any point process in time—not limited to Poisson point processes, the current-based I&F neuronal networks with delta-pulse interactions cannot exhibit chaotic dynamics. In fact, two nearby trajectories of such a network generically coalesce after a finite time regardless of whether the dynamics is in a balanced state or not. Our results demonstrate that in the delta-pulse coupled, current-based I&F system the irregular activity of the balanced state is not a consequence of a chaotic dynamical state. Our proof remains valid when the system possesses only excitatory or inhibitory population. This conclusion extends the previous results (Jin, 2002; Zillmer et al., 2006; Jahnke et al., 2009) that dynamics of the strongly inhibition-dominated networks are stable. By our analysis, stable dynamics with the irregular firing activity can occur in the delta-pulse coupled, current-based I&F neuronal network of any size.
2. Materials and Methods
2.1. I&F Model
We model neurons as integrate-and-fire (I&F) units (Dayan and Abbott, 2001; Newhall et al., 2010; Zhou et al., 2010). The governing equation for the membrane potential of the ith neuron in the kth population is
where gL denotes the leakage conductance, is the resting voltage of the kth population and is the corresponding input current (k = E, I). The voltage evolves according to Equation (1) when , where is the threshold of the kth population. When crosses , the neuron spikes, and then is reset to . Upon resetting, is immediately governed by Equation (1) again. In simulation, corresponds to the membrane time constant 20 ms. The dimensionless values used in simulations are , and , which correspond to the parameters in Vreeswijk and Sompolinsky (1998).
The instantaneous current projecting into the ith neuron in the kth population,
consists of two terms, where is the excitatory input, is the inhibitory input, δ(·) is the Dirac delta function, fk is the strength of the external input, and Jkl is the coupling strength from the lth population to the kth population, k, l = E, I. The coupling constant is an element of the adjacency matrix of the network. It describes the connection from the jth neuron in the lth population to the ith neuron in the kth population. The first term in corresponds to the current arriving from the external input. The spike-time sequence, , corresponds to the external input into the ith neuron in the kth population. At the arrival time of the sth spike, the voltage of the ith neuron in the kth population jumps by the amount of fk. In the simulations below, we use Poisson trains for the external inputs. The second term in and the term in correspond to the presynaptic input from neurons of the excitatory and inhibitory populations in the network, respectively, where is the arrival time of the sth spike from the jth neuron in the kth population for k = E, I.
In our simulation, each neuron in the network has, on average, K excitatory and K inhibitory presynaptic neurons. Because each neuron may receive a large number of synapses in the cortex (Peters, 1987; Braitenberg and Schuz, 1998), and the connection between cortical neurons is often sparse with low connection probability (Holmgren et al., 2003), we choose K to be sufficiently large but much smaller than the total number of neurons in the network. Experimentally, for example, it is observed that cells in the primary visual cortex of adult cats fire much more irregularly than cells in vitro when they are both stimulated by injecting direct current through the electrode. Therefore, there is a substantial influence of fluctuations of synaptic inputs on the irregular activity (Holt et al., 1996). To capture the effect of fluctuations as observed in the experiment, we follow the balanced network theory to set the scaling of the coupling strength Jkl to be of order , thus leading to the scaling of fluctuations in the total synaptic inputs as of order 1. As a consequence, the fluctuations persist in the large-K limit (van Vreeswijk and Sompolinsky, 1996; Vreeswijk and Sompolinsky, 1998; Vogels and Abbott, 2005).
The connection from the jth neuron in the lth population to the ith neuron in the kth population in our simulation follows a Bernoulli distribution, i.e., the probability and , where Nl is the total number of the lth population, k, l = E, I. The parameter values in our simulations are as follows: NE = 32000, NI = 8000, K = 400, , , , , , and the external inputs are Poisson processes with rate νE = νI = ν0K, where ν0 controls the magnitude of Poisson rate.
Simulations of the neuronal network model are carried out to the machine accuracy using the event-driven algorithm (Brette et al., 2007). The event-driven algorithm for a general point process as external inputs proceeds by generating the time of the next external spike. Some neurons in the network receive this spike. Upon receiving an external spike, these neurons' voltages increase by the amount of fE for the excitatory population and fI for the inhibitory population. If some of them have their voltages exceed the threshold, they fire. Their voltages are held at the reset voltage, and the voltages of their postsynaptic neurons are instantaneously increased (for excitatory inputs) or decreased (for inhibitory inputs). It is possible that the voltages of these postsynaptic neurons may now be above threshold as well, then, these neurons also fire. Their voltages are held at the reset voltage, while their postsynaptic neurons' voltages are changed. This process repeats until no new neurons spike. We emphasize that in our dynamics we hold the voltage of the neurons that just fired at the resting potential in order to prevent any of these neurons from firing more than once at any given time. After all the neurons are updated at this time, we release these neurons from the reset voltage to follow the dynamics governed by Equation (1) until the next external spike.
2.2. Analysis of the I&F Network
As is well known, for the balanced state in binary neuronal networks, the population-averaged firing rate is a linear function of the external input. We now address the question of whether the balanced state in the I&F network also possesses the linear response property of the population-averaged firing rate to the external drive.
We consider the mean firing rate mE and mI in the large-K limit. In a balanced network, the firing events of different neurons are nearly independent of one another (Vreeswijk and Sompolinsky, 1998). This is illustrated in Figure 1 in which the cross-correlation between firing events of pairs of neurons is narrowly distributed around zero. According to Equations (1) and (2), it is obvious that, under Poisson-train external inputs, spiking events of a neuron in the network, in general, are not a Poisson train, i.e., and do not follow Poisson point process for a fixed j in Equation (2). However, the input to the ith neuron is a spike train summed over outputs from many neurons in the network. Since the firing events of neurons are statistically independent from one another, the summed spike train from a large number of output spike trains of neurons in the network asymptotically approaches a Poisson spike process (Cinlar, 1972). The recurrent excitatory (inhibitory) input of each neuron can be treated as a Poisson train with rate KmE (KmI), where mE (mI) is the mean firing rate per neuron averaged over the excitatory (inhibitory) population. Each neuron receives three Poisson spike trains in Equation (2). By homogeneity of our networks, the subscript i in Equation (2) will be dropped for the remainder of this discussion.
Figure 1. Defining characteristics of a balanced network. (A) Balanced excitatory and inhibitory inputs into a sample I&F neuron. The amplitude of excitatory (red) and inhibitory (blue) inputs (normalized by gL) deviate significantly from the firing threshold (green solid line) whereas the summed total input (scaled by gL) (black) fluctuates around the threshold while keeping its mean (cyan) below the threshold; (B) The average fano-factor across the network as a function of the time bin size (black line). The values of fano-factor for different time bin sizes (larger than 100 ms) are above unity and nearly constant. The inset is the distribution of the fano-factor calculated from each neurons spike events over the entire network with the time bin size 400 ms. All the values are above unity, again exemplifying the irregularity of the spiking activity; (C) Distribution of neuronal firing rates normalized by the entire population-averaged firing rate. The distribution is broad and skewed as a consequence of the heterogeneity of the single neuron activity; (D) The distribution of cross-correlation coefficient of spike events for each pair of neurons in the network. It has a sharp peak around zero, therefore, the system is nearly uncorrelated. We choose the bin size Δt = 2ms to calculate the cross-correlation; (E) The upper panel is a raster plot of 500 excitatory (red dots) and 500 inhibitory (blue dots) neurons; the lower panel shows the number of the spiking neurons in the excitatory population (red) and inhibitory population (blue) in each time window. The raster plot shows there is no synchrony in the neuronal dynamics. The lower panel exhibits stationarity of the number of firing neurons in the network; (F) The mean firing rate of the excitatory and inhibitory population as a linear function of the external driving rate. For comparison, we vary the value of K here: K = 100 (dashed-dot), K = 400 (solid lines), and K = 3600 (dashed lines); red for excitatory population whereas blue for inhibitory population. The mean firing rates in the K → +∞ limit, and , are also shown (black line, the theoretical gain curves for the excitatory and inhibitory populations overlap). Note that the minor disagreement between the theoretical prediction and the simulated gain curve arises from the finite size effect, and the gain curve converges to the infinite-K limit as K increases. The parameters are NE = 32000, NI = 8000, , , , , , νE = νI = ν0K and ν0 = 30Hz. In (A–E), K = 400.
Since the voltage of neuron in the kth population is reset to after spiking, we consider Equation (1) with initial value for k = E, I. We can readily obtain with , where Mk,l is the total spike count of the lth input to the kth population, described by a Poisson distribution with the average number of successes νk,lt, k = E, I and l = ext, E, I. For each given Mk,l, the spike times of the lth input to the kth population, , s = 1, 2, …, are uniformly distributed on the interval [0, t]. Clearly, the random variable takes value in the interval with the following probability density for . Since vk,l(t) is given by a sum of independent identically distributed random variables for given Mk,l, the average neuronal voltage at time t can be simply expressed as
for k = E, I. As discussed above, , fk, JkE, and JkI are of order ; νk is of order K; and mE and mI are of order 1. Then, the leading order of uk(t) is . Obviously, the membrane potential cannot become infinite as K → +∞. Therefore, the leading order of should vanish and one can obtain
in the large-K limit. Equation (4) describes the linear dependence of the mean firing rate on the external input. This relation is a defining feature of a balanced state, similar to what is obtained for the binary neuronal system (Vreeswijk and Sompolinsky, 1998).
3. Results
3.1. Existence of Balanced States
In Vreeswijk and Sompolinsky (1998), the properties of the balanced state are shown in detail with binary neuronal networks. We use numerical simulations to investigate whether the current-based I&F neuronal network coupled with delta-pulse interactions can exhibit the dynamical characteristics of a balanced state. Our results demonstrate that there indeed exists a balanced state in I&F neuronal networks. Figure 1 summarizes the defining characteristics of the balanced state as exhibited in I&F neuronal networks: balanced inputs (Figure 1A), irregular spiking activity (Figure 1B), heterogeneous firing rate (Figure 1C), weak correlation (Figure 1D), stationary asynchronous dynamics (Figure 1E), and linear response (Figure 1F). The results above confirm the important properties of balanced networks as discussed in previous theoretical work (van Vreeswijk and Sompolinsky, 1996; Vreeswijk and Sompolinsky, 1998; Renart et al., 2010; Litwin-Kumar and Doiron, 2012).
Next, we turn to the investigation of the persistence of the balanced state in the I&F system and derive conditions under which the balanced network state exists. Note that, in our simulation, we choose JIE = JEE. Therefore, from Equation (4), requiring firing rates to be non-negative implies
or
Note that Equation (5) is equivalent to the balance condition derived in Vreeswijk and Sompolinsky (1998) for a binary neuronal system. However, Equation (6) admits a solution for which mE = 0 — that is, only inhibitory neurons fire in the system, the mean firing rate of the inhibitory population is then given by mI = fIν0/JII. Therefore, for k = E in Equation (3), we can find , where C is a positive constant independent of K, which demonstrates that membrane potential is highly negative, and the excitatory neurons' firing activity is suppressed. In the simulation, we choose a range of parameters, in particular, fE, fI, JEI, and JII to examine the competition between the excitatory and inhibitory inputs to a neuron and verify whether the system maintains the characteristic balanced-state properties. If the neuronal network dynamics indeed exhibits the features displayed in Figure 1, that network is classified as balanced. Figure 2A summarizes our results for the parameters fE/fI and JEI/JII that we have scanned. Each dot in the figure can represent a set of parameters fE, fI, JEI, and JII with fixed ratios of fE/fI and JEI/JII. Each red dot in the parameter space indicates the fact that the system with the corresponding parameters can reach a balanced state; each blue dot represents the system with the corresponding parameters that exhibits synchronous dynamics; and each green dot indicates that at those parameter values only inhibitory neurons fire in the system but the inhibitory population still retains the characteristic balanced-state properties in Figure 1.
Figure 2. Dynamical network state as a function of fE/fI and JEI/JII. The range of parameters here is as follows: NE = 32000, NI = 8000, K = 400, 0 < fE < 0.2, 0.04 < fI < 1.0, νEfE = 10, νI = νE, , , . We choose values of the four parameters fE, fI, JEI and JII at random to verify whether the state of the network with these parameters is a balanced state. Note that each dot (C1, C2) in the panel (A) can represent a group of different parameter values as long as they satisfy and , where C1 and C2 are constant. A red dot represents a balanced state in which both excitatory and inhibitory neurons fire, with a sample raster in panel (B) (which has the parameters marked by the black asterisk in the area of red dots, and only 1000 excitatory and 1000 inhibitory neurons are shown here). A green dot represents a state in which only the inhibitory population is firing. This dynamical state also exhibits all the defining characteristics of a balanced state for the inhibitory population, while the excitatory population fires few spikes over a long time. A corresponding sample raster is shown in the panel (C) (which has the parameters marked by the black asterisk in the area of green dots, and only 1000 excitatory and 1000 inhibitory neurons are shown here). Finally, a blue dot represents a state of strong synchrony, for which a sample raster is shown in the panel (D) (which has the parameters marked by the black asterisk in the area of blue dots, and only 1000 excitatory and 1000 inhibitory neurons are shown here). The black dashed lines in the panel (A) are y = x, y = 1 and x = 1. In the raster plot, red dots denote the spike times of the excitatory population and blue dots denote the spike times of the inhibitory population.
The balanced state with only the inhibitory population firing is established by the balance between strong excitatory external input to the inhibitory neurons and strong inhibitory recurrent current. As shown in Figure 2C, in this state, only the inhibitory population exhibits asynchronous firing activity while the excitatory population is silenced due to relatively strong recurrent inhibition as compared with the inhibitory population.
The parameter values for the results of the representative case as presented in Figure 2B are marked by the black asterisk in the red-dot area in Figure 2A. The parameters for the representative case shown in Figure 2C are marked by the black asterisk in the green-dot area in Figure 2A. The broad region covered by the red and green dots demonstrates the existence of the balanced state over a wide range of parameter space. The parameter values marked by the black asterisk in the blue-dot area in Figure 2A with strong fE and JII in comparison to fI and JEI, respectively, render the system robustly synchronized as shown in Figure 2D.
Note that the area represented by Equation (5) is smaller than the red-dot region that supports stable balanced states in our numerical results. We also point out that the area described by Equation (6) contains a small number of stable balanced states in which both excitatory and inhibitory neurons are spiking in addition to balanced states in which excitatory neurons are almost silent while only inhibitory neurons are spiking. As a matter of fact, upon inserting mE = 0 and mI = fIν0/JII into Equation (3), and requiring that uE becomes sufficiently negative as K becomes large so the excitatory population rarely fires, we obtain the condition
which contains Equation (6) and covers a greater area than the parameter regions in Figure 2A filled with green dots. These differences between the numerical results and the theoretical conditions (Equations 5–7) arise from the finite size effect of K.
3.2. Absence of Chaos
As shown above, the I&F networks with delta-pulse coupling can persistently manifest the dynamics of a balanced state. We now address our central question of whether the irregular firing activity of neurons in the balanced network is a consequence of chaotic dynamics of the network. We can mathematically prove that the I&F networks with delta-pulse coupling cannot exhibit chaotic dynamics. Therefore, chaos may not underpin the irregular firing activity of neurons in a balanced state. Below is a proof that there is no chaos in the dynamics of the current-based I&F network coupled with delta-pulses and with a general point process (not limited to Poisson) as external inputs, the details of whose dynamics are described in the section Materials and Methods.
For each reference voltage trajectory, ,, , we consider the perturbed voltage trajectory , with sufficiently small perturbation size ϵ at initial time t = t0, i.e., ϵ = |ṽ(t0) − v(t0)| ≪ 1. The dynamics of the perturbed trajectory ṽ(t) is described by the equation
where is the sth spike of the jth neuron in the kth population along ṽ(t), and k = E, I. The external spike times are the same as those along the reference trajectory v(t). Then the largest Lyapunov exponent,
can be calculated for this I&F network dynamics (Zhou et al., 2009). As is well known, positive λmax measures the average exponential spreading of nearby trajectories, while negative λmax measures the exponential convergence of trajectories onto the attractor (Oseledec, 1968; Ott, 2002; Parker and Chua, 2012). Generically, an attractor is defined to be non-chaotic if λmax is non-positive. In what follows, we show that the largest Lyapunov exponent λmax of this I&F network is always negative and in fact approaches negative infinity for any spike train input.
3.2.1. Spike Train Sorting
For the preparation of the proof of the absence of chaos, we label all the neurons as {1, 2, …, NE + NI}, in which {1, 2, …, NE} labels the excitatory neurons and {NE + 1, NE + 2, …, NE + NI} labels the inhibitory neurons for easy description. That is τjs () is the sth spike of the jth neuron, where j = 1, 2, …, NE stands for an excitatory neuron and j = NE + 1, NE + 2, …, NE + NI stands for an inhibitory neuron.
Then for any fixed finite time T, we sort the spike times of the two trajectories {τpq} and into the increasing lists. Recall the fact that the voltage of both the reference and the perturbed trajectories of any neuron will cross the threshold only upon receiving spikes either from the external or excitatory recurrent input. Hence, the neurons fire precisely at the arrival time of the external or excitatory recurrent input spikes. There may be a group of neurons that fire at a particular time but with a distinct firing sequence amongst these neurons (see our firing dynamics described in the event driven algorithm in section Materials and Methods). Note that these simultaneously firing neurons only spike at the time when an external spike is received by some of the neurons within this group (possibly by all the neurons in the group) for the type of I&F neuronal networks with delta-pulse coupling. When we meet simultaneous firings during the sorting process, suppose these simultaneous firings are τp1q1, τp2q2, ···, τpaqa such that τp1q1 = τp2q2 = ··· = τpaqa, then we perform the following strategy to make the spike timing list unique.
1. p1 is chosen to be the smallest neuron label among the neurons that receive an external input spike in the simultaneously-firing group;
2. the remaining sequence of τp1q1, τp2q2, ···, τpaqa is sorted and reordered such that p2 < p3 < ··· < pa.
The same rule is applied in the perturbed network. Then we can obtain two unique increasing spike time sequences
where τprqr () denotes the th spike of the th neuron in the reference (perturbed) trajectory, stands for an excitatory neuron and NE + 2, …, NE + NI stands for an inhibitory neuron. M and are the total number of spikes in the reference and perturbed trajectories, respectively.
3.2.2. Proposition and Proof
We next turn to the following proposition,
PROPOSITION For each reference voltage trajectory described by Equation (1), and its perturbed voltage trajectory described by Equation (8), with their initial difference ϵ = |ṽ(t0) − v(t0)| at time t = t0, one can obtain the increasing spiking lists (10) for both trajectories over any fixed finite time T according to the sorting process described in the above. If the initial perturbation is sufficiently small, i.e., ϵ ≪ 1, then , , and for any r in the lists (10).
Note that, before the first spike of reference and perturbed trajectories, the time evolution of the perturbation can be obtained from the system of equations
That is, for a sufficiently small initial perturbation size ϵ, the distance between v(t) and ṽ(t) decays exponentially. Then we use mathematical induction to prove this proposition:
1. For the first firing event r = 1, either τp1q1 or is equal to some external input spike time (the first spike must be induced by the external input). Because the reference and perturbed trajectories receive the same external input, for a sufficiently small initial perturbation size ϵ, we have , and ;
2. Suppose we have , and for r ≤ m, we then show that , and ;
3. In fact, before the (m+1)th spike of both the reference and perturbed trajectories, we can find that δv(t) is still governed by Equation (11) because they receive the same external and recurrent input spikes. For a sufficiently small initial perturbation size ϵ, Equation (11) also ensures and . Then, if the qm+1th spike of the pm+1th neuron is caused by the arrival of an excitatory recurrent input spike, this spike can always be traced back to the spike that is caused directly by an external input spike. By the ordering rule of neurons in the simultaneously firing group above, we have τpm+1qm+1 = τpsqs and , where either the spike of τpsqs or the spike of is caused by an external spike to the psth neuron (s ≤ m). Consequently, . If the qm+1th spike of the pm+1th neuron is caused directly by the arrival of an external spike, clearly we have because the reference and perturbed trajectories receive the same external input spike;
4. Finally, we can readily obtain .
Therefore, holds for r = 1, 2, …, M with and . Next, according to the Equation (9), it is obvious that the largest Lyapunov exponent is always negative for the dynamics described by Equation (11). That is, the I&F network with delta-pulse coupling exhibits no chaotic dynamics while receiving pulse-like external input. Note that the voltage of both the reference and perturbed trajectories of any neuron, for example, the ith neuron in the kth population, will be reset to after its first firing. Because both the external input spike time and the synaptic input spike time are the same in the reference and perturbed networks, the two trajectories will no longer separate from each other. This means that the difference between the reference and perturbed trajectories of the ith neuron will converge after its first firing event. Thus, the largest Lyapunov exponent will approach negative infinity. As show in Figure 3, the numerical results for the given fixed network size N and mean degree connectivity K indeed show that the total perturbation of all neurons' voltages always exponentially decay with time and eventually converge to zero in a finite time.
Figure 3. Time evolution of the total perturbation of all neurons' voltages. The total perturbation of all neurons' voltages |δv(t)| in either region described in Figure 2A always exponentially decay with time and eventually converge to zero in a finite time. (A) The parameters are marked by the black asterisk in the area of red dots in Figure 2A; (B) The parameters are marked by the black asterisk in the area of green dots in Figure 2A; (C) The parameters are marked by the black asterisk in the area of blue dots in Figure 2A. The parameters are chosen to make the mean firing rate of spiking neurons around 30 Hz. The total perturbation of all neurons' voltages at initial time in each case is chosen to be |δv(0)| = 5 × 10−4. The slope of the solid part of the curve in each case is −gL. The vertical dashed part of each curve indicates that the difference between the reference and perturbed trajectories vanishes, namely, the largest Lyapunov exponent approaches negative infinity.
As a matter of fact, the proof can be extended to a more generalized system. For example, the membrane potential of the ith neuron in the kth population obeys equation with refractory periods
but the neuron will come into a refractory period after crosses the threshold . That is when , one has for , where is the refractory period of the ith neuron in the kth population. and are the excitatory and the inhibitory synaptic delay of the ith neuron, respectively. In this system, the leakage conductance, the threshold, the resting potential, the refractory period and the synaptic delay of one neuron can be different from that of another neuron. Note that a neuron in this system can generate a spike also only upon receiving spikes either from the external or excitatory recurrent input, and the neuron fires precisely either at the arrival time of the external spike or at the arrival time of the excitatory recurrent spike with a fixed synaptic delay. Therefore, we can extend the proof above directly here to conclude that there is no chaotic dynamics in this generalized system. Despite the demonstration of irregular firing of neurons in the balanced state of this current-based I&F network as above, the irregular activity is not induced by chaos.
4. Discussion
Cortical neurons often exhibit spiking dynamics that are highly irregular. It is believed that the irregular neural spiking activity can be generated from a balance between excitatory and inhibitory inputs to a neuron. In the early work (Vreeswijk and Sompolinsky, 1998), the neuronal network consisting of binary neurons was applied to explain the balanced neuronal state. Meanwhile, the dynamics of the binary neuronal system was found to be chaotic in the balanced state. Therefore, chaos was then often thought to be closely related to the irregular firing activity in a balanced network. Note that the binary neuron is a highly idealized model for describing a neuron.
Here we address the issue of chaotic origin of spiking irregularity using a current-based leaky integrate-and-fire (I&F) model with delta-pulse coupling and pulse-like external input. By examining the defining characteristics of a balanced state (Vreeswijk and Sompolinsky, 1998), and by exploring a wide range of parameter values, we have found that the balanced state persists robustly in this I&F system.
We then apply the analysis of the largest Lyapunov exponent to characterize the dynamics of this I&F system. We mathematically demonstrate that the largest Lyapunov exponent is always negative in the current-based I&F system with delta-pulse coupling and pulse-like external input, in which each neuron can have its own distinct values of leakage conductance, resting potential, spiking threshold, refractory period, and synaptic decay. The reference and perturbed realizations of a neuronal network trajectory can converge in finite time. It is worthwhile to point out that the non-chaotic property in our proof holds for neuronal networks of any size.
To understand under what conditions cortical-like firing irregularity can be generated, several efforts have been made. For instance, an early study (Renart et al., 2007) investigated the possible existence of multiple balanced steady states with persistent activity. It shows that the Poisson-like irregular spiking activity can arise from balanced regimes with sustained persistent activity. The neuronal network system used in the study (Renart et al., 2007) is the same as that in our work, but with different scaling of synaptic connection strengths, e.g., homogeneous or heterogeneous multicolumnar architecture. Note that the proof of non-chaotic mechanism in our system is independent of scaling structure of connection strengths as long as the synaptic interaction is delta-pulse coupled. Therefore, it is expected that the non-chaotic mechanism also works in persistent activity states.
In addition, another recent study found that the I&F neuronal network endowed with probabilistic synaptic transmission can underlie the Poisson-like spiking variability over a wide range of firing rates (Moreno-Bote, 2014). Here, we show that the neuronal network consisting of I&F neurons with delta-pulse couplings is always non-chaotic over a wide range of firing rates. It could be interesting to investigate whether our non-chaotic neuronal network with probabilistic synaptic transmission can achieve Poisson-like variability of the spiking activity as observed in cortex. For example, we can consider that each pair of connected neurons in our system has multiple delta-pulse synaptic contacts, and each contact is activated with certain probability.
As is shown here in our works, the irregular spiking dynamics is always stable in the current-based I&F model with delta-pulse coupling and pulse-like external input. Therefore, the irregular activity cannot simply arise from the underlying chaotic dynamics. We point out that the proof of the non-chaotic dynamics in our system relies on the facts that both the external input and the synaptic interactions are in the form of delta-pulse coupling. Therefore, for any fixed network size N and the mean degree connectivity K, the neuron fires precisely at the arrival time of the external or excitatory recurrent spikes. As a result, the difference between the reference and perturbed trajectories of each neuron will converge after its first firing event. Thus, the largest Lyapunov exponent will approach negative infinity and the dynamics is non-chaotic. As for the case of smooth synaptic coupling, e.g., the current in Equation (2) is an α-like function (Dayan and Abbott, 2001) with the rise time constant τr and the decay time constant τd as , it is expected that the non-chaotic mechanism still holds when the time constants τr and τd are both relatively small. However, if these time constants are not sufficiently small, e.g., greater than 2 ms, it is known that the network dynamics is not always stable and can be chaotic (Zhou et al., 2010; Harish and Hansel, 2015). In general, the phenomenon of chaos is model-dependent (Brette, 2004; Zhou et al., 2009; Sun et al., 2010) and could not be the ultimate source of irregularity in neuronal activity in the brain.
Author Contributions
QG, ZT, GK, DZ, and DC: Conceived and designed the research, Performed experiments and analyzed data, Wrote the paper.
Conflict of Interest Statement
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Acknowledgments
This work is supported by NYU Abu Dhabi Institute G1301 (QG, ZT, DZ, and DC); NSFC-11671259, NSFC-11722107, NSFC-91630208, and Shanghai Rising-Star Program-15QA1402600 (DZ); NSF DMS-1009575, NSFC-31571071 (DC); Shanghai 14JC1403800, 15JC1400104, and SJTU-UM Collaborative Research Program (DC and DZ). The remaining authors dedicate this paper to their late coauthor and mentor DC.
References
Braitenberg, V., and Schuz, A. (1998). “Comparison between synaptic and neuronal density,” in Cortex: Statistics and Geometry of Neuronal Connectivity (Berlin; Heidelberg: Springer), 37–38.
Brette, R. (2004). Dynamics of one-dimensional spiking neuron models. J. Math. Biol. 48, 38–56. doi: 10.1007/s00285-003-0223-9
Brette, R., Rudolph, M., Carnevale, T., Hines, M., Beeman, D., Bower, J. M., et al. (2007). Simulation of networks of spiking neurons: a review of tools and strategies. J. Comput. Neurosci. 23, 349–398. doi: 10.1007/s10827-007-0038-6
Britten, K. H., Shadlen, M. N., Newsome, W. T., and Movshon, J. A. (1993). Responses of neurons in macaque mt to stochastic motion signals. Vis. Neurosci. 10, 1157–1169. doi: 10.1017/S0952523800010269
Cinlar, E. (1972). “Superposition of point processes,” in Stochastic Point Processes: Statistical Analysis, Theory, and Applications, ed P. A. W. Lewis (New York, NY: John Wiley), 549–606.
Compte, A., Constantinidis, C., Tegnér, J., Raghavachari, S., Chafee, M. V., Goldman-Rakic, P. S., et al. (2003). Temporally irregular mnemonic persistent activity in prefrontal neurons of monkeys during a delayed response task. J. Neurophysiol. 90, 3441–3454. doi: 10.1152/jn.00949.2002
Gütig, R., and Sompolinsky, H. (2006). The tempotron: a neuron that learns spike timing–based decisions. Nat. Neurosci. 9, 420–428. doi: 10.1038/nn1643
Haider, B., Duque, A., Hasenstaub, A. R., and McCormick, D. A. (2006). Neocortical network activity in vivo is generated through a dynamic balance of excitation and inhibition. J. Neurosci. 26, 4535–4545. doi: 10.1523/JNEUROSCI.5297-05.2006
Harish, O., and Hansel, D. (2015). Asynchronous rate chaos in spiking neuronal circuits. PLoS Comput. Biol. 11:e1004266. doi: 10.1371/journal.pcbi.1004266
Hertz, J., and Prügel-Bennett, A. (1996). Learning short synfire chains by self-organization*. Network 7, 357–363. doi: 10.1088/0954-898X_7_2_017
Holmgren, C., Harkany, T., Svennenfors, B., and Zilberter, Y. (2003). Pyramidal cell communication within local networks in layer 2/3 of rat neocortex. J. Physiol. 551, 139–153. doi: 10.1113/jphysiol.2003.044784
Holt, G. R., Softky, W. R., Koch, C., and Douglas, R. J. (1996). Comparison of discharge variability in vitro and in vivo in cat visual cortex neurons. J. Neurophysiol. 75, 1806–1814. doi: 10.1152/jn.1996.75.5.1806
Jahnke, S., Memmesheimer, R.-M., and Timme, M. (2008). Stable irregular dynamics in complex neural networks. Phys. Rev. Lett. 100:048102. doi: 10.1103/PhysRevLett.100.048102
Jahnke, S., Memmesheimer, R.-M., and Timme, M. (2009). How chaotic is the balanced state? Front. Comput. Neurosci. 3:13. doi: 10.3389/neuro.10.013.2009
Jin, D. Z. (2002). Fast convergence of spike sequences to periodic patterns in recurrent networks. Phys. Rev. Lett. 89:208102. doi: 10.1103/PhysRevLett.89.208102
Litwin-Kumar, A., and Doiron, B. (2012). Slow dynamics and high variability in balanced cortical networks with clustered connections. Nat. Neurosci. 15, 1498–1505. doi: 10.1038/nn.3220
London, M., Roth, A., Beeren, L., Häusser, M., and Latham, P. E. (2010). Sensitivity to perturbations in vivo implies high noise and suggests rate coding in cortex. Nature 466, 123–127. doi: 10.1038/nature09086
Miura, K., Tsubo, Y., Okada, M., and Fukai, T. (2007). Balanced excitatory and inhibitory inputs to cortical neurons decouple firing irregularity from rate modulations. J. Neurosci. 27, 13802–13812. doi: 10.1523/JNEUROSCI.2452-07.2007
Monteforte, M., and Wolf, F. (2012). Dynamic flux tubes form reservoirs of stability in neuronal circuits. Phys. Rev. X 2:041007. doi: 10.1103/PhysRevX.2.041007
Moreno-Bote, R. (2014). Poisson-like spiking in circuits with probabilistic synapses. PLoS Comput. Biol. 10:e1003522. doi: 10.1371/journal.pcbi.1003522
Newhall, K. A., Kovacic, G., Kramer, P. R., Zhou, D., Rangan, A. V., and Cai, D. (2010). Dynamics of current-based, poisson driven, integrate-and-fire neuronal networks. Commun. Math. Sci. 8, 541–600. doi: 10.4310/CMS.2010.v8.n2.a12
Oseledec, V. I. (1968). A multiplicative ergodic theorem. lyapunov characteristic numbers for dynamical systems. Trans. Moscow Math. Soc. 19, 197–231.
Ostojic, S. (2014). Two types of asynchronous activity in networks of excitatory and inhibitory spiking neurons. Nat. Neurosci. 17, 594–600. doi: 10.1038/nn.3658
Parker, T. S., and Chua, L. (2012). Practical Numerical Algorithms for Chaotic Systems. Springer Science & Business Media.
Peters, A. (1987). “Number of neurons and synapses in primary visual cortex,” in Cerebral Cortex (Springer), 267–294.
Renart, A., De La Rocha, J., Bartho, P., Hollender, L., Parga, N., Reyes, A., et al. (2010). The asynchronous state in cortical circuits. Science 327, 587–590. doi: 10.1126/science.1179850
Renart, A., Moreno-Bote, R., Wang, X.-J., and Parga, N. (2007). Mean-driven and fluctuation-driven persistent activity in recurrent networks. Neural Comput. 19, 1–46. doi: 10.1162/neco.2007.19.1.1
Shadlen, M. N., and Newsome, W. T. (1998). The variable discharge of cortical neurons: implications for connectivity, computation, and information coding. J. Neurosci. 18, 3870–3896. doi: 10.1523/JNEUROSCI.18-10-03870.1998
Shu, Y., Hasenstaub, A., and McCormick, D. A. (2003). Turning on and off recurrent balanced cortical activity. Nature 423, 288–293. doi: 10.1038/nature01616
Sun, Y., Zhou, D., Rangan, A. V., and Cai, D. (2010). Pseudo-lyapunov exponents and predictability of hodgkin-huxley neuronal network dynamics. J. Comput. Neurosci. 28, 247–266. doi: 10.1007/s10827-009-0202-2
Sussillo, D., and Abbott, L. F. (2009). Generating coherent patterns of activity from chaotic neural networks. Neuron 63, 544–557. doi: 10.1016/j.neuron.2009.07.018
Troyer, T. W., and Miller, K. D. (1997). Physiological gain leads to high isi variability in a simple model of a cortical regular spiking cell. Neural Comput. 9, 971–983. doi: 10.1162/neco.1997.9.5.971
van Vreeswijk, C., and Sompolinsky, H. (1996). Chaos in neuronal networks with balanced excitatory and inhibitory activity. Science 274, 1724–1726. doi: 10.1126/science.274.5293.1724
Vogels, T. P., and Abbott, L. F. (2005). Signal propagation and logic gating in networks of integrate-and-fire neurons. J. Neurosci. 25, 10786–10795. doi: 10.1523/JNEUROSCI.3508-05.2005
Vogels, T. P., Rajan, K., and Abbott, L. F. (2005). Neural network dynamics. Annu. Rev. Neurosci. 28, 357–376. doi: 10.1146/annurev.neuro.28.061604.135637
Vreeswijk, C. v., and Sompolinsky, H. (1998). Chaotic balanced state in a model of cortical circuits. Neural Comput. 10, 1321–1371. doi: 10.1162/089976698300017214
Wallace, E., Maei, H. R., and Latham, P. E. (2013). Randomly connected networks have short temporal memory. Neural Comput. 25, 1408–1439. doi: 10.1162/NECO_a_00449
Zhou, D., Rangan, A. V., Sun, Y., and Cai, D. (2009). Network-induced chaos in integrate-and-fire neuronal ensembles. Phys. Rev. E 80:031918. doi: 10.1103/PhysRevE.80.031918
Zhou, D., Sun, Y., Rangan, A. V., and Cai, D. (2010). Spectrum of lyapunov exponents of non-smooth dynamical systems of integrate-and-fire type. J. Comput. Neurosci. 28, 229–245. doi: 10.1007/s10827-009-0201-3
Keywords: balanced state, irregular activity, chaotic dynamics, delta-pulse coupling, largest Lyapunov exponent
Citation: Gu QL, Tian ZK, Kovačič G, Zhou D and Cai D (2018) The Dynamics of Balanced Spiking Neuronal Networks Under Poisson Drive Is Not Chaotic. Front. Comput. Neurosci. 12:47. doi: 10.3389/fncom.2018.00047
Received: 17 August 2017; Accepted: 30 May 2018;
Published: 28 June 2018.
Edited by:
David Hansel, Université Paris Descartes, FranceReviewed by:
Ruben Moreno-Bote, Universidad Pompeu Fabra, SpainMarc Timme, Technische Universität Dresden, Germany
Carl Van Vreeswijk, Centre National de la Recherche Scientifique (CNRS), France
Copyright © 2018 Gu, Tian, Kovačič, Zhou and Cai. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Douglas Zhou, zdz@sjtu.edu.cn