1. Neuroscience
Download icon

Maintained avalanche dynamics during task-induced changes of neuronal activity in nonhuman primates

  1. Shan Yu
  2. Tiago L Ribeiro
  3. Christian Meisel
  4. Samantha Chou
  5. Andrew Mitz
  6. Richard Saunders
  7. Dietmar Plenz  Is a corresponding author
  1. National Institute of Mental Health, United States
Research Article
  • Cited 1
  • Views 578
  • Annotations
Cite this article as: eLife 2017;6:e27119 doi: 10.7554/eLife.27119

Abstract

Sensory events, cognitive processing and motor actions correlate with transient changes in neuronal activity. In cortex, these transients form widespread spatiotemporal patterns with largely unknown statistical regularities. Here, we show that activity associated with behavioral events carry the signature of scale-invariant spatiotemporal clusters, neuronal avalanches. Using high-density microelectrode arrays in nonhuman primates, we recorded extracellular unit activity and the local field potential (LFP) in premotor and prefrontal cortex during motor and cognitive tasks. Unit activity and negative LFP deflections (nLFP) consistently changed in rate at single electrodes during tasks. Accordingly, nLFP clusters on the array deviated from scale-invariance compared to ongoing activity. Scale-invariance was recovered using ‘adaptive binning’, that is identifying clusters at temporal resolution given by task-induced changes in nLFP rate. Measures of LFP synchronization confirmed and computer simulations detailed our findings. We suggest optimization principles identified for avalanches during ongoing activity to apply to cortical information processing during behavior.

https://doi.org/10.7554/eLife.27119.001

Introduction

Neuronal activity in the brain has been traditionally separated into ongoing activity, which lacks a particular sensory stimulus or movement, and evoked activity, which is the response to a well-defined stimulus. Yet, studies have consistently shown that ongoing activity predicts stimulus response (Arieli et al., 1996; Tsodyks et al., 1999; Luczak et al., 2009) and behavioral outcome (Supèr et al., 2003; Womelsdorf et al., 2006) suggesting that both forms of activity are closely related. During the last decade, ongoing activity has been found to organize as neuronal avalanches (Beggs and Plenz, 2003) — scale-invariant activity cascades whose size, duration, waveform and inter-cascade intervals are governed by power laws (for review see Chialvo (2010) and Plenz (2012)). From microscale to macroscale organization, neuronal avalanches describe spontaneous firing of local pyramidal neuron groups in awake rodents in vivo (Bellay et al., 2015), the ongoing local field potential (LFP) and subthreshold population activity in rodents and nonhuman primates (Gireesh and Plenz, 2008; Petermann et al., 2009; Scott et al., 2014), as well as resting activity in humans observed using functional magnetic resonance imaging (fMRI) (Fraiman and Chialvo, 2012; Tagliazucchi et al., 2012; Haimovici et al., 2013), magnetoencephalography (MEG) (Palva et al., 2013; Shriki et al., 2013) and electrocorticography (ECoG) (Priesemann et al., 2014). Neuronal avalanches signify a cortical network at a critical or near-critical state (Chialvo, 2010; Plenz, 2012), which experiments (Shew et al., 2009; Gautam et al., 2015) and simulations (Beggs and Plenz, 2003; Bertschinger and Natschläger, 2004; Haldeman and Beggs, 2005; Kinouchi and Copelli, 2006; Rämö et al., 2007; Nykter et al., 2008) have shown to hold numerous advantages in information processing, such as maximum mnemonic repertoire size (Haldeman and Beggs, 2005; Shew et al., 2011; Fagerholm et al., 2016), information diversity (Nykter et al., 2008) and dynamic range (Kinouchi and Copelli, 2006; Shew et al., 2009; Gautam et al., 2015), optimized computational capabilities (Bertschinger and Natschläger, 2004), information transmission (Beggs and Plenz, 2003; Rämö et al., 2007; Fagerholm et al., 2016), sensory discrimination (Tomen et al., 2014) and learning (de Arcangelis and Herrmann, 2010). Yet, despite these potential advantages, evidence for neuronal avalanches outside the realm of ongoing activity, specifically, during behaviorally relevant motor and cognitive tasks, has been controversial.

Recent analysis of neuronal avalanches during sustained dynamic movie stimulation in the visual cortex of turtle (Shew et al., 2015) and visuo-motor tasks in human fMRI (Fagerholm et al., 2015) and EEG (Arviv et al., 2015) suggest transient dynamics that differ from scale-invariant avalanches. Similarly, neural models implied transient deviations from avalanche dynamics during strong external inputs (Millman et al., 2010; Taylor et al., 2013; Hartley et al., 2014; Stepp et al., 2015; Williams-García et al., 2014). Here, we hypothesize that transient changes in activity during evoked responses could simply reflect a change in the rate of avalanches rather than a transition to a different dynamical regime, for example a change in synchronization. Such a rate change likely requires modifications to standard avalanche analysis, as the standard was originally introduced for stationary rates of activity (Beggs and Plenz, 2003). In fact, by taking changes in activity rate into account, here we show in the behaving nonhuman primate that cortical dynamics during sensory and cognitive processing exhibit clear power law organization in line with the maintenance of neuronal avalanches. In line with these findings, our phase-synchronization analysis demonstrates that synchrony levels are maintained for ongoing and task-related activity. We expand on our experimental results using computational models and explore the sensitivity of this advanced approach to distinguish avalanche dynamics from states with changed synchrony, for example supercritical dynamics. Our results extend previous work on the dynamics of resting state activity and suggest that the optimization principles identified for critical dynamics during resting state activity extend to and sculpt cortical information processing during behavior.

Results

In order to demonstrate avalanche maintenance during evoked activity that is independent of specific brain region or behavioral paradigm, we chose two types of tasks and two cortical areas in two adult macaque monkeys. High-density microelectrode arrays (MEA; 10 × 10; corner electrodes missing; inter electrode distance Δd = 400 µm) were chronically implanted in left prefrontal cortex (PFC) of monkey A and left premotor cortex (PM) of monkey B. We recorded extracellular unit activity (0.3–3 kHz) to confirm involvement of recording sites in behavioral tasks as well as the LFP (1–100 Hz), which is commonly used for avalanche analysis (Beggs and Plenz, 2003; Petermann et al., 2009; Scott et al., 2014).

Avalanche maintenance during a cue-triggered cognitive task

In our first behavioral condition, we asked whether scale-invariant dynamics can be found when strong transient activity changes due to cognitive processing are present. To this end, we analyzed avalanche dynamics in the dorsal-lateral prefrontal cortex (dlPFC; left hemisphere) for a visual-motor mapping task (monkey A; Figures 13). Monkey A applied a learned rule (>95% success rate) to two individually presented visual cues, which instructed the retrieval of food from a left feeder (left trials) or right feeder (right trials).

Neuronal avalanches in prefrontal cortex during a cognitive task.

(A) Firing rate changes for putative single units demonstrate PFC region recorded by the array is involved in task performance. Unit raster (top) separated into right (green) and left (black) trials with corresponding average firing rates (bottom). Cue presentation elicits distinct rate changes for units 1 and 2 during left and right trials, but not for unit 3. Colored area (orange/purple) indicates periods that differ significantly (high/low) from baseline (see Materials and methods). Vertical solid lines: Cue on- and off-set. (B) The population of PFC units does not show a task-related change in rate (top) or Fano Factor (bottom). Top: mean-matched/original rate in solid/broken lines; 38/48% survived matching for left/right trials. Bottom: mean-matched spikes. Shaded areas: 95% confidence interval from linear regression per time window. (C) Distinct changes in the LFP for right (green) and left (black) trials during the task. Grey/light green: trial averages for each electrode. Black/dark green: average over electrodes. (D) The negative LFP (nLFP; −2SD threshold; dots) allows for distinguishing right (green) from left (black) trials. Top: Example nLFP raster (single electrode) separated into right (green) and left (black) trials. Bottom: Corresponding time course in average nLFP rate. Colored areas (orange/purple) indicate significant change (high/low) from baseline. (E) Distinct change in nLFP rate separates a baseline (BASE; 400 ms) period from an early (EARLY; 400 ms) and late (LATE; 400 ms) epoch after cue onset (vertical lines). (F) nLFPs on the array when analyzed for full left and right trial periods show avalanche organization. Power law in size probability densities for nLFP clusters (solid; fixed Δt). Arrow: Cut-off at number of electrodes on the array. Broken lines: Corresponding trial-shuffle controls. Shaded areas: Corresponding 95% confidence interval based on bootstrapping. Dotted line: Visual guide for exponent of −1.5.

https://doi.org/10.7554/eLife.27119.002
Adaptive binning tracks avalanche organization during behaviorally induced transient activity changes in PFC.

(A) Fixed binning defines clusters (grey area) by successive time bins of constant duration Δt with at least one nLFP (diamonds). The temporal resolution Δt does not change for trials or epochs. Size s is defined as the number of nLFPs per cluster. (B) Large avalanches (red) dominate for left trials during LATE. Avalanche raster with size coded by color and dot size at fixed Δt for left (top) and right (bottom) trials for 2 consecutive recording days. (C) Significant increases (orange) and decreases (purple) in average time course for avalanche rate (top) and size (bottom) for left (black) and right (green) trials. (D) Adaptive binning links Δt for each trial i to the average nLFP rate during epochs resulting in three different temporal resolutions for cluster definition: ΔtiBASE,EARLY,LATE (cp. A). (E) Same as in B, but for adaptive binning. Note sparseness of very large avalanches (cp. B). (F) Adaptive binning increases avalanche rate significantly during LATE while reducing avalanche size (cp. C).

https://doi.org/10.7554/eLife.27119.003
Adaptive binning demonstrates maintained power law organization during transient activity epochs in prefrontal cortex.

