Power-law statistics and universal scaling in the absence of criticality
Jonathan Touboul1, 2, ∗ and Alain Destexhe3, 4
1
arXiv:1503.08033v1 [q-bio.NC] 27 Mar 2015
The Mathematical Neuroscience Laboratory, CIRB / Coll`ege de France (CNRS UMR 7241,
INSERM U1050, UPMC ED 158, MEMOLIFE PSL*)
2
MYCENAE Team, INRIA Paris
3
Unit for Neurosciences, Information and Complexity (UNIC), CNRS, Gif sur Yvette, France
4
The European Institute for Theoretical Neuroscience (EITN), Paris
(Dated: March 30, 2015)
Critical states are usually identified experimentally through power-law statistics or universal scaling functions. We show here that such features naturally emerge from networks in self-sustained
irregular regimes away from criticality. Power-law statistics are also seen when the units are replaced
by independent stochastic surrogates, and thus are not sufficient to establish criticality. We rather
suggest that these are universal features of large-scale networks when considered macroscopically.
These results put caution on the interpretation of scaling laws found in nature.
PACS numbers: 02.50.-r, 05.40.-a, 87.18.Tt, 87.19.ll, 87.19.lc, 87.19.lm
Criticality, whether self-organized (SOC) or the result
of fine tuning, was invoked as the fundamental origin of a
number of complex phenomena in natural systems, from
sandpile avalanches to forest fires [1, 2]. It is usually
identified by a typical scale invariance or power-law scaling of specific events, as evidenced in canonical systems,
such as the sandpile [3] or the Ising model [4]. The relationship between power-laws in nature and criticality
became popular after the seminal work of Bak, Tang and
Wiesenfeld on the Abelian sandpile model [3]: there, it is
proposed that systems can self-organize to remain poised
at their critical state and that this may be a universal
explanation for the ubiquity of power-law scalings in nature [1] Following this theory, a considerable number of
publications reported the identification of power-law relations and SOC in many complex systems (see [2, 3] for
reviews). At criticality, such systems show a number of
scale-invariant statistics including the size, duration and
shape of the avalanches. Here, we show that statistical mechanisms provide an alternative natural explanation for the emergence of power-law and scale invariant
statistics in large-scale stochastic systems, which is particularly relevant to neuronal networks.
In neuroscience, the theory that neuronal networks operate at a phase transition would not only be a mathematical curiosity, but could have important consequences
in neural coding: a scale-invariant neuronal spiking activity could reveal the presence of long-range correlations
in the system [5] and optimal encoding of information [6].
This is probably one of the reasons why the characterization of the statistics of collective spike patterns has been
the object of intense research and active debate [7–10].
The first empirical evidence that neuronal networks may
operate at criticality was provided by analysis of the activity in neuronal cultures in vitro, which display bursts
of activity separated by silences. Based on an analogy
with the sandpile model, these bursts were seen as “neuronal avalanches”, and were apparently distributed as a
power-law with a slope close to −3/2, consistent with
the distribution of sandpile avalanches [7]. These experiments used macroscopic measurements of firing events
through local field potentials recordings, and avalanches
were related to peaks in the signals. As this procedure
may lead to spurious power laws related to the stochastic
nature of the signal [10], in-depth exploration of in-vitro
firing of several units were performed in cultured neurons [11]. These researches confirmed the presence of
apparent power-law distributions of avalanche size and
duration with critical exponents (respectively, −3/2 and
−2), as well as the existence of a universal mean temporal
profiles of firing events collapsing under specific scaling
onto a single universal scaling function.
These statistics were sometimes interpreted as the hallmark of criticality in analogy with canonical statistical
physics models. And in addition, a few models of neuronal activity related scale-invariant statistics to criticality [11, 12]. Several models used a branching processes
analogy (that rigorously only applies to networks having
solely excitatory connections) and investigated how the
activity propagates throughout the network: in this case
when a spike induces the firing of more than one spike on
average, there is an explosion of activity in the network,
while if it induces less than one spike, the activity dies
out after a typical time, and it is only at criticality, when
a spike induces on average one spike, that one may obtain
scale-free statistics. But avalanches are more described
through the distribution of the periods of silence, and in
biological neuronal networks, the presence of inhibition is
fundamental in shaping the activity and can be a crucial
element to the termination of avalanches. In this Letter,
we show that networks with inhibition do not need to
be at criticality to have power-law avalanche statistics as
reported in experiments, but is rather a generic macroscopic feature of such networks related to fundamental
properties of large interacting systems with inhibition.
In order to observe scale-free avalanches, neuronal
networks must display periods of collective activity of
broadly distributed duration interspersed by periods of
2
(see raster plot in Fig. 1(B)) [15]. This sharply contrasts
with the SI state, in which we can define avalanches
and find that these events display the same statistics
assumed to reveal criticality in cultures [11]. Fig 1(A)
represents the raster plot with typical avalanches taking
place. Strikingly, both avalanche duration and avalanche
size show excellent fit to a power law, validated by Kolmogorov Smirnov test of maximum likelihood estimator [16]. Moreover, the exponents we find for the fits are
close from the critical exponents of the sandpile model,
as also found in cultures [11]. In particular, we find that
the size s of the avalanches (number of firings during one
avalanche) scales as s−τ with τ = 1.63 (Fig. 1(C)), close
to the theoretical value of 3/2, and to the neuronal culture scaling (1.7 reported in [11]). The distribution of the
duration t of the avalanches scales as t−α with α = 1.86
(Fig. 1(D)), close to the theoretical value of 2 in the sandpile model (1.9 in neuronal cultures [11]). Moreover, let
us note that the typical deviations from the power laws
seen in this spiking model are consistent with experimental data: they are rather qualitatively compatible
with a underlying sub-critical regime, as found experimentally [11]. As a result, the average number of neuron
firing scales very nicely as a power of duration with positive exponent γ = 1.23 (Fig. 1(E)), similar to what was
observed in neuronal cultures (with a similar exponent
1.3), and which is consistent with Sethna’s size-duration
universal scaling relationship [14]
silence: this is typical of the Synchronous Irregular (SI)
regimes (see e.g. [13]), known to reproduce the qualitative
features of spiking in anesthetized animals or cultures.
Although partially ordered and partially disordered, SI
regimes occur in a wide range of parameter values (all in
inhibition-dominated regimes) [13] and are not very sensitive to modifications of biophysical parameters. It is
distinct from the Asynchronous Irregular (AI) regime, in
which neurons fire in a Poisson-like manner (and no period of silence occurs), typical of the awake activity and
where no avalanche can be found. It is also distinct from
the Synchronous Regular Regime (SR), in which collective activity is synchronized and quasi-periodic. In the
SI regime the inhibition dominates the excitation, and
when the activity invades the network, it will naturally
die out due to inhibition. This is what we observe in the
simulations in Fig. 1 in a sparsely-connected network of
excitatory and inhibitory integrate-and-fire neurons [13].
We investigate the statistics of spike units in both cases
(see Fig. 1).
(A)
1000
500
500
0
0
Rate
Neuron index
(B)
1000
1
1000
3
2
1200
4
1400
5
1600
1800
2000
6
2200
7
5
Time (s)
6
10
5
4
4
10
s−1.63
3
10
3
10
2
10
1
−1.86
1
10
10
0
10 0
10
t
2
10
0
1
10
2
10 0
10
3
10
10
Avalanche size
1
10
(F)
1.9
4
〈s〉(t;T)
2
10
3.5
3
2.5
1
10
−1.226
~T
0
10 0
10
Avalanche duration (# bins)
〈s〉(t;T) / T0.23
4.5
α−1
= γ.
τ −1
2
10
10
8.5
10
5
10
7.5
7
(E)
6
10
6.5
Time (s)
Mean Avalanche size
Event Frequency
(C)
6
5.5
(D)
1
10
Avalanche duration (# bins)
(G)
1.8
1.7
1.6
1.5
1.4
1.3
1.2
1.1
2
1
1.5
0
5
10
15
Time t (#of bins)
20
25
0.9
0
0.2
0.4
0.6
0.8
1
t/T
FIG. 1. Avalanche spike statistics in the sparsely connected
integrate-and-fire model [13] with N = 5 000 neurons in the
SI and AI states. (A) Raster plot of the SI together with the
firing rate (below). (B) Raster plot in the AI state: spiking
is asynchronous and faster (notice that the time window is
shorter compared to A for legibility): no silence period arises.
In the SI states, (C) avalanche size and (D) avalanche duration
scale as power law (dashed line is the maximum likelihood
fit), and averaged avalanche size scales as a power law with
the duration avalanche duration according to the universal
scaling law [14] (E). Average avalanche shapes collapse onto
the same curve (F,G) very accurately.
First, in the AI state, the sustained and irregular firing
does not leave room for repeated periods of quiescence at
this network scale, preventing the definition of avalanches
2
10
(1)
But one can go to an even finer scale and consider the averaged shape for an avalanche of fixed size. Typical profiles for three different durations are plotted in Fig. 1(F).
These shapes were shown to actually collapse onto a universal avalanche shape when time is normalized by the
total duration T and the mean number of active neurons is normalized by T γ−1 where γ is defined above.
The avalanches of Brunel’s model in the SI state show
very nicely this scale invariance with the particular coefficients found, as we show in Fig. 1(G) and consistent
with neuronal cultures findings.
These scalings are not specific to the particular choice
of parameters used in our simulations: we consistently
find, in the whole range of parameter values corresponding to the SI state of Brunel’s model, and in particular,
away from phase transitions, similar apparent power-law
scalings with similar scaling exponents. We conclude that
these statistics are not a hallmark of criticality in this
system, but rather correspond to the statistics of synchronous irregular states prominent both in anesthesia
and cultures in which power-laws have been found experimentally [7, 17].
The observation of such scaling relationship in noncritical simple models of neuronal networks suggests that
the observed scaling is due to other mechanisms. The
classical theory of the thermodynamics of interacting par-
ρ˙ t = −αρt + σξt
with (ξt ) a Gaussian white noise. This choice is interesting in that the distribution of excursions and their times
is known in closed form [20]. These distributions are, of
course, not heavy tailed: they have exponential tails with
exponent α. We investigated the collective statistics of
N = 2 000 independent realizations of Poisson processes
with this rate. The resulting raster plot is displayed in
Fig. 2(A). While the firing is an inhomogeneous Poisson process, macroscopic statistics show power-law distributions for the size (τ = 1.3, Fig. 2(B)) and for the
duration (α = 1.7, Fig. 2(C)), both statistically significant [16], and yielding the classical linear relationship
between average avalanche size and duration with a coefficient γ = 1.3 satisfying the classical relationship (1).
Averaged shapes of avalanches of a given duration are
again the same curves rescaled through the universal scaling rules of critical events (Fig. 2(D,E)). Similar results
were found using other choices of irregular rate functions
(not shown).
At this point, we conclude that power-law statistics
of avalanche duration and size, as well as collapse of
avalanche shape, may arise in non-critical systems, they
even appear in networks in which firing times are statistically independent. Therefore, they do not necessarily
reveal how neurons interact, but rather the presence of a
hidden dependence of all neurons on the same fluctuating function. A similar hypothesis was formulated as the
origin of Zipf’s law in nature in high dimensional data
with latent variables [21, 22]. Here, the situation is very
different in that we are not considering the frequency of
data, and moreover because neuronal networks are responsible of generating their own ‘latent variable’ that
naturally emerges from network interactions, while pairs
(of p-uplets) of cells become independent.
500
250
Rate
0
4
3
2
1
0
5
Time (s)
5
6
Event Frequency
10
10
5
10
4
10
~s
3
10
~t
4
10
1.3
1.7
3
10
2
10
2
10
1
10
1
10
0
0
10 0
10
1
2
10
3
10
4
10
10
10 0
10
5
10
2
3
10
4
10
10
2
8
s (t,T) / T
0.4
7
6
5
4
1.8
1.6
1.4
1.2
1
3
2
0
1
10
Avalanche duration (# bins)
Avalanche size
s (t,T)
ticle systems states that in large networks (such as those
of the brain), the correlations between neurons vanishes.
This is also known as Boltzmann molecular chaos hypothesis (Stoßzahlansatz ). In our context, in the limit of large
networks, neurons behave as independent jump processes
with a common rate, which is solution of an implicit
equation, as shown in various network models [18, 19].
This motivates the study of ensembles of neurons made
of statistically independent Poisson processes, with identical rates. The case of constant rates corresponds to AI
regimes, but to generate a stochastic surrogate of the SI
regime, we must consider non-homogeneous rates with
periods of silence. An obvious choice would be to replay
the rates extracted from the SI state, and indeed such
a surrogate generated power-law statistics (not shown),
but in this case we could not rule out whether the powerlaw statistics are encoded in the rate functions. To show
that this is not the case, we generated a surrogate independent of the rate functions, by using a common rate of
firing of the neurons given by the positive part ρ+
t of the
Ornstein-Uhlenbeck process:
Neuron index
3
0.8
5
10
15
20
25
30
Time t
35
0
0.2
0.4
0.6
0.8
1
Normalized time t/T
FIG. 2. Avalanches statistics in the independent Poisson
model with Ornstein-Uhlenbeck firing rate (α = σ = 1). Apparent power-law scalings, together with scale invariance of
avalanches shapes.
This being said, we are facing an apparent contradiction. Indeed, our theory provides an account for the presence of critical statistics in networks in the SI regime, but
discards AI states as possibly having critically distributed
avalanches. Evidences of critical statistics of avalanches
in-vivo in the awake brain have been reported, but are
more scarce and controversial. Unlike neuronal cultures,
the activity in the awake brain does not display bursts
separated by silences, but is sustained. Using a macroscopic measurement of neuronal activity, the Local Field
Potential (LFP), power laws could be shown from the distribution of peaked events [9]. However, the emergence
of this statistics may also be accounted for by the random nature of the signal and the thresholding procedure
used in this analysis, and moreover these power-laws may
not be statistically significant [10]. Double-exponential
distributions seem to provide a better fit [8]. In an
attempt to clarify this, we investigate whether if LFP
peaks distribution can distinguish between critical and
non-critical regimes. We considered the current-based
Vogels and Abbott model [23] which provides synaptic
currents which are more biologically plausible, and from
which we can calculate LFP signals VLF P . Indeed, these
signals are generated primarily by the synaptic currents
and depend on the distance between neurons and electrode, as derived from Coulomb’s law [24]:
VLF P =
Re X Ij
,
4π j rj
where VLF P is the electric potential at the electrode position, Re = 230 Ωcm is the extracellular resistivity of
4
brain tissue, Ij are the synaptic currents of and rj is the
distance between Ij and electrode position. Remarkably,
applying to this more sophisticate model the same procedures as in the original paper [9], we found that the
method cannot distinguish between structured or nonstructured activity: for a fixed firing rate and bin size,
we have been comparing in Figure 3 the LFP statistics
of the avalanche duration in a network of neurons in the
AI regime (A) or in the SI regime (C) that shows partial
order of the firing and saw no difference. We also compared the statistics to those of independent Poisson processes with constant rate (B): the three instances show
power-law scaling of avalanche duration, with the same
exponent, that seem to rather be related to the firing
rate and bin size than to the form of the network activity. The scaling coefficient is again, close from 3/2, and
varies with bin size.
AI
Poisson
SI
phenomenon. This is why we found the same scalings in
stochastic surrogates. This also justifies why power-law
scalings of avalanches were found in a model of spiking
cells when excitation balances inhibition [25], the only
regime in which the cells are in the SI state. But when it
comes to macroscopic measurements, it is very difficult to
conclude: both SI and AI regimes, as well as stochastic
surrogates, display power-law statistics in the distribution of LFP peaks. This shows that systems of weakly
correlated units, or their stochastic surrogates, can generate power-law statistics when considered macroscopically. Therefore, the power-law statistics of macroscopic
variables tells nothing about the critical or non-critical
nature of the underlying system. This potentially reconciles contradictory observations that macroscopic brain
variables display power-law scaling [9], while no sign of
such power-law scaling was found in the units [8, 26].
More generally, these results also put caution on the interpretation of power-law relations found in nature.
Acknowledgments
A.D. was supported by the CNRS and grants from the
European Community (BrainScales and Human Brain
Project). We warmly thank Quan Shi and Roberto
Zu˜
niga Valladares for interesting discussions and analyses.
3
10
2
10
1
10
0
10 0
10
1
10
2
10
0
10
1
10
2
10
0
10
1
10
2
10
FIG. 3. Avalanche analysis defined from the macroscopic variable VLF P of a network of integrate-and-fire neurons with exponential synapses. Networks displaying an AI state (left), its
stochastic surrogate (middle) or a SI state (right), all show
similar macroscopic power-law statistics with a universal exponent close to -3/2, and plots are indistinguishable.
We conclude that networks of neurons within synchronous irregular (SI) regimes may naturally display
power-law statistics and universal scaling functions due
to the combination of (i) an emergent irregular collective
activity in a large-scale system, and (ii) the power-law
statistics of such large-scale systems, due to molecular
chaos.] In agreement with this, the same statistics are
reproduced by a surrogate stochastic process where the
periods of firing and silences are themselves generated by
another stochastic process. This observation is evocative
of recent results showing that Zipf laws can arise from the
presence of hidden latent variable in high dimensional
systems [21, 22]. Here, we showed that neuronal networks can naturally adopt a regime in which neurons are
close to be independent, with similar irregular statistics.
In such an analogy, the common rate could be seen as
a latent variable, and both independence and irregularity build up only from the interactions between cells. In
networks with irregular rates, power-laws and data collapse are therefore as universal as the molecular chaos
∗
[email protected]
[1] P. Bak, How nature works (Oxford university press Oxford, 1997).
[2] H. J. Jensen, Self-organized criticality: emergent complex behavior in physical and biological systems, Vol. 10
(Cambridge university press, 1998).
[3] P. Bak, C. Tang, K. Wiesenfeld, et al., Physical review
letters 59, 381 (1987).
[4] M. Kardar, Statistical physics of fields (Cambridge University Press, 2007).
[5] W. L. Shew and D. Plenz, The neuroscientist 19, 88
(2013).
[6] D. R. Chialvo, Nature physics 2, 301 (2006).
[7] J. M. Beggs and D. Plenz, The Journal of neuroscience
23, 11167 (2003).
[8] N. Dehghani, N. G. Hatsopoulos, Z. D. Haga, R. A.
Parker, B. Greger, E. Halgren, S. S. Cash, and A. Destexhe, Frontiers in physiology 3 (2012).
[9] T. Petermann, T. C. Thiagarajan, M. A. Lebedev, M. A.
Nicolelis, D. R. Chialvo, and D. Plenz, Proceedings of
the National Academy of Sciences 106, 15921 (2009).
[10] J. Touboul and A. Destexhe, PloS one 5, e8982 (2010).
[11] N. Friedman, S. Ito, B. A. Brinkman, M. Shimono, R. L.
DeVille, K. A. Dahmen, J. M. Beggs, and T. C. Butler,
Physical review letters 108, 208102 (2012).
[12] A. Levina, J. M. Herrmann, and T. Geisel, Physical
review letters 102, 118110 (2009).
[13] N. Brunel, Journal of computational neuroscience 8, 183
(2000).
5
[14] J. P. Sethna, K. A. Dahmen, and C. R. Myers, Nature
410, 242 (2001).
[15] Sub-sampling of the network may lead to find periods of
quiescence, whose statistics will depend on the network
size. However, these are not robust statistical quantities
to define a critical regime.
[16] A. Clauset, C. R. Shalizi, and M. E. Newman, SIAM
review 51, 661 (2009).
[17] G. Hahn, T. Petermann, M. N. Havenith, S. Yu,
W. Singer, D. Plenz, and D. Nikoli´c, Journal of neurophysiology 104, 3312 (2010).
[18] P. Robert and J. D. Touboul, arXiv preprint
arXiv:1410.4072 (2014).
[19] A. Sznitman, Ecole d’Et´e de Probabilit´es de Saint-Flour
XIX , 165 (1989).
[20] L. Alili, P. Patie, and J. L. Pedersen, Stochastic Models
21, 967 (2005).
[21] D. J. Schwab, I. Nemenman, and P. Mehta, Physical
review letters 113, 068102 (2014).
[22] L. Aitchison, N. Corradi, and P. E. Latham, arXiv
preprint arXiv:1407.7135 (2014).
[23] T. P. Vogels and L. F. Abbott, The Journal of neuroscience 25, 10786 (2005).
[24] P. L. Nunez and R. Srinivasan, Electric fields of the brain:
the neurophysics of EEG (Oxford university press, 2006).
[25] M. Benayoun, J. D. Cowan, W. van Drongelen, and
E. Wallace, PLoS computational biology 6, e1000846
(2010).
[26] C. Bedard, H. Kroeger, and A. Destexhe, Physical review
letters 97, 118102 (2006).