Network Events on Multiple Space and Time Scales in Cultured Neural Networks and in a Stochastic Rate Model

Abstract

Cortical networks, in-vitro as well as in-vivo, can spontaneously generate a variety of collective dynamical events such as network spikes, UP and DOWN states, global oscillations, and avalanches. Though each of them has been variously recognized in previous works as expression of the excitability of the cortical tissue and the associated nonlinear dynamics, a unified picture of the determinant factors (dynamical and architectural) is desirable and not yet available. Progress has also been partially hindered by the use of a variety of statistical measures to define the network events of interest. We propose here a common probabilistic definition of network events that, applied to the firing activity of cultured neural networks, highlights the co-occurrence of network spikes, power-law distributed avalanches, and exponentially distributed 'quasi-orbits', which offer a third type of collective behavior. A rate model, including synaptic excitation and inhibition with no imposed topology, synaptic short-term depression, and finite-size noise, accounts for all these different, coexisting phenomena. We find that their emergence is largely regulated by the proximity to an oscillatory instability of the dynamics, where the non-linear excitable behavior leads to a self-amplification of activity fluctuations over a wide range of scales in space and time. In this sense, the cultured network dynamics is compatible with an excitation-inhibition balance corresponding to a slightly sub-critical regime. Finally, we propose and test a method to infer the characteristic time of the fatigue process, from the observed time course of the network's firing rate. Unlike the model, possessing a single fatigue mechanism, the cultured network appears to show multiple time scales, signalling the possible coexistence of different fatigue mechanisms.

Fig 1. Time course of the network firing rate.

Panel A: noisy mean-field simulations; panel B: ex-vivo data. Random large excursions of the firing rate (network spikes and quasi-orbits) are clearly visible in both cases.

Fig 2. Inter-network-spike interval (INSI) statistics in the noisy mean-field model, for varying levels of excitation (w_exc) and inhibition (w_inh).

Panel A: 〈INSI〉 (the scale is in seconds); panel B: coefficient of variation of INSI (CV_INSI). For high net excitation (bottom-right quadrant) short-term depression plays a determinant role in generating frequent and regular (low CV_INSI) NSs; for weak excitability (upper-left quadrant) random fluctuations are essential for the generation of rare, quasi-Poissonian NSs (CV_INSI ≃ 1). The solid lines are isolines of the real part ℜλ of the dominant eigenvalue of the mean-field dynamics’ Jacobian; white line: ℜλ = 0 Hz; red line: ℜλ = 3.5 Hz; black line: ℜλ = −3.5 Hz. Note how such lines roughly follow isolines of 〈INSI〉 and CV_INSI.

Fig 3. Stability analysis of the linearized dynamics captures most of the variability in the inter-network-spike interval (INSI) statistics.

〈INSI〉 (panel A) and CV_INSI (panel B) vs the real part ℜλ of the dominant eigenvalue of the Jacobian of the linearized dynamics, for two networks that are pointwise identical in the excitation-inhibition plane, except for their size (circles: 200 neurons, as in Fig 2; squares: 8000 neurons). The data points almost collapse on 1-D curves when plotted as functions of ℜλ, leading effectively to a “quasi one-dimensional” representation of the INSI statistics in the (w_exc, w_inh)-plane. The region in which the INSIs are neither regular (CV_INSI ∼ 0) nor completely random (CV_INSI ≃ 1), as typically observed in experimental data, shrinks for larger networks. The filled circles mark a null imaginary part ℑλ.

Fig 4. Algorithms for network events detection.

Panel A: total network activity from simulation (blue line) with detected NS/quasi-orbits (green line) and avalanches (red line). Four large events (green line) are visible; the first and third are statistically classified as network spikes; the other smaller two as quasi-orbits. Note how network spikes and quasi-orbits are typically included inside a single avalanche. Panel B: a zoom over 0.5 seconds of activity, with discretization time-step 0.25 ms, illustrates avalanches structure (red line).

Fig 5. A broad spectrum of synchronous network events: simulations vs ex-vivo data.