(A) Avalanche size probability densities for LATE (orange) significantly deviate from the power law observed during BASE (black/green) and EARLY (purple) for left (left panel), but not right trials (right panel). Area: 95% confidence interval. (B) Size distributions normalized by a power law with exponent −1.5 reveal good agreement with theory (horizontal curve with a sharp cut-off), except during LATE for left trials, where a large deviation can be observed. Inset: Area under the curve of the normalized distributions for avalanches of size > 20 quantify differences between epochs. Error bars: SD from 1000 bootstraps. (C) Same as in A, but for adaptive binning, which collapses size distributions for BASE, EARLY and LATE. (D) Adaptive binning leads to a decrease in the area difference between LATE and the other epochs. (E) Distributions of the difference between the mean avalanche size of LATE and BASE for each trial (positive values indicate avalanches in LATE were larger than those in BASE for a given trial) for fixed (grey) and adaptive binning (black/green). Significant differences between distributions are marked by an asterisk (Kolmogorov-Smirnov test; p0.05). For left trials, LATE is significantly biased toward smaller avalanches for adaptive binning compared to fixed binning.

https://doi.org/10.7554/eLife.27119.004

It is known that the dlPFC is involved in such visual-motor mapping task (Asaad et al., 1998; Puig and Miller, 2012). Consistently, in all of 29 putative single-units recorded in the dlPFC (1 day, 1 hr recording session), unit activity was variable with transient increases or decreases in firing rate around cue on- and off-set that depended on left or right trials (Figure 1A; three examples shown). No significant change in Fano Factor (Churchland et al., 2010) was found during task execution (Figure 1B; see Materials and methods). On the other hand, the average LFP changed consistently for each electrode exhibiting two waveforms that differed between left and right trials (Figure 1C). We extracted peak times of the negative LFP (nLFP) at a threshold thr = –2 SD for further analysis (see Materials and methods), and observed a rapid ~8-fold increase in nLFP rate ~0.5 s after cue presentation, when the monkey was preparing to reach for the left feeder (Figure 1D, single electrode; Figure 1E, all electrodes). A much weaker, more delayed increase was observed for right trials. Based on the delayed increase in nLFP rate for left trials, we divided trials into a baseline period before cue-onset (BASE; –0.4–0 s), and an early (EARLY; 0–0.4 s) and late epoch (LATE; 0.4–0.8 s) after cue-onset.

To study neuronal avalanches, two 1 hr recording sessions from 2 consecutive days resulting in 322 trials (day 1: 139, day 2: 183) were combined. In short, for each trial successive nLFPs on the array were concatenated into contiguous clusters at temporal resolution Δt (Beggs and Plenz, 2003). A cluster started with the transition from an empty time bin of width Δt to a time bin with at least one nLFP on the array. The concatenation continued until a time bin with no nLFP was encountered (see also Figure 2A). The temporal resolution Δt was defined by the inverse of the average nLFP rate on the array and is not a free parameter (Beggs and Plenz, 2003). In line with standard avalanche analysis, we kept Δt fixed (‘fixed binning’) and determined its value from the full recording to be Δt = 4.75 ms for monkey A. We found that the probability density function of cluster size s, that is the number of nLFPs in a cluster, followed a power law, for both left and right trials, with slope α close to –1.5 and a cut-off for s > 96, the maximal number of MEA electrodes (Klaus et al., 2011; Yu et al., 2014) (Figure 1F, solid black and green lines for left and right trials, respectively; exponential vs. power law: p<105; see Materials and methods). Importantly, these hallmarks of avalanche dynamics were destroyed by trial shuffling, which removes spatial trial-to-trial correlations while maintaining average rate changes (Figure 1F, broken lines; exponential vs. power law: 1p<105). This control demonstrates that the power law in cluster sizes, indicative of avalanche dynamics, does not simply arise from transient activity changes during behaviorally relevant periods.

By combining clusters from different behavioral epochs into one single distribution, however, transient deviations from scale-invariance could cancel each other, giving the appearance of maintained avalanche dynamics during behavior. Indeed, when separating clusters according to trial epochs BASE, EARLY and LATE, systematic differences in the corresponding avalanche raster were revealed. Specifically, the increase in nLFP rate during LATE for left trials yielded significantly larger avalanches during the period prior to reach (Figure 2B, red dots; Figure 2C, orange and purple segments indicate significantly high and low values, respectively; p0.05; see Materials and methods). Accordingly, a preponderance of large avalanches was visible in the corresponding size distribution during LATE for left trials (Figure 3A and B, left), whereas size distributions remained close to a power law during BASE and EARLY and for right trials (Figure 3A and B; p<105, all cases). In the case of right trials, where the transient changes in activity rate are much less pronounced, size distributions are more similar among all epochs, with a small tendency for larger avalanches during the LATE epoch (Figure 3A and B, right; p<105, all cases).

These deviations from scale-invariance compared to BASE, on the other hand, were abated when the temporal resolution Δt to define clusters was obtained from the average nLFP rate for each trial i and epoch separately (ΔtiBASE,EARLY,LATE; Figure 2D–F). We define this approach as ‘adaptive binning’. It decreases Δt during periods of high nLFP rate, while increasing Δt for low activity periods (Figure 2D). Indeed, adaptive binning obtained scale-invariant power laws in avalanche sizes that were similar for each epoch (Figure 3C,D; p=0.06 (0.08) and p=0.11 (0.1) for EARLY (LATE) compared to BASE, for left and right trials respectively; see Materials and methods), with a change in exponent consistent with previous results (Beggs and Plenz, 2003; Petermann et al., 2009; Priesemann et al., 2014) and a cut-off for s > 96 (Klaus et al., 2011; Yu et al., 2014). Comparing BASE and LATE avalanches for each trial, the bias toward larger ones in the latter epoch seen when employing fixed binning is significantly decreased after adaptive binning for left trials (Figure 3E, left), while during right trials, in which scale-free avalanches were obtained even for fixed binning, the bias does not change after adaptive binning. As can be seen in the corresponding raster plots and trial averages (Figure 2E,F), avalanche size decreased in LATE concomitant with an increase in avalanche rate (more prominently for left trials). We conclude that the apparent over-representations of large avalanches during behavioral epochs were due to a mismatched temporal resolution Δt by not accounting for systematic changes in nLFP rate during different epochs.

We further explored if nLFP rate within an epoch correlates with behavioral outcome. Indeed, reaction time, defined as the time it took the monkey to reach the reward after the visual cue was switched off, did not correlate with nLFP rate (left trials: R = 0.036, p=0.73; right trials: R = −0.039, p=0.71; p-values indicate the likelihood of obtaining a correlation at least as high by chance), avalanche rate (left trials: R = 0.041, p=0.68 for fixed Δt and R = −0.022, p=0.79 for adaptive Δt; right trials: R = −0.033, p=0.74 for fixed Δt and R = 0.012, p=0.87 for adaptive Δt) or avalanche size (left trials: R = 0.085, p=0.52 for fixed Δt and R = 0.062, p=0.66 for adaptive Δt; right trials: R = 0.058, p=0.63 for fixed Δt and R = −0.025, p=0.72 for adaptive Δt), in line with our finding of a scale-free regime regardless of rates for each trial. Given the monkey’s high performance (>95% successful trials), rates could not be analyzed for correct/incorrect trial outcome.

Avalanche maintenance during self-initiated movements

In monkey B, we studied activity changes during self-initiated movements. The monkey, without further cue or prompting, initiated the touching of a pad with her right hand, after which a reward was provided. Similar to monkey A, two 1 hr sessions were recorded in 2 consecutive days. Since results were highly consistent between different days, the data were pooled together, resulting in 663 trials analyzed (day 1: 319; day 2: 344). As the recording site was within the arm representation area within the PM, recorded activity clearly reflected task execution (Figure 4). About 200 ms before touching the pad, 18/43 (~42%) putative single-units increased in firing rate, whereas 5/43 units (~12%) decreased their firing rate around the time of touching (Figure 4A,B). Fano Factor of unit activity decreased in line with previous reports (Figure 4C) (Churchland et al., 2010). Similarly, the LFP transiently turned negative ~200 ms before touching for all electrodes (Figure 4D), together with a ~4-fold increase in the nLFP rate and a slight decrease thereafter (Figure 4E), which readily separated activity into a baseline epoch (BASE; –0.8 to –0.4 s relative to the touch), followed by a pre-touch (PRE; –0.4 to 0 s) and post-touch epoch (POST; 0 to 0.4 s) (Figure 4F).

Neuronal avalanches in premotor cortex during self-initiated motor task.

(A) Change in firing rate for two putative single units demonstrate PM region recorded by the array is involved in task performance. Unit one increased whereas unit two decreased firing around the time of self-initiated touch (grey arrow). Unit raster (top) trial-aligned to touch (grey line, arrow) and corresponding average firing rate time course (bottom). Orange/purple area indicate significant increase/decrease from baseline. (B) Most units show transient excitation before touching. Average firing rate change for all putative single units normalized to baseline. (C) Average (top, broken) and mean-matched (top, solid) firing rate across all units and trials. Fano Factor (bottom) decreases significantly (purple) before touching (70% of data survived mean-matching). Shaded area: 95% confidence interval. (D) Negative deflection in the LFP, averaged over trials, prior to touch. Grey: trial average for individual electrodes. Black: average over electrodes. (E) Single electrode nLFP raster and corresponding average nLFP rate. Orange/purple areas indicate significantly high/low periods from baseline. (F) Change in nLFP rate defines baseline (BASE), pre-touch (PRE) and post-touch (POST) epochs (vertical lines). (G) nLFPs on the array when analyzed for full trial periods show avalanche organization. Power law in size probability densities for nLFP clusters (solid; fixed Δt). Arrow: Cut-off at number of electrodes on the array. Broken lines: Corresponding trial-shuffle control. Shaded areas: Corresponding 95% confidence interval based on bootstrapping. Dotted line: Visual guide for exponent of −1.5.

https://doi.org/10.7554/eLife.27119.005

When all data across epochs was analyzed using fixed binning, that is Δt based on the full recording (Δt = 10.3 ms), avalanche sizes followed a power law with –1.5 exponent and cut-off at ~96 electrodes (Figure 4G; exponential vs. power law: p<105), which was destroyed by trial shuffling (exponential vs. power law: 1p<105), again demonstrating that scale-invariant LFP clusters were based on the trial-by-trial correlations between cortical sites and did not simply arise from transient changes in nLFP rate.

