L3 · Chapter 6

Neural avalanches from branching spikes, seen through a partial recording

Sandbox Physics · 2026-10-31

The chapters so far looked at one patch of membrane. Now connect many neurons: one spike makes a few neighbours fire in the next step, and they trigger the next batch. In a recording these cascades show up as bursts of activity called neuronal avalanches. This chapter first uses a branching process to work out how large the cascades are, then asks how the statistics of a recording drift away from the network itself when an electrode array sees only part of it.

How many spikes does one spike trigger

The Lab network is a 40 by 25 sheet of 1000 units. Each unit connects to 10 neighbours chosen at random within 5 sites. When a unit fires, each of its neighbours fires in the next step with probability $\sigma/10$ , and a unit that has just fired cannot fire in the next step. So the number of direct offspring of one spike follows the binomial law $\mathrm{Bin}(10,\sigma/10)$ with mean $\sigma$ , called the branching ratio.

If every generation reproduces independently in this way, generation $k$ has mean size $\sigma^k$ and the whole cascade has a mean size given by a geometric series. For $\sigma<1$ every cascade dies out; for $\sigma>1$ some grow for ever; $\sigma=1$ sits exactly between the two. That is the critical point.

We tested this by kicking one unit of a quiet network, with no external input, 20000 times. At $\sigma=0.5$ the first generation has mean $0.4986\pm0.0049$ ; 0.6006 of kicks give only the kicked spike, against the exact $0.95^{10}=0.5987$ ; the mean cascade size is $1.989\pm0.014$ against 2. At $\sigma=0.8$ the second generation has mean $0.630\pm0.007$ , below the 0.64 of independent reproduction, and the mean size is $4.83\pm0.06$ , below 5.

The shortfall comes from the network itself. Two spikes that hit the same neighbour in the same step make it fire only once, and a unit that has just fired cannot fire again at once. The first generation comes from 10 distinct neighbours of a single unit and cannot collide, so it matches the binomial law exactly; later generations collide more and fall short of independent reproduction.

\langle A_k\rangle=\sigma^{k},\qquad \langle s\rangle=\sum_{k\ge0}\sigma^{k}=\frac{1}{1-\sigma}\quad(\sigma<1)

The mean number of spikes in generation

k

and the mean cascade size

\langle s\rangle

, for a cascade that starts from one spike and whose generations reproduce independently.

Source: S. Zapperi · J. Wilting and V. Priesemann

At the critical point sizes fall as minus three halves

Zapperi, Lauritsen and Stanley studied the simplest case: each active site makes two new active sites with probability $p$ and otherwise stops. The branching ratio is $2p$ and the critical point is $p=1/2$ . Cascade sizes are then odd, and the exact distribution follows from the Catalan numbers: $P(1)=0.5$ , $P(3)=0.125$ , $P(5)=0.0625$ , $P(7)=0.0390625$ .

For large sizes their Eq. (8) approximates this by $s^{-3/2}$ times an exponential cutoff. The cutoff size $s_c$ diverges at the critical point: 199 at $p=0.45$ and 4999 at $p=0.49$ . We checked the approximation against the exact Catalan formula. At $p=1/2$ the ratio of the two is 0.99963 at $s=2001$ and 0.99265 at $s=101$ ; at $p=0.45$ and $s=401$ it is 0.99814.

They also find cascade durations falling as $t^{-2}$ . The exponent 3/2 is not special to binary branching: we simulated 20000 trees whose offspring follow a Poisson law with mean 1 and fitted them with the maximum-likelihood method of the next chapter, obtaining $1.4981\pm0.0045$ .

The converse does not hold. The next chapter shows systems with no critical cascades at all that give an exponent very close to 3/2.

P(2k+1)=C_k\,p^{k}(1-p)^{k+1},\quad C_k=\frac{1}{k+1}\binom{2k}{k};\qquad P(s)\simeq\sqrt{\frac{2(1-p)}{\pi p}}\;s^{-3/2}e^{-s/s_c},\quad s_c=-\frac{2}{\ln[4p(1-p)]}

Left: the exact cascade-size distribution of the binary branching process. Right: Eq. (8) of Zapperi and co-workers, valid for

1\ll s

, with cutoff size

s_c

Source: S. Zapperi

From spikes to avalanches: the bin width decides the cut

An experiment does not see the family tree of a cascade, only a stream of time-stamped spikes. The usual recipe pools the events from all electrodes and counts them in time bins of equal width. A run of consecutive non-empty bins, with an empty bin on each side, is one avalanche. Its size is the number of events in those bins and its duration is the number of bins.