Panels A-C: experimental distributions of network events. Panels A and B: ∼40-minute recordings from a very long recording, for the same network; panel C: ∼1-hour recording from another cultured network. Panels D-F: distributions from simulations of networks corresponding to the points in Fig 2 ((w_exc, w_inh) = (0.82, 0.7), (w_exc, w_inh) = (0.82, 0.55), (w_exc, w_inh) = (0.88, 0.55)). The three networks of panels D-F have increasing levels of subcritical excitability. Note the logarithmic scale on the y-axis. The solid lines are fits of the theoretical distribution of event sizes, a sum of an exponential (for quasi-orbits) and a Gaussian (for NS) distribution (see Models and Analysis, Eq 18). The vertical lines mark the probabilistic threshold separating NS and quasi-orbits.

Fig 6. Avalanche size distribution: simulations vs ex-vivo data.

Panels A-C: mean-field simulations, with fixed inhibition w_inh = 1. and increasing excitation (w_exc = 0.9, 0.94, 1). The distributions are well fitted by power-laws; panel B and C clearly show the buildup of ‘bumps’ in the high-size tails, reflecting the increasing contribution from network spikes and quasi-orbits in that region of the distribution. Panels D-F from ex-vivo data, different ∼40-minute segments from one long recording; power-laws are again observed, although fitted exponents cover a smaller range; in panels E and F, bumps are visible, similar to model findings. The similarity between the theoretical and experimental distributions could reflect changes of excitatory/inhibitory balance in time in the experimental preparation. Since all the three simulations lay on the left of or just on the bifurcation line (white line in Fig 2), the shown results are compatible with the experimental network operating in a slightly sub-critical regime.

Fig 7. Slow time-scales inference procedure: test on simulation data.

Panel A: correlation between low-pass filtered network activity f (see Eq 19) and the size of the immediately subsequent network spike plotted against the time-scale τ* of the low-pass integrator (continuous line). The correlation presents a clear (negative) peak for an ‘optimal’ value ${\tau }_{\text{optim}}^{*}=0.58$ s of the low-pass integrator; such value is interpreted as the effective time-scale of the putative slow self-inhibitory mechanism underlying the statistics of network events—in this case, short-term synaptic depression (STD); as a reference, the dotted line marks the value computed for surrogate data (see text). Panel B: for each point in the (w_exc, w_inh)-plane (see Fig 2), ${\tau }_{\text{optim}}^{*}$ vs average network activity; the continuous line is the best fit of the theoretical expectation for STD’s effective time-scale (Eq 20); the fitted values for the STD parameters τ _STD and u _STD are consistent with the actual values used in simulation (τ _STD = 0.8 s, u _STD = 0.2).

Fig 8. Slow time-scales inference procedure on ex-vivo data.

Correlation between low-pass filtered network activity f (see Eq 19) and the size of the immediately subsequent network spike plotted against the time-scale τ* of the low-pass integrator for two experimental datasets (different periods—about 40 minutes each—in a long recording). The plot in panel A is qualitatively similar to the simulation result shown in panel A of Fig 7: a peak, although broader and of smaller maximum (absolute) value, is clearly identified and statistically significant (with respect to surrogate data, dotted line). Panel B shows two significant peaks in the correlation plot, a possible signature of two concurrently active fatigue processes, with time scales differing by roughly an order of magnitude. Panel A: same data as Fig 5, panel B.

Fig 9. Slow time-scales inference procedure on simulation data with STD and spike-frequency adaptation.

Correlation between low-pass filtered network activity f (see Eq 19) and the size of the immediately subsequent network spike plotted against the time-scale τ* of the low-pass integrator. In this case, the mean-field model includes, besides short-term depression (STD), a mechanism mimicking spike-frequency adaptation. Panel A: spike-frequency adaptation with characteristic time τ _SFA = 15 s. Panel B: τ _SFA = 30 s. In both cases the correlation presents a STD-related peak at around τ* ≃ 500 ms (τ _STD = 800 ms), consistently with Fig 7. The peaks at higher τ*s, found respectively at 12 and 22 s, in accordance with what is reported in the plot legends and in the main text, roughly preserve the ratio of the corresponding τ _SFA values.