Evaluation of avalanche raster plots based on epochs revealed systematic differences in the corresponding avalanche raster (Figure 5A,B). Specifically, the increase in nLFP rate during PRE yielded significantly more avalanches and of larger size during that period, whereas the converse was true during POST when nLFP rate decreased (Figure 5A; colored dots; Figure 5B; orange and purple segments respectively; p0.05; see Materials and methods). The corresponding size distributions exhibited deviations from scale-invariance obtained for BASE (Figure 6A and B, left) with an excess of large avalanches for PRE and a deficit of large avalanches for POST, despite all distributions remaining close to power laws (exponential vs. power law: p<105, all three epochs).

Adaptive binning tracks avalanche organization during self-initiated motor task in premotor cortex.

(A) Avalanche raster with size coded by color and dot size for fixed Δt. Large avalanches (red) emerge during PRE and fewer avalanches are found during POST. (B) Significant increases (orange) and decreases (purple) in average time course for avalanche rate (top) and size (bottom) with respect to BASE. (C) Same as in A, but for adaptive binning. Note absence of very large avalanches during PRE. (D) Adaptive binning reduces avalanche size while increasing avalanche rate (cp. B) during PRE.

https://doi.org/10.7554/eLife.27119.006
Adaptive binning demonstrates maintained power law organization during transient activity epochs in premotor cortex.

(A) Avalanche size distributions for PRE and POST deviated from the power law observed during BASE when fixed binning was employed (left), but collapse together after adaptive binning (right). Area: 95% confidence interval. (B) Size distributions normalized by a power law (exponent −1.5 for fixed and −1.3 for adaptive binning). Inset: Area under the curve for the normalized distributions, considering avalanches of size > 20, emphasize the bias toward larger (smaller) avalanches during PRE (POST), in line with higher (lower) activity rates observed during that epoch compared to BASE for fixed binning (left). The bias is significantly decreased after adaptive binning (right). Error bars: SD from 1000 bootstraps. (C) Distributions of the difference between the mean avalanche size of PRE and BASE during each trial for fixed (grey) and adaptive binning (black). Large differences occur much more often for fixed binning (left: linear-linear plot; right: log-log plot). Significant differences between distributions are marked by an asterisk (Kolmogorov-Smirnov test; p0.05).

https://doi.org/10.7554/eLife.27119.007

As in the case for the cognitive task, these systematic deviations in size distributions were abated when calculating Δt from the average nLFP rate for each trial i and epoch (ΔtiBASE,PRE,POST; Figure 6A and B, right). Through adaptive binning we once again obtained scale-invariant power laws in avalanche size that were similar for each epoch (p=0.177, PRE vs. BASE; p=0.05, POST vs. BASE), with the bias towards larger avalanches during PRE being significantly reduced (Figure 6C). As can be seen in the corresponding raster plots and trial averages (Figure 5C,D), avalanche size decreased in PRE concomitant with an increase in avalanche rate confirming our findings in monkey A. During POST, under representation of large avalanches due to mismatched temporal resolution Δt can also be partially compensated for by adaptive binning.

Impact of ‘adaptive binning’ in capturing rate change vs. change in dynamics

We next illustrate the expected deviation from scale-invariance when rates change and its remediation using ‘adaptive binning’ in a simple schematic. In the hypothetical example shown in Figure 7A, we consider a baseline period with nLFP rate of 3 Hz per electrode, followed by a period with the rate increased to 9 Hz (Figure 7A, left). Fixed binning, which sets Δt according to the average rate of 6 Hz, underestimates the temporal resolution for baseline activity, whereas it overestimates Δt for the high activity period. The resulting power law size distributions will be steeper for the baseline period and shallower for the high activity period (Beggs and Plenz, 2003). The combined distribution will show an upward bend at the cross-over of the individual distributions (Figure 7A; middle, solid line) which seems to suggest a deviation from avalanche dynamics with an overabundance of large clusters similar to what would be expected for supercritical dynamics. Adaptive binning, that is calculating Δt for each period separately approximates the exponents for both power laws, resulting in a final distribution that matches a power law with an intermediate slope (Figure 7A, right) (Beggs and Plenz, 2003; Petermann et al., 2009; Priesemann et al., 2014). We note that simply using shorter time bins for all data would retain different slopes for different regimes and thus is not an alternative to adaptive binning.

Adaptive binning recovers the power law in the face of consistent rate changes when dynamics remain critical.

(A) Schematic impact of avalanche analysis from non-stationary event rates. Left: A low event rate period (purple) followed by a high event rate period (orange) results in an intermediate Δt (fixed binning) based on the mean event rate (broken line). Middle: The superposition of two power laws with different slopes from their corresponding rate regime (broken lines) results in an avalanche size distribution (grey) that deviates from a power law, with a characteristic up-ward bend at the cross-over point. Right: Adaptive binning steepens/reduces the slope for the high/low rate period respectively resulting in a distribution collapse at an intermediate slope. (B) Left: Average nLFP rate for monkey B with trials separated into high (black) and low (green) nLFP rate during PRE (shaded area). Middle: Corresponding distributions obtained with fixed binning increasingly deviate with rate (arrow) from baseline (broken, purple). Right: Adaptive binning collapses all distributions. (C) Simulations using a transient external Poisson drive match experimental findings (cp. B). (D) In simulations with transiently switching from critical to supercritical dynamics (left), adaptive binning fails to compensate for overabundance of large avalanche sizes from supercritical dynamics (middle/right). (E) Simulations using transient Poisson noise. Distributions obtained with fixed binning do not follow power laws, even in the low-noise regime (middle). Distributions obtained with adaptive binning do not have a clear cut-off (right).

https://doi.org/10.7554/eLife.27119.008

Although we illustrated this principle using two sequential periods with different rates, the same holds for repeat trials with different rates. In Figure 7B, we separated trials of monkey B into low- and high-rate responses based on PRE epochs (Figure 7B, left). Size distributions obtained from high-rate trials at fixed binning exhibit an overabundance of large avalanches (Figure 7B, middle, arrow). With adaptive binning, deviations from the baseline distribution decreased and all distributions approached power laws with a –1.2 exponent (Figure 7B, right). This demonstrates that trials in monkey B exhibit similar dynamics as those observed for monkey A or, more generally, that deviations from scale-invariance observed in the two experimental regimes can be explained by nLFP rate changes alone.

Neuronal simulations support ‘adaptive’ binning to recover scale-invariance for critical dynamics

Using simulations, we explored the conditions under which adaptive binning recovers scale-invariant size distributions. Simulations were carried out on a 10 × 10 cellular automaton network, the approximate size of our MEAs, with cascades unfolding according to a critical branching process (Zapperi et al., 1995) (see Materials and methods). Quiescent times between cascades were randomly drawn from the experimentally obtained quiescent time distribution of monkey B. Simulated baseline activity, by randomly initiating avalanches in the network, was superimposed after 0.4 s by a second process, which produced a transient increase in rate mimicking movement initiation. We tested three different processes for two levels of rate increase and obtained corresponding size distributions for fixed and adaptive binning. Increasing rate by adding Poisson inputs that could trigger avalanches (Figure 7C, ‘Input driven’) closely matched our experimental results. For fixed binning, an overabundance of large avalanches close to the cut-off was found in the size distribution (Figure 7C, middle, arrow), which originated from the concatenation of spontaneous with input-triggered avalanches. We note that this overabundance manifests while the model interactions remained tuned to the critical point in the limit of zero input. Accordingly, the power law in the size distribution was recovered by means of adaptive binning (Figure 7C, right). Next, a high rate was achieved by increasing the likelihood of nodes exciting each other, that is changing to a supercritical branching process with increased synchronization (Figure 7D, ‘Supercritical’). The corresponding size distribution exhibits an overabundance of large avalanches that (Figure 7D, middle) reflects explosive growth as activity propagates in the network. Adaptive binning, which is tied to the overall increase in event rate, was insufficient to compensate for the increase in synchrony. It failed to collapse the cut-off of the distributions, instead further separating non-global cascades from synchronized, global activity (Figure 7D, right, arrows). Finally, event rate was increased by adding uncorrelated activity, that is Poisson inputs that do not trigger avalanches (Figure 7E, ‘Poisson noise’). The resulting size distribution exhibits a preponderance of intermediate size avalanches that reflects the mean rate of uncorrelated events introduced (Figure 7E, middle). Here, even low-noise levels destroy the power-law regime. While distributions tend to get closer to a power law with adaptive binning, their cut-off at system size is lost in the process (Figure 7E, right).

‘Adaptive’ thresholding as an alternative approach to ‘adaptive’ binning

Adaptive binning reduces Δt for periods of high rates, thereby reducing the number of nLFPs per time bin. Alternatively, one can increase the threshold for nLFP detection during these periods, likewise reducing the number of nLFPs encountered per time bin (see also [Petermann et al., 2009]). The reverse argument holds when periods of low rates are encountered. We tracked the mean and standard deviation of the LFP in successive windows of 100 ms width and used a moving threshold thr = mean – 2SD to identify nLFPs within each window (Figure 8A,C). Thus, activity modulations induced by the task were mainly reflected in nLFP amplitude, but not rate. This approach, which trades temporal resolution for local sensitivity, collapsed distributions from different periods (grey and orange areas in Figure 8B,C) in monkey A into a similar power law, indicative of sustained avalanche dynamics during task performance (Figure 8D,E).

Adaptive thresholding as an alternative means to collapse size distributions during rate changes.