The count sequence in the figure is cut into three avalanches of sizes 6, 1 and 6 and durations 3, 1 and 2. Double the bin width and neighbouring bins merge, so the last two, separated by a single empty bin, become one. The size distribution of avalanches therefore depends on the bin width as well as on the network.

The Lab’s default electrode array covers 608 of the 1000 units and records 2979 avalanches in 300 000 steps, with mean size 114.5 and largest 7461. Averaged over 10 random seeds, bin widths of 1, 4 and 8 steps give fitted exponents of 1.461, 1.393 and 1.387, a small change. An uncoupled mimic whose input rate jumps among a few levels gives 6.16, 3.96 and 1.87 for the same three widths. Priesemann and Shriki point to this sensitivity to the bin width as one way to tell the two apart.

Cutting avalanches from binned counts. Bars are event counts per bin; each amber bracket marks one avalanche, labelled with its size. Top: twelve bins of one step; runs of non-empty bins between empty ones give avalanches of sizes 6, 1 and 6. Bottom: the same events in bins twice as wide; the last two avalanches merge into one of size 7. The counts are an invented example, not a recording.

Source: V. Priesemann and O. Shriki · S. Zapperi

Recording only part of the network makes the branching ratio look smaller

The most direct estimate regresses the count in the next bin on the count in this bin and reads the slope as the branching ratio; this is the conventional estimate. Wilting and Priesemann show that when only part of the units is recorded, the regression slope $r_k$ at every lag $k$ is multiplied by the same factor $b$ . The one-step slope is biased low, but the rate at which the slopes decay with $k$ is unchanged, so fitting the whole decay recovers $m$ . This is multistep regression. In their example a process with $m=0.9$ recorded at 10% gives a conventional estimate of only 0.312.

If each spike is kept independently with probability 10%, the one-step slope after thinning follows exactly from the mean, variance and one-step covariance of the full record. For one Lab record the prediction is 0.601 and the measurement $0.6015\pm0.0104$ ; across 12 records every deviation is within 2.4 standard errors.

At the Lab defaults the array’s one-step slope is 0.965 and multistep regression gives $0.9848\pm0.0007$ ; multistep regression on all 1000 units gives 0.9855. The nominal branching ratio is 1.05 but the effective one is only 0.985, again because of collisions and because a unit that has just fired cannot fire at once. Kicking one unit of a quiet network triggers 155 spikes on average, which corresponds to $1-1/155=0.9936$ .

r_k=\frac{\operatorname{Cov}(a_{t+k},a_t)}{\operatorname{Var}a_t}=b\,m^{k},\qquad \hat m_{\mathrm{C}}=r_1=b\,m,\qquad \tau=-\frac{\Delta t}{\ln m}

a_t

is the recorded count in bin

t

. Partial recording multiplies every

r_k

by the same

b

; the conventional estimate uses only

r_1

, while multistep regression fits the whole curve.

\tau

is the correlation time of the activity.

Source: J. Wilting and V. Priesemann

Worked example: read the branching ratio from 20 electrodes

Problem. From the same record of the Lab’s default network, use only a 5 by 4 block of 20 electrodes in one corner. The regression slopes at lags one, ten and twenty are $r_1=0.4908$ , $r_{10}=0.2738$ and $r_{20}=0.2175$ . What is the conventional estimate? Use the slopes at ten and twenty steps to estimate $m$ , $b$ and the correlation time, and compare with 0.9855 for the whole network.

Model. Use $r_k=b\,m^k$ . Dividing the two slopes removes $b$ : $m=(r_{20}/r_{10})^{1/10}$ .

Numbers. The conventional estimate is just $r_1=0.491$ ; taken literally, a cascade would hold only $1/(1-0.491)\approx2$ spikes on average. The two-lag estimate is $m=(0.2175/0.2738)^{1/10}=0.9772$ , with $b=0.2738/0.9772^{10}=0.345$ and a correlation time $\tau=-1/\ln0.9772=43$ steps. The whole network’s 0.9855 corresponds to 68 steps.

Explanation. The two-lag estimate lifts 0.491 back to 0.977: most of the bias comes from recording only part of the network and can be removed. The rest comes from where the electrodes sit. Scattering the 20 electrodes at random over the sheet in the same record, five placements give two-lag estimates between 0.981 and 0.988 and multistep estimates between 0.9845 and 0.9864, in line with the whole network. A likely reason is that electrodes packed into one corner see the same local cascade together, and the cascade soon leaves the corner. That adds a fast early decay to the slopes which is not the network’s branching ratio.