(A) Example LFP on single electrode in monkey A. Solid grey lines: Cue presentation. nLFPs detected at fixed (broken, purple) or adaptive threshold (solid, purple) vary in size and rate accordingly. Triangles: nLFPs obtained with fixed (empty) or adaptive thresholding (filled). (B) nLFP raster (top) and corresponding average nLFP rate (bottom) from A at fixed threshold. Color code: nLFP amplitude in μV. (C) Adaptive thresholding produces a relatively constant event rate (cp. B). (D) Adaptive thresholding (right) successfully collapses avalanche size distributions obtained for fixed thresholding (left; cp. Figure 3A and C, left). (E) Distributions normalized by a power law of −1.5 (fixed thresholding, left) or −1.3 (adaptive thresholding, right) together with area differences are in line with the results from adaptive binning (cp. Figure 3B and D, left). Error bars: SD from 1000 bootstraps.

https://doi.org/10.7554/eLife.27119.009

Phase synchrony analysis supports recovered scale-invariance through ‘adaptive’ binning

Experiments and theory have consistently shown that size distributions of clustered activity are sensitive to the degree of synchrony in a system (Shew et al., 2009; Larremore et al., 2011; Shew et al., 2011; Yang et al., 2012; Bellay et al., 2015; Shew et al., 2015). Our demonstration of scale-invariant size distributions suggests avalanche dynamics are maintained throughout task performance, which should mirror a sustained average synchrony among cortical sites. Previous studies have shown overall phase-synchronization between sites sensitively captures fast dynamical changes in a network (Kitzbichler et al., 2009; Yang et al., 2012), in particular at higher frequencies (Meisel et al., 2016). In line with our results from adaptive binning and thresholding, levels of LFP synchrony R (see Materials and methods) were not different across epochs in monkey A (Left: 0.56 ± 0.015, 0.56 ± 0.01 and 0.56 ± 0.01; Right: 0.56 ± 0.01, 0.55 ± 0.01 and 0.56 ± 0.004; mean ± standard deviation for BASE, EARLY and LATE respectively; p>0.05; two-tailed paired student’s t-test) as well as in monkey B (BASE, 0.450 ± 0.003; mean ± standard deviation; vs. PRE, 0.452 ± 0.003; vs. POST, 0.455 ± 0.003; p>0.05).

Discussion

Our findings suggest that cortical avalanches are maintained in nonhuman primates during information processing, thus adding avalanche dynamics to the many links found between ongoing and evoked activities in the brain (Arieli et al., 1996; Tsodyks et al., 1999; Womelsdorf et al., 2006). We suggest that optimal information processing capabilities of avalanche dynamics (Beggs and Plenz, 2003; Bertschinger and Natschläger, 2004; Haldeman and Beggs, 2005; Kinouchi and Copelli, 2006; Rämö et al., 2007; Nykter et al., 2008; Shew et al., 2009; Tomen et al., 2014; Gautam et al., 2015) might guide sensory, motor and cognitive processing in the brain.

During ongoing activity, avalanche dynamics have been demonstrated in various brain areas including the primary visual (Hahn et al., 2010), primary motor (Petermann et al., 2009), somatosensory (Gireesh and Plenz, 2008), premotor and prefrontal (Yu et al., 2014) cortices, reflecting a common dynamic feature regardless of specific cortical areas. Therefore, we aimed to extract a similar rule that govern avalanches during evoked activity. To achieve this goal while minimizing the number of nonhuman primates used, we chose to record activity from different cortical sites and different behavioral tasks in two monkeys, while the common condition is that the recorded region was strongly activated by the corresponding task. Monkey A performed a visual-motor mapping task, in which the correct actions have to be applied to specific cues (Wise and Murray, 2000). This task is known to involve the dorsal-lateral prefrontal cortex (Asaad et al., 1998; Puig and Miller, 2012), in line with our finding of differential spiking as well as LFP responses for the left and right cues. Conversely, the target area in our self-initiated motor task for monkey B was the premotor cortex, which is involved in movement planning (Weinrich et al., 1984; Brasted and Wise, 2004) and self-initiated movement (Hoffstaedter et al., 2013). Consistent with previous studies, we found that, at both the single neuron and population level, the premotor cortex exhibited strong activation a few hundred milliseconds before the self-initiated touch. Our highly consistent results obtained for these different cortical areas and behavioral tasks suggest that avalanche dynamics might be a general principle governing cortical dynamics underlying various behaviors.

Our trial-shuffling control demonstrated that the power law in avalanche sizes reflects correlations that are maintained during behavior on a trial-by-trial basis and that this statistical feature did not simply arise from non-stationary rates. Recently, power laws, for example Zipf’s law, have been shown to arise when the range of event probabilities becomes enlarged, for example from underlying common correlations due to common inputs typically captured in latent variables (Schwab et al., 2014; Aitchison et al., 2016). Because trial-shuffling maintains such latent variables as captured in the average evoked response, yet, abolishes scale-invariance, our findings suggest that the scale-invariance encountered in our behavioral data reflects intrinsic correlation structure of trial-to-trial variability and not common correlations.

Originally introduced for resting or spontaneous activity (Beggs and Plenz, 2003), a fixed temporal resolution Δt assumes stationary rates, which is typically violated during sensory epochs or periods of movement when event rates change systematically. For avalanche analysis, the temporal resolution Δt is based on the inter-event interval distribution of the population, for example successive nLFP events on the array. This temporal scale Δt is not a free parameter. It further depends on the spatial scale, that is inter-electrode distance Δd, as well as local minimal event threshold thr. Importantly, all three parameters Δt, Δd and thr have been experimentally linked and shown to affect the slope of the power law in a predictable manner (Beggs and Plenz, 2003; Petermann et al., 2009). In our initial approach, that is adaptive binning, we kept the threshold thr and inter electrode distance Δd constant while changing Δt in accordance with the change in event rate, which in our experiments increased up to a factor of 10. This recovered power law statistics, which was confirmed by our simulation in which rate changes were created by increasing avalanche rate. These findings are in line with reported changes of the power law slope with changes in Δt (Beggs and Plenz, 2003; Petermann et al., 2009; Priesemann et al., 2014).

Our results highlight the importance of the principle of ‘separation of time scales’ (Bak et al., 1988; Vespignani and Zapperi, 1998; Plenz, 2012; Priesemann et al., 2014) for measuring avalanche dynamics. That is, the time scale at which cascading events unfold within individual avalanches, characterized by the avalanche duration, should be significantly shorter than the time scale at which different avalanches emerge, characterized by the inter-avalanche interval. Estimating the temporal resolution Δt from the inverse of the event rate provides a reasonable separation of time scales, because it balances two potential errors – false concatenation of separate avalanches and false separation of a single avalanche (Beggs and Plenz, 2003; Plenz, 2012).

In an alternative approach, we kept spatial resolution Δd and temporal resolution Δt fixed while changing threshold thr. This approach does not affect the power law organization, as demonstrated in nonhuman primates (Petermann et al., 2009) and human fMRI (Tagliazucchi et al., 2012). Particularly when Δt is difficult to change, changing thr might provide an alternative means to reestablish a separation of time scales, since large local events are less common than small local events. Our recovery of a power law for large local events during high activity periods demonstrates the analogous roles time and local event size play in avalanche dynamics. Our separate demonstration that synchrony does not change between BASE, EARLY/PRE and LATE/POST conditions further supports our finding that avalanches are maintained despite changes in event rates. Conversely, if thresholds are calculated individually for different temporal segments corresponding to different conditions of interest and the resulting distributions deviate from power laws along with concurrent changes in synchrony measures, then one has strong indication that the underlying dynamics deviate from criticality (Meisel et al., 2013).

It is important to emphasize that adaptive binning and thresholding are only tools to reveal the signatures of scale-free dynamics after they become blurred by increased drive and absence of separation of time-scales. They are not meant to be a perfect theoretical solution to the problem and have obvious limitations. For instance, there is no clear answer as to how large the window for calculating adaptive binning or thresholding should be. If it is too small, the inherent fluctuations of a critical system will be regarded as rate changes and therefore no power law dynamics can be observed. If it is too large, the rate changes can be averaged out inside them, and therefore the problem of the non-stationarity would remain. In the case of our experiments, we had a clear indication of when we should expect a rate change, which gave us a good parameter to choose the window size (~400 ms). Another problem arrives for extreme rate regimes. For very high drive, in which Δt would be much shorter than the time it takes for spikes to propagate from a neuron to its neighbors, it is unlikely that a power law can be recovered since activity becomes uncorrelated at that temporal resolution. On the other extreme, if there is almost no activity to be observed, a very large bin may require unrealistic long recordings to produce power laws. These two extreme cases may explain some of the residual error we found for left trials during PRE for monkey A (slightly larger avalanches even after adaptive binning) and during POST for monkey B (slightly smaller avalanches even after adaptive binning).

Another important point is that the advantages of criticality within a certain area allow it to optimally process information, regardless of the nature of the input received by that area. In fact, one of the advantages of criticality is maximal dynamic range (Kinouchi and Copelli, 2006; Shew et al., 2009; Gautam et al., 2015), which means that the downstream network which employs criticality will be able to respond properly for a larger range of different inputs. Nevertheless, it is possible downstream neurons employ a sort of adaptive binning mechanism. It has been suggested in simulations that the activity of an upstream network can influence how individual downstream neurons perform spatial and temporal integration at the synaptic level (Bernander et al., 1991; Rapp et al., 1992). The more input a cell receives, the more synchronized these incoming spikes have to be in order to produce a spike on the post-synaptic cell, effectively increasing the temporal resolution by which this downstream cell process its input (Bernander et al., 1991), similarly to how our proposed adaptive binning method works.

Neuronal avalanche dynamics are predominantly located in cortical layers 2/3 (Stewart and Plenz, 2006; Petermann et al., 2009) and exhibit strong non-linear components involving the interaction between cortical sites, that is the interactions of neurons or neuronal groups (Thiagarajan et al., 2010; Yu et al., 2011; Plenz, 2012). Accordingly, maintained avalanche dynamics during self-initiated and cue-related responses suggest non-linear interactions in superficial layers during cortical processing, which numerous studies have shown to be indeed the case. For example, evoked responses in superficial layers of visual cortex in the nonhuman primate have been found to be non-linear, potentially involving the local recurrent network (Williams and Shapley, 2007; Xing et al., 2012). Similarly, functional connectivity based on the interaction from neuronal firing in layer 2/3 of nonhuman primates has been found to contribute to motor coding (Hatsopoulos et al., 1998) and to be as robust or even superior over tuning curves in predicting motor outcome (Stevenson et al., 2012). With respect to their spatiotemporal spread, avalanches reveal a fractal dimension when analyzed individually (Plenz, 2012), with nearest-neighbor relationships in the average (Gireesh and Plenz, 2008; Yu et al., 2014), in line with compact spatiotemporal spreading of average evoked responses reported, for example, for premotor cortex (Rubino et al., 2006).

Our current findings may retroactively explain previous reports that found deviations from avalanche dynamics during stimulus presentation. Arviv et al. (2015) found power law avalanche size distribution when considering the full 1 s window for a visual detection task in human MEG recordings. This is consistent with our results when considering the entire period of task performance (~2 s). They also observed a trend toward supercritical dynamics during transient activity periods, just as we found predominantly large avalanches during high rate periods. Our analysis, though, shows adaptive binning as one potential approach to distinguish true supercritical dynamics from avalanche dynamics. Fagerholm et al. (2015) reported subcritical dynamics during focused attention in human EEG, in line with the observation of reduced correlations during such episodes (Cohen and Maunsell, 2009; Mitchell et al., 2009; Harris and Thiele, 2011). However, changes in activity levels were not reported in that study, and size distributions were evaluated outside their cut-off (Yu et al., 2014). Finally, early stimulus responses in the turtle’s visual system (Shew et al., 2015) have been reported to deviate from avalanches, a finding which might reflect the statistics of the stimulus used rather than intrinsic cortical dynamics.

Numerical simulations of systems under external drive or with non-zero spontaneous neuronal activity can never be truly critical, since the quiescent phase (and with it the transition connecting it to an active phase) disappears under these conditions (Taylor et al., 2013; Hartley et al., 2014; Williams-García et al., 2014). This effect is increasingly pronounced in simulations with small network size and large activity drive. Given the relatively small change in firing rate observed during behavior and the large size of the cortical network, our findings are in line with simulations that show weakly driven critical system to exhibit power law behavior when the temporal resolution is matched to the activity rate (Hartley et al., 2014).

In conclusion, our findings provide strong evidence that the scale-free organization of neuronal avalanches is maintained during evoked responses in premotor and prefrontal cortex, despite systematic fluctuations in firing rates. We therefore suggest that neuronal information processing during tasks might capitalize on various functional benefits shown for critical dynamics. This calls for further investigation into how critical dynamics may explain the execution of specific brain functions. Of equal importance, our results also shed new, methodological light on how to identify true deviations from avalanche dynamics under pathological conditions.

Materials and methods

Behavioral training and electrophysiological setup