Check. The randomly scattered electrodes have one-step slopes of only 0.30 to 0.34, even lower than the corner’s 0.491, yet they give the better branching ratio. The branching ratio sits in how fast the slopes decay with lag; the size of the one-step slope mostly reflects how much of the network the array records.

Source: J. Wilting and V. Priesemann

Try it in the Lab

Critical and subcritical cascades

Predict first: lower the coupling from $\sigma=1.05$ to 0.8. How will the tail of the size distribution and the branching estimate change?
Start with Cascades near the critical point and read the fitted exponent and the multistep estimate. Then choose Cascades that die out and compare the size chart and the two readouts.
Explain: the branching ratio drops to about 0.77, cascades die quickly, the largest avalanche holds little more than a hundred events and the power-law fit is rejected. Kicking one unit now gives only a few spikes.

Bin width and electrode placement

Predict first: in the near-critical network, drag the analysis bin from 1 step to 8. How much will the fitted exponent change? And for Drifting input, no coupling?
Drag the bin in each of the two presets and note the exponent. Then choose A few electrodes see little and drag the array or its corner while watching the branching estimate.
Explain: the coupled network’s exponent moves only from about 1.46 to about 1.39, while the drifting mimic changes far more. A small, compact array also biases the multistep estimate low; enlarging the array brings it back near 0.985.

Exercises with explained answers

At $\sigma=0.5$ , what is the probability that a spike has no offspring at all? What if the offspring count is approximated by a Poisson law with mean 0.5?
Show the answer
The binomial law gives $(1-0.05)^{10}=0.599$ and the Poisson approximation $e^{-0.5}=0.607$ . The simulation gives 0.6006, matching the binomial law; with only 10 neighbours the Poisson approximation is slightly high.
How close to 1/2 must $p$ be for the cutoff size of the binary branching process to exceed 7461, the largest avalanche in the Lab’s default record?
Show the answer
From $-2/\ln[4p(1-p)]=7461$ we get $4p(1-p)=e^{-2/7461}$ , so $p=0.4918$ , less than 0.01 from 1/2. In a finite record, slightly subcritical and exactly critical are hard to tell apart.
With a bin width of 1, the count sequence 0, 2, 3, 1, 0, 0, 1, 0, 4, 2, 0, 0 gives three avalanches of sizes 6, 1 and 6. What sizes and durations do you get after doubling the bin width, merging each pair of neighbouring bins?
Show the answer
The merged sequence is 2, 4, 0, 1, 6, 0, which holds only two avalanches: size 6 lasting 2 bins and size 7 lasting 2 bins. The avalanche of size 1 has merged into the next one.
Partial recording multiplies every slope by $b=0.6$ and the true branching ratio is $m=0.985$ . What are the slopes at one and thirty steps?
Show the answer
$r_1=0.6\times0.985=0.591$ and $r_{30}=0.6\times0.985^{30}=0.381$ . Looking only at $r_1$ would wrongly call the network clearly subcritical.

Sources and what we recomputed

S. Zapperi, K. B. Lauritsen and H. E. Stanley, Self-organized branching processes: mean-field theory for avalanches, arXiv:cond-mat/9509161v1; Physical Review Letters 75, 4071 (1995)
Version 1, read in full. The binary branching process with critical probability one half. Eqs. (6a) and (7): the generating function of avalanche sizes. Eq. (8): the size distribution with exponent 3/2 and cutoff size −2/ln[4p(1−p)]. Eq. (11): the critical distribution. Lifetimes with exponent 2. We recomputed the exact small sizes from the Catalan numbers and checked the asymptotic prefactor.
J. Wilting and V. Priesemann, Inferring collective dynamical states from widely unobserved systems, arXiv:1608.07035v3; Nature Communications 9, 2325 (2018)
Version 3, main text pp. 1 to 3; the supplementary methods were not read. The branching process with drive, the conventional estimator as the slope of the next count against the current count, the result that subsampling multiplies every lagged slope by the same factor, and the multistep regression estimator that fits that exponential across lags. Their examples: a process with m = 0.9 sampled at 10 % gives a conventional estimate of 0.312.
V. Priesemann and O. Shriki, Can a time varying external drive give rise to apparent criticality in neural systems?, PLoS Computational Biology 14, e1006081 (2018)
Introduction, Results and Discussion read; Methods not read. How avalanches are defined from binned events. A Poisson process with a constant rate gives no power law; one whose rate switches among a few levels gives sizes close to a power law with exponent near 3/2 over a limited range. The measures that tell the two apart, and their conclusion that passive analysis cannot rule out a drive that mimics critical activity, so a small controlled perturbation is needed.