All procedures followed the Institute of Laboratory Animal Research (part of the National Research Council of the National Academy of Sciences) guidelines and were approved by the NIMH Animal Care and Use Committee (protocol #LSN-11). Two adult rhesus monkeys (Macaca mulatta) were surgically implanted with a titanium head post each under sterile conditions while under isoflurane anesthesia. After recovery, the monkeys were trained to sit head-fixed in a primate chair for behavioral performance. In the cue-initiated task, monkey A (male, 9 years old, 8 kg) had to press a bar in front of the chair upon presentation of the ‘trial-initiation’ cue on a computer screen (a grey square in the center of the screen). After ~2 s, the initiation cue was followed by an ‘instruction’ cue, either ‘green cross’ or ‘red circle’, for the duration of 1 s. Upon cue disappearance, monkey A had to release the bar and reach with his right arm to one of two specialized feeders (Mitz et al., 2001). The ‘green cross’ instructed the monkey to reach to the left feeder for reward (‘left trials’; contralateral to the reaching hand); a ‘red circle’ instructed reaching for the right feeder (‘right trials’; ipsilateral to the reaching hand). Approaching the incorrect feeder rapidly triggered a proximity sensor to sequester the food rewards in both feeders, which prevented the monkey from obtaining a reward on that trial. The inter trial interval was 3–5 s. In the self-initiated motor task, monkey B (female, 8 years old, 7 kg) had to move her right arm to touch a pad placed ~30 cm in front of the monkey chair after which a food reward was given. Pad touching was self-initiated: no cue was presented. After the monkeys learned their respective tasks, a multi-electrode array (MEA; 96 channels - 10 × 10 without corners, inter-electrode distance: 400 μm; BlackRock Microsystems) was chronically implanted in the left prefrontal area (area 46, monkey A; electrode length: 0.55 mm) or the arm representative region of the left premotor cortex (monkey B; electrode length: 1 mm – PM is thicker than PFC, therefore longer shanks were employed). After recovering from the implantation surgery (~1 week), behavioral training resumed. The LFP (1–100 Hz band pass filtered; 2 kHz sampling frequency) and extracellular unit activity (0.3–3 kHz band pass filtered; 30 kHz sampling frequency) were simultaneously obtained from the implanted MEA. Electrophysiological signals as well as the timing of behaviorally relevant events, for example touching the pad, presentation of visual cues, etc., were stored for off-line analysis.

Unit analysis

Extracellular unit activity was extracted offline by spike sorting (Offline Sorter, Plexon Inc.). Putative single-unit activity was identified whenever a clear clustering and separation of waveforms could be identified in at least one feature space. In total, 29 and 43 putative single-units were identified in monkeys A and B, respectively. For calculating the peristimulus time histogram (PSTH), firing rates were calculated for a 250 ms window moving in steps of 50 ms. Statistical significance of firing rate changes was based on surrogate spike trains from individual units obtained by randomly permuting inter-event intervals. Specifically, a firing rate change was considered significantly above (below) expectation when it was among the top (bottom) 2.5% of the firing rate distribution obtained from the shuffled data. For this analysis, left and right trials were treated separately for monkey A.

Fano Factor

We employed the method described by Churchland et al. (2010) to calculate the Fano Factor, defined as the variance (over trials) of the spike count divided by the mean, in order to evaluate how the variability of the neuronal population studied evolved with time. We computed the variance and mean spike count for each time window (100 ms sliding window moving in steps of 40 ms) for each neuron separately. After that, a linear regression was performed on the scatterplot of the variance versus the mean spike count in which each point represents a single neuron. The slope of this regression measures the Fano Factor for the relevant time window, and the estimated error in the slope calculation provides a 95% confidence interval. To account for the effect from the increase in firing rate on the Fano Factor, a mean-match procedure was employed as described previously (Churchland et al., 2010). In short, this procedure removes units from the analysis until the mean firing rate of the remaining units does not change significantly between the periods studied (see Figures 1B and 4C; the amount of data surviving the respective mean-matching is given in the legends). Statistical significance of changes in the Fano Factor due to the task was assessed by a p-value computed from the probability of observing a given slope based on the confidence interval calculated from the baseline level.

LFP analysis

Negative deflections in the LFP (nLFPs) were detected by applying a threshold at –2 (monkey A) or –2.5 (monkey B) SD of the continuous LFP fluctuations estimated for each electrode separately (fixed and adaptive binning; for adaptive thresholding see below). The time stamps of nLFPs, determined by the data points with the largest negative amplitude, were then employed for avalanche analysis. The same procedure employed to assess the significant of firing rates was also employed for nLFP rates.

Avalanche analysis using fixed and adaptive temporal binning

Avalanches were identified by concatenating nLFPs occurring in successive time bins of width Δt on the array into temporally contiguous spatiotemporal clusters as described originally (Beggs and Plenz, 2003). The avalanche size was defined as the number of nLFPs within the avalanche. The temporal resolution Δt for the avalanche analysis was defined by the inverse of the average nLFP rate on the array over a given period and thus is not a free parameter. Two approaches for binning the data were employed. The first one, fixed binning, was obtained by calculating the average nLFP rate across all epochs (see Figures 1 and 4 for different epochs). The second one, adaptive binning, was obtained by calculating the average nLFP rate over different epochs separately. Therefore, this second method resulted in Δt that varied between epochs according to changes in the activity rate. The adaptive binning method calculates the temporal resolution on a trial-by-trial basis, taking into account highly variable responses across trials.

Significance of changes in average avalanche rate and size

In order to test for significant changes in average avalanche rate and size during task performance, we compared the values from each point in time to those obtained from 100 shuffled sets, providing a p-value. Significance was defined at p0.05. For size analysis, shuffled sets were constructed by randomly permuting avalanche sizes, while keeping avalanche occurrence time (thus also keeping the avalanche rate unchanged). For rate analysis, shuffled sets were constructed by randomly assigning a new occurrence time (within the same trial, drawn from a uniform distribution) for each avalanche. For monkey A, left and right trials were treated separately.

Avalanche size distribution analysis

In order to fit power laws and exponentials to the probability densities of avalanche size, and obtain the corresponding exponents, we employed the Maximum Likelihood Estimation method, as previously described (Clauset et al., 2009; Klaus et al., 2011). The distribution used for the power law fit had a sharp cut-off: ps=Cs-αexp[-(s/s0)γ], where the normalization constant is given by C=1/(s=sminsmaxs-αexp[-(s/s0)γ]). The three parameters to be determined were the exponent of the power law α, the cut-off s0 and the strength of the decay beyond the cut-off γ. In order to evaluate the uncertainty in the size distributions, we employed a bootstrap procedure (Efron, 1979) to estimate the 95% confidence interval for the probability density at each avalanche size from 1000 resampled data sets, which were obtained by randomly sampling the same number of avalanches found in each case from the original distribution. In the case of model distributions (see below), 100 different instances of the network were simulated, from which the confidence intervals were obtained.

The significance of the different fits to the data was obtained by employing the likelihood test: D= 2[ln(Lpl)ln(Lalt)], where D is the test statistics and Lpl is the likelihood for the power-law fit. A p-value is obtained using the calculated statistics from the chi-squared distribution. When comparing power laws to exponentials Lalt is the likelihood for the exponential fit, and the p-value computes the chances that the fit with lower likelihood value is actually better. When comparing different distributions Lalt is the likelihood for the alternative power-law fit (e.g. when comparing PRE distribution to BASE in Figure 6A right, the alternative fit is the one obtained for the BASE distribution). In this case, the p-value computes the chances that the compared distribution cannot be explained by the alternative fit or, in other words, how similar they are. Note that a p-value above 0.05 in this case does not imply that the distributions are significantly different (whereas p-values of 0.05 or lower indicate that the distributions are significantly similar).

Trial shuffling

In order to assess the influence of rate modulation introduced by the task on avalanche size distributions, we compared the original data to data obtained by randomly shuffling the trial order for each electrode independently, that is, the shift predictor (Gerstein and Perkel, 1969; 1972). This procedure created datasets that preserved the average change in nLFP rate for all electrodes and trials, but destroyed spatio-temporal correlations within individual trials.

Cellular automaton models

To study the interplay between avalanches and activity rate, we employed a network of 10 × 10 cellular automata (Kinouchi and Copelli, 2006; Ribeiro et al., 2014). Each site, which represents the activity of a single channel from the experimental data, cycles through its 11 states: xit=0 if the ith site is quiescent at time t, xit=1 if it is active, and xit=2, , 10 if it is refractory. A quiescent site at time t can become active at t + 1 if any of its pre-synaptic neighbors is active at t and transmits successfully, each connection independently with probability P. Once a site is active, its state is incremented according to the following equation until it is back to quiescence: xit+1= [xit+1] mod 11 (deterministic refractory period). Each site sends k=16 connections to randomly chosen post-synaptic sites, thus forming a random network. Each connection has a probability P of transmitting a spike, which was tuned according to P=1/k in order to achieve critical dynamics in the network. For avalanche initiation, we randomly chose a site to become active in an otherwise quiescent network. After the propagation of activity ended, a time interval drawn from the inter-avalanche interval distribution obtained from experimental data was imposed before another avalanche was initiated.

In this model, the event, that is nLFP, rate was adjusted according to the baseline level in the experimental data by employing a time step of 2 ms (implying a refractory period of 18 ms). Four different approaches to introduce a transient increase in the activity rate, as observed in the recordings from both monkeys, were studied. In the Input driven model, we added independent Poisson inputs that triggered avalanches for each site, with a rate that increases and then decreases linearly from trial time –0.35 s to –0.05 s, peaking at two different possible levels: in the low drive regime, hmax=1 s-1 and in the high drive regime, hmax=6.5 s-1. In the supercritical model, we increased the probability of transmission P by multiplying it by a modifier with the same temporal profile described above for the drive, with the modifier ranging from 1 (no change) to a maximum of either 1.5 (50% increase in the probability of spike transmission) or 2.5 (150% increase). Note that this led to both increased activity rate and supercritical dynamics. In the Poisson noise model, we added independent Poisson inputs that could not trigger avalanches to each site, with a rate similar to the one employed for the Input driven model, but peaking at either 17.5 s−1 (low-noise regime) or 35 s−1 (high-noise regime).

For avalanche analysis, the methods used for experimental data were also employed for the model. In order to reduce artificial separations of single avalanches during adaptive binning, that is when Δt becomes smaller than the time step of the model, we introduced a small jitter in event times (≤ half of the 2 ms time-step). This procedure effectively transformed the discrete temporal dynamics of the model into a continuous one. Jittering did not change our results for experimental data or simulations (data not shown).

Adaptive thresholding

For each electrode, the mean LFP amplitude and SD were calculated for consecutive windows of duration 100 ms. For each window, the adaptive threshold thr was set to mean – 2SD.

Synchronization measures

We derived estimates of mean phase synchronization for band-pass filtered data. After filtering the data (50–100 Hz frequency band; phase neutral filter by applying a second-order Butterworth filter in both directions), we first obtained a phase trace θi(t) from each LFP channel Fi(t) by applying its Hilbert transform H[Fi(t)]: θit= tan-1HFit/Fit. Next, we quantified the mean synchrony R in each segment by R= rt= L-1t=1Lrt, where L is the length of the data segment in samples and r(t) is the Kuramoto order parameter: rt= N-1j=1Neiθjt, which was used as a time-dependent measure of phase synchrony. Here, N is the number of channels in the data segment. The length of the segment in samples L is the product of the time segment considered (0.4 s) and the sampling frequency (2000 Hz), that is L=800.

References

  1. 1
  2. 2
  3. 3
  4. 4
  5. 5
  6. 6
    Neuronal avalanches in neocortical circuits
    1. JM Beggs
    2. D Plenz
    (2003)
    Journal of Neuroscience 23:11167–11177.
  7. 7
  8. 8
  9. 9
  10. 10
  11. 11
  12. 12
  13. 13
  14. 14
  15. 15
  16. 16
  17. 17
  18. 18
  19. 19
  20. 20
  21. 21
  22. 22
  23. 23
  24. 24
  25. 25
  26. 26
  27. 27
  28. 28
  29. 29
  30. 30
  31. 31
  32. 32
  33. 33
  34. 34
  35. 35
  36. 36
  37. 37
  38. 38
  39. 39
  40. 40
  41. 41
  42. 42
  43. 43
  44. 44
    Neuronal avalanches and coherence potentials
    1. D Plenz
    (2012)
    The European Physical Journal Special Topics 205:259–301.
    https://doi.org/10.1140/epjst/e2012-01575-5
  45. 45
  46. 46
  47. 47
  48. 48
  49. 49
  50. 50
  51. 51
  52. 52
  53. 53
  54. 54
  55. 55
  56. 56
  57. 57
  58. 58
  59. 59
  60. 60
    Internal state of monkey primary visual cortex (V1) predicts figure-ground perception
    1. H Supèr
    2. C van der Togt
    3. H Spekreijse
    4. VA Lamme
    (2003)
    Journal of Neuroscience 23:3407–3414.
  61. 61
  62. 62
  63. 63
  64. 64
  65. 65
  66. 66
  67. 67
  68. 68
  69. 69
  70. 70
  71. 71
  72. 72
  73. 73
  74. 74
  75. 75
  76. 76

Decision letter

  1. Cynthia Chestek
    Reviewing Editor; University of Michigan

In the interests of transparency, eLife includes the editorial decision letter and accompanying author responses. A lightly edited version of the letter sent to the authors after peer review is shown, indicating the most substantive concerns; minor comments are not usually included.

Thank you for submitting your article "Maintained avalanche dynamics during task-induced changes of neuronal activity in nonhuman primates" for consideration by eLife. Your article has been favorably evaluated by a Senior Editor and three reviewers, one of whom is a member of our Board of Reviewing Editors. The following individual involved in review of your submission has agreed to reveal his identity: Woodrow Shew (Reviewer #2). We have worked together to provide you with a consensus decision letter to help you prepare a revision.

This manuscript sets out to explore the hypothesis that task-related population dynamics in the cerebral cortex of behaving monkeys are organized according to a special type of scale-invariant spatiotemporal statistics called neuronal avalanches. This is an important hypothesis because previous work has shown that cortical states with neuronal avalanches are associated with functional advantages. However, until this paper, no one has carefully examined how and whether avalanche dynamics manifest during the execution of a well-defined behavioral task. The paper's primary, and most important finding confirms that task-related cortical activity is indeed organized as neuronal avalanches. Secondary to this, the authors have developed a new method of analyzing nonstationary data (as all task-related activity is) in the framework of neuronal avalanches.

Major Issues for Revision:

1) Issues related to the statistical robustness of the results

1a) The primary point of the paper is about the existence of neuronal avalanches during task-evoked cortical activity. However, most of the figures and work presented in the paper are focused on their new method of avalanche analysis using "adaptive binning". The authors claim that without adaptive binning, their measured avalanche size distributions deviate from power-law (subsection “Avalanche maintenance during self-initiated movements”, third paragraph and subsection “Transient rate changes during cue-triggered cognitive tasks”, last paragraph). These are some of the most important, but also weakest claims in the paper given that it is not clear they are supported by quantitative statistical measures. To a naive reader, the BASE and PRE distributions in Figure 2D appear to be extremely similar to each other. Without quantitative statistical support, the reader is left to suspect that the differences are not actually statistically significant. Perhaps they are slightly less similar to each other than those in Figure 2H, but this seems to be a very subtle difference. Without some quantitative, statistically rigorous support for this claim, much of the following work seems poorly motivated. In other words, if there was no significant deviation from power-law to begin with, then adaptive binning would seem to be unnecessary.

1b) Are the avalanches shown in Figure 1D really invariant in timescale? Is trial averaging the electrodes obscuring a lot of variance in deflection size? And is light green line the overall average? "All electrodes" could also imply that should be the more variable set of traces showing individuals.

1c) The authors suggest that the maintenance of scale invariant avalanches could also maintain the functional advantages associated with neuronal avalanches. It would greatly improve their paper if they could directly test whether function is actually maintained. Perhaps one way to do this would be to follow up on the approach taken in Figure 5B. For monkey B, they might compare some aspect of task performance (reaction time, success rate, etc.) for low rate vs. high rate trials. If function is not significantly different between these two groups of trials, it might strengthen the claim that function is maintained in spite of rate changes in avalanches. Or, if function is different for high rate vs. low rate trials, this might cast some doubt on the hypothesis that function is maintained throughout the changes in rate.

1d) In Figure 4E, there is, indeed a shift toward increased rates and smaller changes in size compared to B however there are still statistical changes in avalanche size unlike Figure 2. How do you account for this?

2) Conceptual issues:

2a) Although the dataset of Yu et al. shows rather subtle (if any) deviations from power-law avalanche distributions, other datasets, like that of Shew et al. 2015, seem to show bigger deviations from power-laws at the onset of visual stimulation, for example. This raises the question: would adaptive binning rescue the power-law distribution for all levels of increased rate, or just for sufficiently small increases in rate? This is important given that a sufficiently strong stimulus can evoke a saturated response, i.e. all or nearly all electrodes exhibiting response in unison. Can such an extremely saturated response really be scale invariant? It has one scale – the system size. We are not suggesting that the monkeys in the experiments presented here are in this saturated regime. But, in principle, this could happen, particularly in primary sensory cortex. This seems to be an important limiting case, in which adaptive binning would be expected to fail. This should be mentioned as a limitation of the adaptive binning method. This case should also be discussed as an example in which sufficiently strong input could create cortical dynamics that are not expected to be scale invariant even if the underlying network is 'truly' critical.

This is not just a conceptual point about whether cortex always supports scale invariant dynamics; it is also likely to matter for function. For example, if we are subjected to a strong enough increase in sensory input, there will be a brief transient period during which the details of the stimulus are difficult to perceive, i.e. function is suboptimal for a brief transient period. For example, leaving darkness and entering a very bright space, it takes some time (~ few hundred ms) to adjust before you can see well. In other words, not only should we expect a deviation from neuronal avalanches during such an intense input, we should also expect degraded function. As suggested by Shew et al. 2015 Nat Phys and many others, it seems that adaptation (tuning of synapses) is necessary to overcome the poor sensory processing that comes with a strongly saturated response.

2b) It is confusing when the authors reference the hard-core physics notion that criticality cannot really exist in a driven system. If that's the case, it raises the question: why do any of this analysis? Perhaps, it would be helpful to mention that the hard-core physics notion of true criticality may not be relevant given that every real neural system has some external input. The functional advantages seem to be robust to some amount of external input.

2c) It has previously been suggested that networks depart from criticality during task-evoked activity. However, given the mathematical advantages of critical dynamics, it would be curious if networks deviated precisely at the time when they should be processing task-relevant information. This paper presents a compelling explanation, which is that increases in average activity correspond to a higher avalanche frequency, and using inappropriate time bins effectively lumps multiple cascades together making the network appear to exhibit supercritical dynamics. Indeed, with shorter time bins during task-evoked increases in activity the number of large avalanches is reduced. Given all of this, it's not clear why this shorter time bin shouldn't be an appropriate size for the entire data set. In other words, if the argument is that it's only the rate and not the size of the avalanches that change, then shouldn't this smaller time bin still capture all avalanches and just show that there are fewer of them in the baseline and post-event epochs? Or is it simply more generalizable to recommend estimating the bin size based on local statistics rather than searching the data set for the optimal bin size?

2d) Related to this, adaptive binning resulted in longer time bins in the post-event epoch, but this did not normalize the avalanche statistics as it did in the pre-event epoch. At this time there were still fewer and smaller avalanches than baseline (Figure 2G). How is this explained? Is this evidence against criticality?

3) Issues concerning interpretation

3a) Reinterpreting large avalanches as multi-avalanches is a reasonable and interesting thing to try, but was guaranteed to shift size into rate as shown. Therefore, the main result is really the tighter power law curves. However, we would like to see a more careful consideration of alternative hypotheses?

3b) Somewhere in the Introduction or Discussion, it would be beneficial to more clearly explain the difference between a system with 'true' criticality but without scale invariant dynamics versus a truly non-critical system. Both would lack scale invariant dynamics. The hypothesis proposed by the authors is that their experiments are a case of the former (true criticality without power-law avalanches). However, the alternative, competing hypothesis is also quite plausible. That is, there are many well-known ways that interactions among neurons can change (synaptic depression, facilitation, neuromodulators, etc.) and therefore true deviations from criticality might result. I think the paper would be stronger if this alternative was explained.

3c) From the perspective of a researcher, the adaptive binning method seems to be a useful tool for learning something about the state of the cortex in a way that is less dependent on activity rate nonstationarities. However, adaptive binning may be irrelevant from the perspective of a downstream neuron that receives input from a region with changing rate. Unless the downstream neuron has some way of adapting its own synaptic integration timescales, it seems like it would see input that deviates from scale-free avalanches.

3d) The authors suggest that the functional benefits associated with neuronal avalanches might be preserved, but it seems that, functionally, real neurons do not have the ability to do the adaptive binning they would need to get the functional benefits of neuronal avalanches. Moreover, rate coding is prominent in many neural systems. If neurons do implement some kind of adaptive binning, it seems they would be ignoring rate coded information. Please discuss this in the Discussion section.

https://doi.org/10.7554/eLife.27119.011

Author response

Major Issues for Revision:

1) Issues related to the statistical robustness of the results

1a) The primary point of the paper is about the existence of neuronal avalanches during task-evoked cortical activity. However, most of the figures and work presented in the paper are focused on their new method of avalanche analysis using "adaptive binning". The authors claim that without adaptive binning, their measured avalanche size distributions deviate from power-law (subsection “Avalanche maintenance during self-initiated movements”, third paragraph and subsection “Transient rate changes during cue-triggered cognitive tasks”, last paragraph). These are some of the most important, but also weakest claims in the paper given that it is not clear they are supported by quantitative statistical measures. To a naive reader, the BASE and PRE distributions in Figure 2D appear to be extremely similar to each other. Without quantitative statistical support, the reader is left to suspect that the differences are not actually statistically significant. Perhaps they are slightly less similar to each other than those in Figure 2H, but this seems to be a very subtle difference. Without some quantitative, statistically rigorous support for this claim, much of the following work seems poorly motivated. In other words, if there was no significant deviation from power-law to begin with, then adaptive binning would seem to be unnecessary.

We thank the reviewers for these comments. In fact, significance was demonstrated in our original submission using a boot-strap method. This provided shaded confidence intervals, which unfortunately, made differences visually less distinctive. In response, we now provide separate figures for each monkey that detail differences in size distributions and we added additional quantifications of differences (new Figures 3 and 6). These figures now show each size distribution normalized by the expected power law. Thus, the expected power law translates into a horizontal curve, which allows power-law deviations to be visualized clearly. The area under each normalized distribution now quantifies how much LATE and PRE distributions are biased towards larger avalanches. Finally, we present the actual distribution of differences in the mean avalanche size of LATE/PRE and BASE, for each trial, using fixed and adaptive binning. Comparing these distributions demonstrate the significant effect in adaptive binning in recovering a power law. As a minor point, we now present the distributions of the monkey performing the cognitive task first.

1b) Are the avalanches shown in Figure 1D really invariant in timescale? Is trial averaging the electrodes obscuring a lot of variance in deflection size? And is light green line the overall average? "All electrodes" could also imply that should be the more variable set of traces showing individuals.

We apologize for the confusion. Figure 1D (now Figure 4C) displayed the average, touch-pad aligned LFP waveform for each single electrode superimposed and the corresponding average (center line). These data support the extraction of negative LFP threshold crossings (nLFPs) to identify nLFP rate changes that identify BASE, PRE and POST periods. We have reworded the figure legends accordingly to avoid this confusion.

1c) The authors suggest that the maintenance of scale invariant avalanches could also maintain the functional advantages associated with neuronal avalanches. It would greatly improve their paper if they could directly test whether function is actually maintained. Perhaps one way to do this would be to follow up on the approach taken in Figure 5B. For monkey B, they might compare some aspect of task performance (reaction time, success rate, etc.) for low rate vs. high rate trials. If function is not significantly different between these two groups of trials, it might strengthen the claim that function is maintained in spite of rate changes in avalanches. Or, if function is different for high rate vs. low rate trials, this might cast some doubt on the hypothesis that function is maintained throughout the changes in rate.

We thank the reviewer for this important question given the performance task we analyzed. In fact, we test both trial outcome and reaction time. Regarding trial outcome, the monkey’s success was > 95% leaving us with too few false trials to establish statistical significance. We then extracted reaction time as the period from the visual cue being switched off to reaching the feeder. There was no correlation between the number of nLFPs on the array and reaction time. Similarly, no correlation was found for avalanche size and rate using adaptive or fixed binning. These results, though negative, demonstrate that adaptive binning does introduce new correlations and have now been added to the manuscript (subsection “Avalanche maintenance during a cue-triggered cognitive task”, last paragraph).

We have also added new paragraphs in the Discussion (seventh and eighth paragraphs) which elaborates on the issue that the dynamical state (as identified by power law statistics) is maintained in the monkeys independent from task-induced activity rate changes. It’s important that we emphasize that our message is that the dynamical state is maintained regardless of rate changes. That doesn’t necessarily imply that the performance doesn’t change with rate: the fact that some trials have decreased rate could be linked to decreased input, which could be, for instance, a reflection of lack of attention, which would lead to worse performance, regardless of the dynamical state of the network in pre-frontal cortex. The fact that this network is critical means it can do the best possible processing given the input it received. The only true way of testing the link between performance and dynamical state is by disturbing the dynamical state of the network we’re investigating.

1d) In Figure 4E, there is, indeed a shift toward increased rates and smaller changes in size compared to B however there are still statistical changes in avalanche size unlike Figure 2. How do you account for this?

In our new Figure 3E we now show the statistics of significantly reduced bias towards larger avalanches for adaptive binning. This figure also quantifies the residual error this referee pointed out. We now discuss potential sources of this residual error in the manuscript. We believe one major source originates from the necessity to parse the behavioral continuous time course into three episodes from average trial activity for which adaptive bins are calculated. This error is inherent to any approach of estimating rate from a non-zero time period. A second error arises when there are too few events, for e.g. during the POST epoch, when activity is bounded by a zero-activity floor. These limitations are now discussed in a separate paragraph in the Discussion (seventh paragraph).

2) Conceptual issues:

2a) Although the dataset of Yu et al. shows rather subtle (if any) deviations from power-law avalanche distributions, other datasets, like that of Shew et al. 2015, seem to show bigger deviations from power-laws at the onset of visual stimulation, for example. This raises the question: would adaptive binning rescue the power-law distribution for all levels of increased rate, or just for sufficiently small increases in rate? This is important given that a sufficiently strong stimulus can evoke a saturated response, i.e. all or nearly all electrodes exhibiting response in unison. Can such an extremely saturated response really be scale invariant? It has one scale – the system size. We are not suggesting that the monkeys in the experiments presented here are in this saturated regime. But, in principle, this could happen, particularly in primary sensory cortex. This seems to be an important limiting case, in which adaptive binning would be expected to fail. This should be mentioned as a limitation of the adaptive binning method. This case should also be discussed as an example in which sufficiently strong input could create cortical dynamics that are not expected to be scale invariant even if the underlying network is 'truly' critical.

This is not just a conceptual point about whether cortex always supports scale invariant dynamics; it is also likely to matter for function. For example, if we are subjected to a strong enough increase in sensory input, there will be a brief transient period during which the details of the stimulus are difficult to perceive, i.e. function is suboptimal for a brief transient period. For example, leaving darkness and entering a very bright space, it takes some time (~ few hundred ms) to adjust before you can see well. In other words, not only should we expect a deviation from neuronal avalanches during such an intense input, we should also expect degraded function. As suggested by Shew et al. 2015 Nat Phys and many others, it seems that adaptation (tuning of synapses) is necessary to overcome the poor sensory processing that comes with a strongly saturated response.

We agree with the reviewer and have included a discussion of saturated responses to the method in the appropriate section (Discussion, seventh paragraph).

2b) It is confusing when the authors reference the hard-core physics notion that criticality cannot really exist in a driven system. If that's the case, it raises the question: why do any of this analysis? Perhaps, it would be helpful to mention that the hard-core physics notion of true criticality may not be relevant given that every real neural system has some external input. The functional advantages seem to be robust to some amount of external input.

The reviewer is correct and we have modified the main text according to his/her suggestions (Discussion, eleventh paragraph).

2c) It has previously been suggested that networks depart from criticality during task-evoked activity. However, given the mathematical advantages of critical dynamics, it would be curious if networks deviated precisely at the time when they should be processing task-relevant information. This paper presents a compelling explanation, which is that increases in average activity correspond to a higher avalanche frequency, and using inappropriate time bins effectively lumps multiple cascades together making the network appear to exhibit supercritical dynamics. Indeed, with shorter time bins during task-evoked increases in activity the number of large avalanches is reduced. Given all of this, it's not clear why this shorter time bin shouldn't be an appropriate size for the entire data set. In other words, if the argument is that it's only the rate and not the size of the avalanches that change, then shouldn't this smaller time bin still capture all avalanches and just show that there are fewer of them in the baseline and post-event epochs? Or is it simply more generalizable to recommend estimating the bin size based on local statistics rather than searching the data set for the optimal bin size?

The issue with using a single bin to resolve all data is the non-stationarities, and the explanation is in (old) Figure 5A. Using a small bin for all the data would be equivalent of setting the bin at the orange level, let’s assume. That would result in power laws with different exponents which when combined would lead to bumps as observed in the grey curve. It’s important to notice that the bin is calculated for every trial, and not just one for each epoch. We added this explanation to the main text (subsection “Impact of ‘adaptive binning’ in capturing rate change vs. change in dynamics”, first paragraph).

2d) Related to this, adaptive binning resulted in longer time bins in the post-event epoch, but this did not normalize the avalanche statistics as it did in the pre-event epoch. At this time there were still fewer and smaller avalanches than baseline (Figure 2G). How is this explained? Is this evidence against criticality?

We agree with this reviewer and have now added an extensive paragraph in the Discussion about the limitation of adaptive binning. In general, saturation, when responses cannot be discriminated anymore and very low-activity, when even long time bins can’t concatenate events into clusters for the recording duration given, provide upper and lower limits respectively (see Discussion, seventh paragraph).

3) Issues concerning interpretation

3a) Reinterpreting large avalanches as multi-avalanches is a reasonable and interesting thing to try, but was guaranteed to shift size into rate as shown. Therefore, the main result is really the tighter power law curves. However, we would like to see a more careful consideration of alternative hypotheses?

In fact, our simulations considered three alternative (non-critical) hypotheses (Figure 7): supercritical, added noise, and faster propagation (not shown in the manuscript due to non-physiologically high propagation rates to obtain observed rate changes). All of them failed to reproduce the observed experimental results. We also provided two additional methods (Figure 8; end of Results section): adaptive thresholding and synchronization measure, which are completely separated and independent of temporal binning. All are in line with the results observed for adaptive binning supporting that the increased drive is the most likely explanation for the data. We believe that the revision of the manuscript and additional analysis provided have increased the transparency of our work such that these alternative considerations stand out more.

3b) Somewhere in the Introduction or Discussion, it would be beneficial to more clearly explain the difference between a system with 'true' criticality but without scale invariant dynamics versus a truly non-critical system. Both would lack scale invariant dynamics. The hypothesis proposed by the authors is that their experiments are a case of the former (true criticality without power-law avalanches). However, the alternative, competing hypothesis is also quite plausible. That is, there are many well-known ways that interactions among neurons can change (synaptic depression, facilitation, neuromodulators, etc.) and therefore true deviations from criticality might result. I think the paper would be stronger if this alternative was explained.

We agree with the reviewer that is fundamentally important to discuss the alternative that a true deviation from criticality is the reason for the deviation towards larger avalanches in the data. I believe we addressed this point by providing 8 references in the Introduction (last paragraph) before we present our hypothesis as alternative view. We used the Results section to detail our simulations in which true deviations from criticality are explored further. Accordingly, we show in these simulations that in a truly supercritical system, our approach of adaptive binning leads to even larger deviations from power law distributions, which is in absolute contrast with the experimental data. See point 3a above. We believe our current layout of addressing non-critical dynamics and our suggested approach works well given the space requirements and flow of the manuscript.

3c) From the perspective of a researcher, the adaptive binning method seems to be a useful tool for learning something about the state of the cortex in a way that is less dependent on activity rate nonstationarities. However, adaptive binning may be irrelevant from the perspective of a downstream neuron that receives input from a region with changing rate. Unless the downstream neuron has some way of adapting its own synaptic integration timescales, it seems like it would see input that deviates from scale-free avalanches.

We thank the reviewer for this insightful comment. We now reference work in which neurons have been suggested to reduce their integration time constant as they receive more input, which essentially mimics adaptive binning as a physiological process. We have added this point to the Discussion section (eighth paragraph).

3d) The authors suggest that the functional benefits associated with neuronal avalanches might be preserved, but it seems that, functionally, real neurons do not have the ability to do the adaptive binning they would need to get the functional benefits of neuronal avalanches. Moreover, rate coding is prominent in many neural systems. If neurons do implement some kind of adaptive binning, it seems they would be ignoring rate coded information. Please discuss this in the Discussion section.

We thank the reviewer for this comment. We added that to the Discussion section (see point 3c above).

https://doi.org/10.7554/eLife.27119.012

Article and author information

Author details

  1. Shan Yu

    Section on Critical Brain Dynamics, National Institute of Mental Health, Bethesda, United States
    Present address
    Institute of Automation, and Center for Excellence in Brain Science and Intelligence Technology, Chinese Academy of Sciences, Beijing, China
    Contribution
    Conceptualization, Data curation, Software, Formal analysis, Supervision, Validation, Investigation, Visualization, Methodology, Writing—original draft, Writing—review and editing
    Contributed equally with
    Tiago L Ribeiro
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-9008-6658
  2. Tiago L Ribeiro

    Section on Critical Brain Dynamics, National Institute of Mental Health, Bethesda, United States
    Contribution
    Software, Formal analysis, Investigation, Visualization, Methodology, Writing—original draft, Writing—review and editing
    Contributed equally with
    Shan Yu
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0003-3195-9284
  3. Christian Meisel

    Section on Critical Brain Dynamics, National Institute of Mental Health, Bethesda, United States
    Contribution
    Software, Formal analysis, Visualization, Methodology, Writing—review and editing
    Competing interests
    No competing interests declared
  4. Samantha Chou

    Section on Critical Brain Dynamics, National Institute of Mental Health, Bethesda, United States
    Contribution
    Data curation, Writing—review and editing
    Competing interests
    No competing interests declared
  5. Andrew Mitz

    Laboratory of Neuropsychology, National Institute of Mental Health, Bethesda, United States
    Contribution
    Resources, Data curation, Methodology
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-8045-0970
  6. Richard Saunders

    Laboratory of Neuropsychology, National Institute of Mental Health, Bethesda, United States
    Contribution
    Resources, Data curation, Methodology
    Competing interests
    No competing interests declared
  7. Dietmar Plenz

    Section on Critical Brain Dynamics, National Institute of Mental Health, Bethesda, United States
    Contribution
    Conceptualization, Resources, Software, Supervision, Funding acquisition, Validation, Investigation, Visualization, Methodology, Writing—original draft, Project administration, Writing—review and editing
    For correspondence
    plenzd@mail.nih.gov
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-0008-3657

Funding

National Institute of Mental Health

  • Shan Yu
  • Tiago L Ribeiro
  • Christian Meisel
  • Samantha Chou
  • Andrew Mitz
  • Richard Saunders
  • Dietmar Plenz

Conselho Nacional de Desenvolvimento Científico e Tecnológico

  • Tiago L Ribeiro

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

We thank the members of the Plenz lab for helpful discussion and John Hylton for help in spike sorting. This work was supported by the Intramural Research Program of the National Institute of Mental Health (NIMH) and Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq).

Ethics

Animal experimentation: All procedures followed the Institute of Laboratory Animal Research (part of the National Research Council of the National Academy of Sciences) guidelines and were approved by the NIMH Animal Care and Use Committee.(protocol #LSN-11).

Reviewing Editor

  1. Cynthia Chestek, University of Michigan

Publication history

  1. Received: March 23, 2017
  2. Accepted: October 28, 2017
  3. Version of Record published: November 8, 2017 (version 1)

Copyright

This is an open-access article, free of all copyright, and may be freely reproduced, distributed, transmitted, modified, built upon, or otherwise used by anyone for any lawful purpose. The work is made available under the Creative Commons CC0 public domain dedication.

Metrics

  • 578
    Page views
  • 92
    Downloads
  • 1
    Citations

Article citation count generated by polling the highest count across the following sources: Crossref, PubMed Central, Scopus.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Download citations (links to download the citations from this article in formats compatible with various reference manager tools)

Open citations (links to open the citations from this article in various online reference manager services)