Neuronal Avalanches Across the Rat Somatosensory Barrel Cortex and the Effect of Single Whisker Stimulation

1. Introduction

The cortex operates in a state of restless activity, whose meaning and functionality are not yet understood. The critical brain hypothesis suggests that this is the result of the brain operating in the vicinity of the critical point of a phase transition, leading to a rich and variable dynamics at rest. In general, it has been argued that criticality provides biological systems with an optimal balance between robustness against perturbations and the flexibility to adapt to changing conditions. In the case of the brain, this would confer optimal computational capabilities [e.g., by optimizing the correlation length and the dynamic range, leading to the existence of large dynamical repertoires accompanied by maximal transmission and storage of information (Kinouchi and Copelli, 2006; Shew et al., 2009; Shew and Plenz, 2013; Hidalgo et al., 2014; Muñoz, 2018)]. In this context, Hidalgo et al. (2014) have shown that complex adaptive systems that have to cope with a great variety of stimuli are much more efficient when operating in the vicinity of a critical point, and thus they benefit from dynamically tuning themselves to that point.

By analyzing LFPs of cortical neurons in culture, the seminal work of Beggs and Plenz (2003) provided the first evidence of power law distributed neuronal avalanches, i.e., cascades of activity interspersed by periods of quiescence typical of critical systems. In particular, the exponents of these power laws were remarkably close to the ones of a critical branching process, hence suggesting a scenario in which neuronal networks are characterized by a marginal propagation of the activity at the critical point between an active and an absorbing phase. Since then, such power laws have been observed repeatedly in different experimental settings (Mazzoni et al., 2007; Gireesh and Plenz, 2008; Pasquale et al., 2008; Petermann et al., 2009; Hahn et al., 2010; Yu et al., 2011), thus strengthening the critical brain hypothesis. Despite that, ambiguities, and open questions remain.

A first problem concerns the experimental definition of avalanche. In order to estimate avalanches, one needs to define discrete events. While neuronal spikes and multi-unit activities (MUAs) are events by nature, the conversion of coarse-sampled brain signals such as LFPs into a discrete form (e.g., by simply applying a threshold) is more ambiguous and difficult to interpret in terms of neural correlates. As the relation between events and actual underlying neural activity becomes more uncertain, the definition of neural avalanches becomes fuzzy (Touboul and Destexhe, 2010, 2017; Dalla Porta and Copelli, 2019; Villegas et al., 2019).

An additional limiting factor is spatial sampling of experimental recordings. Studies using coarse-sampled activity like LFPs typically yielded power-law distributions both in vivo and in vitro, but several experiments relying on spikes in awake animals did not (Ribeiro et al., 2010; Dehghani et al., 2012; Priesemann et al., 2014). One possible cause is insufficient spatial sampling of recording in awake conditions. Avalanches are a population phenomenon, but spikes are sparsely recorded in these experiments and reflect a subpopulation of neurons, in fact missing a consistent fraction of the real activity (Ribeiro et al., 2010). On the other hand, LFPs are average and composite signals that reflect neuronal populations but are difficult to interpret in terms of single neurons activity (Einevoll et al., 2013; Wilting and Priesemann, 2019; Neto et al., 2020). This makes the event-extraction process less straightforward and threshold dependent (Touboul and Destexhe, 2010), even though nLPFs are known to correlate with synchronous spiking activity (Beggs and Plenz, 2003; Petermann et al., 2009; Yu et al., 2011). For all these reasons, it is informative to include spikes, MUAs and LFPs when measuring neural avalanches (see for example Hahn et al., 2010).

Up to now, most of the work on neuronal avalanches has focused on spontaneous activity, while much remains to be understood about their behavior after perturbations caused by incoming inputs. Investigating the response to sensory stimuli in primary cortical areas is a clear-cut strategy to address this point. First, these brain regions can be expected to benefit from operating around a critical point to encode the sensory stimuli themselves. Moreover, sensory inputs, which are under direct experimental control, propagate to cortical networks across a limited number of well characterized processing stages (contrary to, e.g., associative or motor areas). Avalanches were studied in the turtle visual cortex ex-vivo while the retina was exposed to a movie acting as a continuous visual stimulus (Shew et al., 2015; Clawson et al., 2017). Based on neural avalanches extracted from LFPs, results suggested that the cortical network self-adapts to a critical state after a short period from the stimulus. However, in these experiments, LFPs were sparsely measured and reflected populations of neurons scattered across the visual processing cascade and downstream to important processing structures including the retina. Work on the primary auditory cortex hinted, instead, at a critical behavior both in resting and post-stimulus conditions (Bowen et al., 2019). Noteworthy, measurements were confined to either layers 2/3 or 4 and limited by the slow dynamics of calcium imaging to obtain an indirect estimate of neural activity.

In this work we contribute to verify the critical brain hypothesis with a systematic study of neuronal avalanches in the rat barrel cortex (the region of the rat primary somatosensory cortex that encodes tactile sensory inputs from the whiskers). We run our measurements across cortical layers in single barrel columns of the rat anesthetized with tiletamine. This common preparation for electrophysiology (Sorrenti et al., 2021) displays rich cortical spontaneous activity, including UP and DOWN states and oscillations that have been linked to avalanches and criticality (Scarpetta and de Candia, 2014; di Santo et al., 2018). Given the current challenges to test the critical brain hypothesis and the various nature of recorded signals (Ribeiro et al., 2010; Dehghani et al., 2012; Priesemann et al., 2014), we explored activity across a wide frequency range, covering spikes, MUAs and LFPs (i.e., up to 3,000 Hz). In particular, we investigated both spontaneous activity and the evoked activity that followed a single whisker deflection and analyzed neural avalanches through a protocol based on state-of-the-art maximum likelihood statistical method (Gerlach and Altmann, 2019).

2. Materials and Methods

2.1. Electrophysiological Recordings and Surgical Procedures

Extracellular spikes and MUAs and LFPs recordings. Spikes and MUAs were recorded using a neural probe with a linear array of 32 Iridium Oxide (IrOx) microelectrodes with 65 μm inter-electrode distance and 2015 μm array length (E32+R-65-S1-L6 NT; Atlas Neuroengineering) (Figure 1). Raw signals were acquired by an Open Ephys Acquisition Board (OEps Tech, Lisbon, Portugal) with a 32-channels head stage (RHD2000, Intan Technologies) and an SPI cable at 25 KHz sampling frequency and band-pass filtered (300–3,000 Hz). Single unit activity (spikes) represented only about one percent of the activity recorded which was, therefore, mostly represented by MUAs. The maximum number of neurons contributing to a MUA is of a few tens per microelectrode as it can be estimated assuming a spherical volume of tissue centered around the recording site, an average density of neurons in the rat barrel cortex ~ 55,000/mm3 from Ren et al., 1992) and and approximate distance of one hundred micrometers at which the amplitude of an extracellular action potential is above the background noise (~ 10 microvolts RMS in the 0.1–3 kHz bandwidth). LFPs were recorded at high density using a CMOS based neural probe with an array of 256 microelectrodes (7.4 μm in diameter size and organized in four vertical columns and sixty-four horizontal rows) (Schroder et al., 2015) with 32 μm inter-electrode distance along both the horizontal and vertical axis and 2023.4 μm array length (Figure 1). Raw multiplexed signals were acquired through a NI PXIe-6358 (National Instruments) board (sampling frequency 1.25 MS/s at 16 bit) and demultiplexed using a custom-made LabVIEW software. The resulting whole-array LFP signal was sampled at 976.56 Hz and band-pass filtered (2–300 Hz). Once inserted in the barrel column, both arrays were spanning across all the six cortical layers (from 0 to 1,800 μm).

FIGURE 1

Figure 1. (A) Experimental setting for recording with the anesthetized rat immobilized on the stereotaxic apparatus. The Atlas probe for spikes and MUAs recording is inserted in the barrel cortex through the dedicated cranial window. The cortical surface is bathed by Krebs' solution that is grounded through the immersed Ag/AgCl reference electrode. The piezoelectric bender with cannula used to control single whisker deflection is visible on the left. An identical arrangement, but with the custom high-density probe, was adopted for LFPs recording. (B) Custom 2D 64 × 4 array used for LFPs and commercial 1D array used for MUAs spanning across cortical layers as during recording (schematic drawing). The actual number of electrodes spanning the barrel cortex was 220 for LFPs (organized in a 55 × 4 matrix form) and 27 for MUAs, with the other electrodes of the arrays remaining above the cortical surface bathed by the grounded electrolyte. Examples of MUAs (C) and LFPs (D) traces (from three representative channels each).

Surgical implantation and single whisker stimulation. Wistar rats were maintained under standard environmental conditions in the animal research facility of the Department of Biomedical Sciences - University of Padova. All the procedures were approved by the local Animal Care Committee (O.P.B.A.) and the Italian Ministry of Health (authorization number 522/2018-PR). Rats of both sexes, aged 36–50 days (P36-P50) and weighting between 150 and 230 g, were anesthetized with an intra-peritoneal induction mixture of tiletamine-xylazine (2 mg and 1.4 g/100 g body weight, respectively), followed by additional doses (0.5 mg and 0.5 g/100 g body weight) every hour. The anesthesia level was constantly monitored by testing the absence of eye and hind-limb reflexes and whiskers' spontaneous movements. Before starting with surgery, the rat was fixed on a stereotaxic apparatus by teeth and ear bars. The body temperature was monitored continuously with a rectal probe and maintained at 37°C by a heating pad. The skull was exposed through an anterior-posterior opening of the skin in the center of the head and a window was drilled over the right somatosensory barrel cortex at stereotaxic coordinates from −1 to −4 AP, from +4 to +8 ML referred to bregma (Swanson, 2003). A slit in the meninges was made with dedicated fine forceps at coordinates −2.5 AP, +6 ML for the subsequent insertion of the recording probe, and the brain was constantly bathed in Krebs' solution (in mM: NaCl 120, KCl 1.99, NaHCO₃ 25.56, KH₂PO₄ 136.09, CaCl₂ 2, MgSO₄ 1.2, glucose 11).

The recording probe was fixed to a dedicated holder connected to a Patchstar micromanipulator (Scientifica Ltd, East Sussex, UK), which was used for inserting the probe into the cortex orthogonal to the cortical surface. The depth was set at 0 μm when the electrode proximal to the chip tip touched the cortical surface. An Ag/AgCl electrode bathed in Krebs' solution in proximity of the probe was used as reference.

Contralateral whiskers were trimmed at around 10 mm from the mystacial pad. To control deflection, single whiskers were inserted for 8 mm inside a cannula glued to a piezoelectric bender with integrated strain gauges (P-871.122; Physik Instrumente (PI) GmbH & Co. KG) and driven by a custom-made digital closed-loop control system. Each stimulus, delivered by a waveform generator (Agilent 33250A 80 MHz, Agilent Technologies Inc., Colorado, USA), was consisting of a voltage pulse of 5 ms duration and 100 μs rise/fall time applied to the piezoelectric bender. The principal (maximally responding) whisker identified on the basis of the amplitude of evoked LFP responses was selected for the recording session.

2.2. Avalanches Analysis

MUAs and high-density LFPs recordings were performed separately in four and five rats, respectively. The minimum time interval between consecutive whisker stimuli was set to 2 s to avoid receptors and central adaptation phenomena. Accordingly, 2 s of recording after the stimuli were excluded from the analysis of basal activity. Due to the different duration of the stimulus-evoked avalanches in the LFPs and MUAs domains (Figure 2), post-stimulus intervals were set for analysis at 2 s or 500 ms for LFPs and MUAs, respectively. In total, 2 min long recordings of LFPs basal activity and 5 min long recordings of MUAs basal activity were analyzed for each rat. Forty stimulations of the whisker were considered for each rat in the analysis of LFPs evoked activity, while in MUAs data the recordings of each rat included at least 60 stimulations of the whisker.

FIGURE 2

Figure 2. Evoked response by whisker stimulation: LFPs and MUAs. Example of LFP (A) and spike (B) traces with the 3 standard deviations (SD) thresholds used for events detection (indicated by dots). The corresponding raster plots are shown in (C,D). The 2D array in (C) was rearranged in one dimension (groups of four sensors in a row are reported consecutively along the vertical axis). The early response to the stimulus occurred within few tens of milliseconds and was characterized by a high degree of events synchronization across channels both for LFPs and MUAs. The typical LFP evoked response consisted in a first negative peak followed by a slower positive wave. Both were detected as events which is reflected in the raster plot. After the early response large and frequent avalanches typically followed, taking the form of oscillations of synchronous events across channels fading after a few hundreds of milliseconds from the stimulus. Correspondingly, epochs of high frequency firing (bursts) were observed in the MUAs domain (C,D).

For the detection of LFP events, the standard deviation (SD) and the mean of the signal was computed for each channel. In order to distinguish real events from noise, a three SD threshold was chosen basing on the distribution of the signal amplitudes which significantly deviated from a Gaussian best fit above that threshold (see Supplementary Figure 1). Both negative and positive LFPs (i.e., nLPFs adn pLFPs, respectively) were considered as events in accordance with previous work (Shew et al., 2015; Clawson et al., 2017). Within our specific experimental settings, one reason is that across the depth of the cortex there are polarity changes of the LFP signal because of compensatory capacitive ionic currents particularly along dendrites of pyramidal cells (Buzsaki et al., 2012; Einevoll et al., 2013). Since in our experiments electrodes span multiple cortical layers, both nLFPs and pLFPs were found and detected. Moreover, alternatively, pLFPs can be related to activation of populations of inhibitory neurons leading to inhibitory outward postsynaptic currents, which also justifies their inclusion in the events count. For detection, each deflection was considered terminated only after it crossed the mean of the signal. For completeness, we report analyses on nLFPs only, and also considering separately the superficial layers of the cortex and the bottom layers (Supplementary Figures 8–13).

For extracellular spikes and MUAs the detection threshold was set at three SD of the noise (Quiroga et al., 2004). In fact, this threshold excludes the Gaussian “white noise” component from the signal (Supplementary Figure 2). Events recorded at the same time frame (temporal resolution: 0.04 ms) by different microelectrodes were ascribed to the same neuron and thus counted as one event, although avalanches results were not significantly affected by this correction.

For both post stimulus and basal (resting state) activity, an average inter event interval (〈IEI〉) was calculated and used for temporal binning to estimate avalanches. Since the 〈IEI〉 is calculated in the periods post-stimuli separately from the periods of resting state, our procedure is analogous to adaptive binning (Yu et al., 2017). Avalanches were defined as sequences of 〈IEI〉 time bins presenting activity in the form of events, with the end of the avalanche identified by the first empty bin. The number of events in each avalanche accounted for its size, while the duration was the number of temporal bins comprising the avalanche.

2.3. Power Law Fitting and Statistical Testing

The avalanches sizes and durations distributions are fitted using the maximum likelihood method. The fitting function for both avalanche sizes and duration is a discrete power-law:

$\begin{array}{ll}p(y;\alpha )=\frac{y^{-\alpha }}{{\displaystyle {\sum }_{x=x_{min}}^{x=x_{max}}}{x}^{-\alpha }}.& (1)\end{array}$

The parameter x_max is set to the maximum observed size or duration. Then the tails of the distributions are fitted by selecting as parameter x_min the one that minimizes the Kolmogorov-Smirnov distance (KS), following the method proposed by Clauset et al. (2009):

$\begin{array}{ll}KS={\displaystyle \underset{y\ge x_{min}}{max}}|S(y)-\widehat{P}(y)|& (2)\end{array}$

where S(y) is the cumulative distribution function (CDF) of the data and $\widehat{P}(y)$ is the CDF of the theoretical distribution fitted with the parameter that best fits the data for y ≥ x_min.

After finding the best-fit power law, to assess goodness-of-fit we compared the experimental data against 1,000 surrogate datasets drawn from the best-fit power law distribution with the same number of samples as the experimental dataset. The deviation between the surrogate datasets and a perfect power law was quantified with the KS statistic. The p-value of the power-law fit was defined as the fraction of these surrogate KS statistics which were greater than the KS statistic for the experimental data. Note that the data were considered power law distributed if the null hypothesis could not be rejected, namely if the the p-value turned out to be greater than the significance level, which was set to a conservative value of 0.1.

However, when estimating the parameters and evaluating the p-value, we take into consideration another aspect that has been recently pointed out in Gerlach and Altmann (2019). A point often ignored is that maximum likelihood methods rely on two assumptions:

1. the observations y are distributed as p(y; α), where α is the power law exponent;

2. the empirical observations y_i, i = 1, …, N, are independent.

While the first assumption corresponds to our choice of a statistical law, statistical tests rely on the second one, which for instance is implicitly assumed when the log-likelihood is computed as $\sum _{i=1}^{i=N}logp(y_i)$. However, complex systems are often characterized by strong temporal and spatial inter-dependencies, thus often violating the independence assumption. This may lead to false rejections of the statistical laws and to over-optimistic uncertainties of the estimated parameters. The authors of Gerlach and Altmann (2019) propose a method to distinguish between these assumptions, and we exploit it here to estimate the parameters and evaluate the goodness-of-fit. Briefly, we take the timeseries of sizes or durations of consecutive avalanches, and we estimate the time τ^* after which two observations (e.g., the avalanche sizes) are independent from each other. In practice, τ^* is obtained by computing the time at which the autocorrelation of the timeseries reaches an interval around zero (1-percentile of the random realization). Then, the original sequence of length N is randomized, N^* = N/τ^* observations are selected and the standard statistical analysis is applied to the new sample. This guarantees that the new sample of dimension N^* < N comprises only uncorrelated avalanches.

Indeed, in our timeseries of sizes and durations we find non negligible values of τ^* (see Supplementary Figure 4), and we verify that the acceptance rate of the statistical law increases for uncorrelated avalanches. Because of the variability of the different realizations of the subsampling procedure, the exponents and the p-values shown in Tables 1, 2 are obtained averaging over 20 repetitions of this subsampling.

TABLE 1

Table 1. Exponents before and after subsampling in LFP data, with the corresponding p-values.

TABLE 2

Table 2. Exponents before ad after subsampling in MUAs data, with the corresponding p-values.

3. Results

3.1. Preliminary Considerations

The standard approach to infer criticality is to search for neuronal avalanches whose sizes and durations follow scale-free distributions in resting state (i.e., unperturbed) conditions. Indeed, at criticality, it is expected that such distributions scale as the power laws

$\begin{array}{ll}P(S)~S^{-\tau },& (3)\end{array}$

$\begin{array}{ll}P(T)~T^{-{\tau }_t},& (4)\end{array}$

where S is the number of events in an avalanche (i.e., the size), T is its duration (also called avalanche lifetime) and τ and τ_t are the related critical exponents. Power laws, however, can also stem from non critical systems and generative mechanisms. Thus, a more robust test of criticality is to verify whether the so-called crackling noise relation holds. This scaling relation was first developed in the context of crackling noise (Sethna et al., 2001), hence the name, but nonetheless it is expected to hold in general in all systems close to their critical point (Friedman et al., 2012), and in particular in systems with absorbing states (di Santo et al., 2017). The relation predicts that the critical exponent δ, which relates the duration of an avalanche to its mean size as

$\begin{array}{ll}\langle S\rangle (T)~T^{{\delta }_{\text{fit}}}& (5)\end{array}$

obeys the scaling relation

$\begin{array}{ll}{\delta }_{\text{pred}}=\frac{{\tau }_t-1}{\tau -1}.& (6)\end{array}$

We estimated δ in two independent ways, as δ_pred and as δ_fit, i.e., the slope of the least square fit of the average sizes given their durations. In principle, if these two estimates are compatible, then the system is compatible with criticality. Proving such relation is however challenging. First, it is sensitive to the fitting methods of the distributions of avalanche sizes and lifetimes, as it has been recently shown (Destexhe and Touboul, 2021). Second, in the case of LFPs, the range of avalanche lifetimes typically extends over one order of magnitude only, which undermines reliability of power law fitting.

A complementary approach is to look for criticality by perturbing the system, that is by shifting activity to a subcritical or supercritical regime and then measure the distance from a critical state (Shew et al., 2009; Meisel, 2020). Thus, when measuring neural activity across a cortical barrel column of the rat brain, we also provided sensory stimuli consisting of impulsive deflections of the corresponding whisker and therefore representing strong, well defined, and accurately reproducible perturbations.

Neural avalanches sizes and durations were fitted with a discrete power law through the maximum likelihood method (see section 2.3). As mentioned in the Materials and Methods section, we found non negligible values of the correlation time τ^* both for spikes/MUAs and LFPs, suggesting that events were not independent. Therefore, to test if avalanches were power law distributed, we corrected for dependencies between avalanches by subsampling the data of sizes and durations as described in section 2.3. Further results are reported in Supplementary Tables 1–8.

3.2. Avalanches in LFP Data

We first analyzed avalanches in LFP recordings and found similar results across five animals (Table 1). First of all, we focused on spontaneous (i.e., resting, non stimulated) activity. We estimated τ^* and subsampled the sizes and the lifetimes accordingly to ensure that individual observations were uncorrelated (thus not violating the assumption of the maximum likelihood method). Following this correction, avalanches resulted power law distributed across the five rats, except for one case for the avalanches durations (see Figure 3 and Table 1).

FIGURE 3

Figure 3. Exponents obtained after subsampling in Local Field Potential data (τ is the exponent for the sizes, τ_t for the durations, δ_fit is the fitted exponent from $\text{S}(\text{T})~T^{{\delta }_{\text{fit}}}$). The exponents predicted by the critical branching process are shown in the plot as gray dashed lines. The line δ = 1.28 is also plotted, which is reported in Fontenele et al. (2019) as a universal exponent found in many different experiments. It is noteworthy that, as seen in the plot on the right, the exponents δ_fit are always far from the hallmark of the branching process, except for one case after stimulation, where the presence of large bumps in the size distribution makes the slope of the line in the (〈S(T)〉, T) plane steeper.

Then, we analyzed stimulus-evoked responses over 2 s after the stimulus and confirmed power law scaling (accepted by statistical tests). However, at a more careful look, the size distribution was altered by the presence of a bump and similar to the heap observed in previous work (di Santo et al., 2018) (Figure 4). Clearly, the bump derived from an excess of large size avalanches (i.e., involving a large number of microelectrodes). As suggested by Figure 2, these avalanches were not only large, but also composed by highly synchronous events across cortical layers. In fact the bump vanished in the duration distribution, confirming that the size increase was not accompanied by a corresponding increase of lifetime and that events remained concentrated within a few time bins.

FIGURE 4

Figure 4. LFP avalanches. Representative probability densities of avalanches sizes in one rat after the sensory stimulus (A), or in resting state conditions (B). A bump deviating from the power law can be recognized in the post-stimulus distribution in correspondence of large avalanche sizes, pointed to by an arrow. Distributions of the avalanches' durations do not display clear alterations instead (C,D). Durations are expressed in ms by multiplying the number of bins by 〈IEI〉. In the randomized dataset the exponential distribution provides a better fit in all cases (dashed gray lines). The shuffling procedure consists in randomizing the occurrence times of the events of each channel, so that the events rate of each channel is preserved. (E,F) The crackling noise relation is verified in both cases within the experimental errors (i.e., δ_pred is compatible with δ_fit). Once more, notice in (E) the presence of a bump (indicated by an arrow) in the post-stimulus regime.

In general, we found that the crackling noise relation was verified both at resting and during the post-stimulus (Figure 4 and see Supplementary Figure 5 for the avalanches results on all the rats). However, after the stimulus, a localized deviation from the expected trend was observed in correspondence of the bump found in the size distribution, and therefore also caused by large and synchronized waves of activity unleashed by whisker deflection. As anticipated previously, avalanches durations, and consequently also 〈S〉(T), extended over about two orders of magnitude. In fact, since LFPs are average signals that integrate over space and time single neuron events, it is expected that the maximum size of the avalanches extracted from LFPs is of the order of the array size (see Supplementary Figure 3 for an analysis of finite size effects in LFPs avalanches) (Beggs and Plenz, 2003). Specifically, a cutoff in P(S) around a value N_C ≈ N_E, with N_E the number of the recording electrodes, is commonly observed experimentally, meaning that during an avalanche each electrode is typically activated just once. This cutoff in P(S) implies that 〈S〉(T) < N_C, and, from 〈S〉(T) ~ T^δ, it also implies $T<N_C^{\frac{1}{\delta }}$ (Neto et al., 2020). Thus, if δ > 1, the cut-off in P(S) causes a much earlier cut-off in both P(T) and 〈S〉(T).

3.3. Avalanches in Spikes and MUAs Data

Spikes and MUAs avalanches (largely dominated by MUAs as reported in the Methods) were analyzed for four rats. Albeit similar to LFPs, the results unveiled important differences. First of all, MUAs were much less affected by the finite size of the recording array. In fact, despite the low number of microelectrodes (twenty-seven sites spanning across the cortex in the vertical direction), avalanches could be observed with size greater than 10² events. The reason is that, in MUAs, the same avalanche can reach an electrode repeatedly and in quick succession, contrary to LFPs where single neuron contributions are integrated over space and time within the brain tissue and coalesce to generate the recorded signal (Neto et al., 2020). As a consequence, MUAs provided a much better estimate of the avalanche duration distributions.

Similarly to LFPs, we found a power law distribution of avalanches both in terms of size and duration, except for one rat deviating from this general trend with respect to duration alone (Figure 6 and Table 2). Moreover, we confirmed the emergence of the bump of activity in the post-stimulus size distribution caused by the abundant number of large-sized avalanches (Figure 5). The exponents τ and τ_t were greater than the ones found in LFPs (Table 2), but also greater than the ones predicted for a critical branching process (Figure 6). Interestingly, the exponent δ was consistently close to the value δ ≈ 1.28 which was found in Fontenele et al. (2019) to be universal, i.e., to hold across different experimental conditions, from cultured slices to freely moving or anesthetized mammals.

FIGURE 5

Figure 5. Avalanches in MUAs data. Avalanches sizes probability density in one rat after a stimulus (A), and during resting state (B). Notice the bump in the post-stimulus distribution. (C,D) The same, for avalanche durations. Durations are expressed in ms by multiplying the number of bins times 〈IEI〉. In the randomized datasets the exponential distribution provides a better fit in all cases. The shuffling procedure consists in randomizing the occurrence times of the events of each channel, so that the events rate of each channel is preserved. (E,F) The crackling noise relation is verified in both cases within the experimental errors. Once more, notice in (E) the presence of a bump in the post-stimulus regime.

FIGURE 6

Figure 6. Exponents obtained after subsampling in multi-unit activities data (τ exponent for the sizes, τ_t for the durations, δ_fit the fitted exponents of the crackling-noise relation). The exponents predicted by the critical branching process are also plotted. The line δ_fit = 1.28 is also plotted, which is found in Fontenele et al. (2019), and which turns out to be remarkably close to our results.

Finally, as avalanche size and duration were less affected by finite size effects, MUAs allowed us to reliably test the crackling noise relation, which was verified in all dataset (Figure 5, see also Supplementary Figure 6). In particular, as MUAs data are believed to be more robust in this respect (Friedman et al., 2012), we concluded that post-stimulus avalanches statistics is also compatible with the results recently reported by Fontenele and collaborators (Fontenele et al., 2019). Nevertheless, also in this case, the waves of large and synchronized activity triggered by whisker deflection generated a local bump deviating from the crackling noise relation.

4. Discussion and Conclusions

In this work we investigated criticality in the rat barrel cortex, which offers advantages over other sensory systems. First, there is a clear and well characterized somatotopic representation of the whiskers in this primary sensory cortex, where single whiskers are mapped to single cortical columns with a one-to-one correspondence. This differentiates the barrel cortex from, e.g., the primary visual cortex, where functional columns are not as clearly segregated by septa in a cylindrical structure perpendicular to the cortex surface as in the barrel field, which facilitates the recording of single column activity by a neural probe (Horton and Adams, 2005). Second, whisker receptors, which are transducers placed around the follicle, directly activate primary sensory neurons, that therefore encode whisker deflections without interposed processing. In other circumstances, such as in the visual system, the transduced stimulus is subject to extensive processing already at the periphery (e.g., by the retina network) which makes it difficult to disentangle the dynamics and the contributions of cortical and pre-cortical networks in response to sensory inputs (Diamond et al., 2008; Feldmeyer et al., 2013). Third, single-whisker deflection can be controlled with high accuracy through closed-loop piezoelectric systems enabling a tight experimental control over repeated sensory stimuli.

As for other sensory cortical areas, it can be hypothesized that the barrel cortex takes advantage of criticality to efficiently map tactile stimuli. This holds also for the single barrel column, which faces the severe challenge to represent the parameters related to deflection of its corresponding whisker (e.g., amplitude, direction and velocity of displacement), as transduced by the follicle receptors, in an efficient, noise-tolerant manner and in real-time. Along a similar line of reasoning, previous work has addressed the barrel cortex (Gautam et al., 2015) using an air puff on the snout for stimulation and with avalanches analysis that was tailored to UP/DOWN states dominating activity under urethane anesthesia. According to a simplified general model of the processing in the barrel column, it is believed that layer IV acts as main input stage of sensory information propagating from the thalamus, whereas layer V is the main output. However, the few thousands neurons composing a single barrel form complex microcircuits of excitatory and inhibitory connections across layers. These microcircuits generate a rich dynamics that has been shown, both in-vitro and in-vivo, to include avalanches or synchronization states, such as UP and DOWN states and oscillations, and whose significance in terms of tactile information mapping and processing is far from being understood (Feldmeyer, 2012; Petersen, 2019).

We run our experiments under tiletamine anesthesia by precisely controlling single whisker deflection by a closed-loop controlled piezoelectric actuator. With tiletamine, that contrary to urethane is also used in clinics, spontaneous activity remains rich, with a mix of oscillations and UP and DOWN states, and with the activity evoked by whisker stimulation that shows clear similarities to the one in the awake animal (Rojas et al., 2006). In the context of our study, another advantage of the anesthetized animal with respect to the awake condition was that, by reducing cortico-cortical communication, the anesthetic was insulating the somatosensory barrel region from external “contaminating” waves of activity propagating from other brain areas (Aronoff et al., 2010; Voss et al., 2019; Sorrenti et al., 2021).

We measured neural activity across the six cortical layers of a single barrel both in the domain of spikes/MUAs (although MUAs were largely overwhelming) and LFPs, thus extending the avalanches analysis over a wide frequency range and covering both single neuron and population dynamics. Recently, some first attempts were made to characterize criticality in neural systems in a broader sense (Agrawal et al., 2019; Meshulam et al., 2019; Nicoletti et al., 2020) but, despite some well-known limitations (Beggs and Timme, 2012), as of now scale-free avalanches are still the most employed approach to test the critical brain hypothesis. For LFPs, we used a high-density 2D array of microelectrodes developed at the purpose to monitor at high spatial resolution the electrical potential within a planar section of the barrel column (Thewes et al., 2016) whereas for spikes and MUAs we used a 32 channels linear array probe.

On the one hand, and contrary to previous works based on imaging methods (Scott et al., 2014; Bellay et al., 2015), we found that in the anesthetized animal avalanche sizes and lifetimes followed power laws during basal activity, in agreement with findings on LFPs and spikes (Hahn et al., 2010; Ribeiro et al., 2010). Notably, we observe inter-rats variability in the power-law exponents (Figures 3, 6), but this phenomenon is not unexpected as it was observed in previous experiments (Shew et al., 2015) and a theoretical explanation has also been recently proposed in terms of quasi-criticality (Fosque et al., 2021). In particular, in Fontenele et al. (2019) a statistical meta-analysis of a large number of experiments suggests that, even if the avalanche exponents vary across individuals and across experimental settings, they all lie along a scaling line consistent with the crackling noise relation and an exponent δ ≈ 1.28. We do find that the exponent δ consistently converge to the same value (Figure 6). Intriguingly, and in contrast to the classical view of avalanches seen as generated by a quiescent-to-active phase transition—even though a different scenario involving a disorder-induced transition has been recently proposed as well (Ponce-Alvarez et al., 2018)—the analysis of Fontenele et al. (2019) suggests that the critical transition occurs at the edge of synchronization (di Santo et al., 2018) and this non-trivial exponent might emerge in different contexts (Dalla Porta and Copelli, 2019; Neto et al., 2020; Carvalho et al., 2021). In this context, we also note that if we assume as a generative mechanism a model with absorbing states, such as the branching process, it would imply no correlation between avalanches (Dalla Porta and Copelli, 2019). However, in our data, we found non-negligible values of τ^* hinting at the existence of temporal correlations between subsequent avalanches. These latter results are in agreement with previous works in which long range temporal correlations emerged through a detrended fluctuations analysis (DFA) of the signals (Linkenkaer-Hansen et al., 2001; Ribeiro et al., 2010; Hardstone et al., 2012). Although the nature of such correlations remains unknown, we speculate that they may derive from a variety of neuronal and circuital mechanisms including, e.g., a neuromodulation of the barrel column, acting at time scales that are large compared to avalanches lifetimes.

After the stimulus, power laws displayed localized bump-like alterations. The associated Kuramoto order parameters (see Supplementary Table 9) point at the emergence of a high, but brief, synchronization across layers. A control test of the avalanches' distributions without the first 200 ms after the stimulus shows indeed that the bumps in the sizes' distributions disappear both in LFPs and MUAs (see Supplementary Figure 7). Importantly, bumps are not just a spurious effect of the bin chosen, as we adopted a similar procedure to the one proposed by Yu et al. (2017) by evaluating 〈IEI〉 in the periods after stimuli separately from the periods of resting state. However, at the same time, one should note that bumps are also influenced by the finite size of the system (de Arcangelis, 2012). Indeed, being bounded by the system size, an excess of large avalanches can create a peak at a characteristic scale in the distributions. This is precisely the case of LFP avalanches, that are known to be constrained by the size of the array. In the case of MUAs, instead, finite size constraints are more relaxed, bumps are genuine although less evident in some MUAs distributions (see Supplementary Figure 6). Hence, in summary, the presence of bumps marks the occurrence of avalanches with large sizes that appear with a higher probability than what predicted by the power law scaling and, in the case of LFPs, these sizes are of the order of the array size.

This seemingly simple observation has quite profound implications. These avalanches, characterized by a highly synchronous activity across microelectrodes (see Figure 2), are related to a strong response to the stimulus that is thought to emerge prominently in layer IV first, and then to quickly spread among the other layers eliciting a global activity in the barrel. This suggests that synchronization waves together with oscillations play a fundamental role in shaping the neural activity during the mapping of the tactile stimulus. Indeed, the LFPs evoked response is also characterized by oscillations at a frequency close to 6 Hz (Supplementary Figure 14), whose power is particularly high in the central layers. These oscillations are not present in the spontaneous activity phase. Let us note that, most likely, these small-amplitude oscillations in the LFP mainly reflect subthreshold synaptic potentials as the 6 Hz component is not visible in the spikes domain (see Supplementary Figure 14). Moreover, being below the signal detection threshold for LFPs, these oscillations do not overall contribute to LFP avalanches distributions. In recent years, some modeling efforts taking into account synchronizations and oscillations phenomena were made (Poil et al., 2012; di Santo et al., 2018; Fontenele et al., 2019) but a comprehensive model of this kind remains elusive and further experimental studies will be needed to address this crucial point.

In perspective, in addition to further investigating the relation between resting state and stimuli with respect to brain criticality, an interesting route could be the investigation of other markers of criticality beyond neuronal avalanches, such as scale-free correlations (Mariani et al., 2021) or other patterns of dynamic functional connectivity, e.g., obtained by computing the instantaneous phase difference between signals at different locations (Cabral et al., 2017). On the other hand, it is likely that sound models of the possible self-organizing mechanisms of the cortex will have to take into account its oscillatory behavior, whose signature in the present work emerges through both the bumps in the avalanche size distribution and the presence of a characteristic peak in the power spectrum in the periods after stimuli. In perspective, avalanches and oscillations may represent two intertwined and functionally relevant faces of the cortical brain dynamics, that will have to be considered together in future theoretical and experimental investigations on the brain's sensory coding.

Data Availability Statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Ethics Statement

The animal study was reviewed and approved by Italian Ministry of Health.

Author Contributions

SV and SS designed the study. SV, SS, and MB supervised the research. MM and AL conducted the experiments. RO designed the whisker control system. BM analyzed the data. BM, GN, MB, MM, SS, and SV interpreted the results and wrote the manuscript. BM and GN prepared the figures. All authors contributed to the article and approved the submitted version.

Funding

SS acknowledges DFA and UNIPD for SUWE_BIRD2020_01 grant, and INFN for LINCOLN grant. This work was supported by the grant SYNCH to SV (European Commission, FET Proactive, GA N. 824162).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We thank A. Maritan for fruitful discussion.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fnsys.2021.709677/full#supplementary-material

References

Agrawal, V., Chakraborty, S., Knpfel, T., and Shew, W. L. (2019). Scale-change symmetry in the rules governing neural systems. iScience 12, 121–131. doi: 10.1016/j.isci.2019.01.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Aronoff, R., Matyas, F., Mateo, C., Ciron, C., Schneider, B., and Petersen, C. C. (2010). Long-range connectivity of mouse primary somatosensory barrel cortex. Eur. J. Neurosci. 31, 2221–2233. doi: 10.1111/j.1460-9568.2010.07264.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Beggs, J. M., and Plenz, D. (2003). Neuronal avalanches in neocortical circuits. J. Neurosci. 23, 11167–11177. doi: 10.1523/JNEUROSCI.23-35-11167.2003

CrossRef Full Text | Google Scholar

Beggs, J. M., and Timme, N. (2012). Being critical of criticality in the brain. Front. Physiol. 3:163. doi: 10.3389/fphys.2012.00163

PubMed Abstract | CrossRef Full Text | Google Scholar

Bellay, T., Klaus, A., Seshadri, S., and Plenz, D. (2015). Irregular spiking of pyramidal neurons organizes as scale-invariant neuronal avalanches in the awake state. Elife 4:e07224. doi: 10.7554/eLife.07224

PubMed Abstract | CrossRef Full Text | Google Scholar

Bowen, Z., Winkowski, D. E., Seshadri, S., Plenz, D., and Kanold, P. O. (2019). Neuronal avalanches in input and associative layers of auditory cortex. Front. Syst. Neurosci. 13:45. doi: 10.3389/fnsys.2019.00045

PubMed Abstract | CrossRef Full Text | Google Scholar

Buzsaki, G., Anastassiou, C., and Koch, C. (2012). The origin of extracellular fields and currents–EEG, ECoG, LFP and spikes. Nat. Rev. Neurosci. 13, 407–420. doi: 10.1038/nrn3241

PubMed Abstract | CrossRef Full Text | Google Scholar

Cabral, J., Vidaurre, D., Marques, P., Magalh aes, R., Silva Moreira, P., Miguel Soares, J., et al. (2017). Cognitive performance in healthy older adults relates to spontaneous switching between states of functional connectivity during rest. Sci. Rep. 7:5135. doi: 10.1038/s41598-017-05425-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Carvalho, T. T. A., Fontenele, A. J., Girardi-Schappo, M., Feliciano, T., Aguiar, L. A. A., Silva, T. P. L., et al. (2021). Subsampled directed-percolation models explain scaling relations experimentally observed in the brain. Front. Neural Circ. 14, 83. doi: 10.3389/fncir.2020.576727

PubMed Abstract | CrossRef Full Text | Google Scholar

Clauset, A., Shalizi, C. R., and Newman, M. E. J. (2009). Power-law distributions in empirical data. SIAM Rev. 51, 661–703. doi: 10.1137/070710111

CrossRef Full Text | Google Scholar

Clawson, W. P., Wright, N. C., Wessel, R., and Shew, W. L. (2017). Adaptation towards scale-free dynamics improves cortical stimulus discrimination at the cost of reduced detection. PLoS Comput. Biol. 13:e1005574. doi: 10.1371/journal.pcbi.1005574

PubMed Abstract | CrossRef Full Text | Google Scholar

Dalla Porta, L., and Copelli, M. (2019). Modeling neuronal avalanches and long-range temporal correlations at the emergence of collective oscillations: continuously varying exponents mimic MEEG results. PLoS Comput. Biol. 15:e1006924. doi: 10.1371/journal.pcbi.1006924

PubMed Abstract | CrossRef Full Text | Google Scholar

de Arcangelis, L. (2012). Are dragon-king neuronal avalanches dungeons for self-organized brain activity? Eur. Phys. J. Special Top. 205, 243–257. doi: 10.1140/epjst/e2012-01574-6

CrossRef Full Text | Google Scholar

Dehghani, N., Hatsopoulos, N. G., Haga, Z. G., Parker, R., Greger, B., Halgren, E., et al. (2012). Avalanche analysis from multielectrode ensemble recordings in cat, monkey, and human cerebral cortex during wakefulness and sleep. Front. Physiol. 3:302. doi: 10.3389/fphys.2012.00302

PubMed Abstract | CrossRef Full Text | Google Scholar

Destexhe, A., and Touboul, J. D. (2021). Is there sufficient evidence for criticality in cortical systems? eNeuro 8:ENEURO.0551-20.2021. doi: 10.1523/ENEURO.0551-20.2021

PubMed Abstract | CrossRef Full Text | Google Scholar

di Santo, S., Villegas, P., Burioni, R., and Mu noz, M. A. (2017). Simple unified view of branching process statistics: Random walks in balanced logarithmic potentials. Phys. Rev. E 95, 032115. doi: 10.1103/PhysRevE.95.032115

PubMed Abstract | CrossRef Full Text | Google Scholar

di Santo, S., Villegas, P., Burioni, R., and Mu noz, M. A. (2018). Landau–Ginzburg theory of cortex dynamics: Scale-free avalanches emerge at the edge of synchronization. Proc. Natl. Acad. Sci. U.S.A. 115, E1356–E1365. doi: 10.1073/pnas.1712989115

PubMed Abstract | CrossRef Full Text | Google Scholar

Diamond, M. E., von Heimendahl, M., Knutsen, P. M., Kleinfeld, D., and Ahissar, E. (2008). “where” and “what” in the whisker sensorimotor system. Nat. Rev. Neurosci. 9, 601–612. doi: 10.1038/nrn2411

PubMed Abstract | CrossRef Full Text | Google Scholar

Einevoll, G., Kayser, C., Logothetis, N., and Panzeri, S. (2013). Modelling and analysis of local field potentials for studying the function of cortical circuits. Nat. Rev. Neurosci. 14, 770–785. doi: 10.1038/nrn3599

PubMed Abstract | CrossRef Full Text | Google Scholar

Feldmeyer, D. (2012). Excitatory neuronal connectivity in the barrel cortex. Front. Neuroanat. 6:24. doi: 10.3389/fnana.2012.00024

PubMed Abstract | CrossRef Full Text | Google Scholar

Feldmeyer, D., Brecht, M., Helmchen, F., Petersen, C. C., Poulet, J. F., Staiger, J. F., et al. (2013). Barrel cortex function. Progr. Neurobiol. 103, 3–27. doi: 10.1016/j.pneurobio.2012.11.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Fontenele, A. J., de Vasconcelos, N. A. P., Feliciano, T., Aguiar, L. A. A., Soares-Cunha, C., Coimbra, B., et al. (2019). Criticality between cortical states. Phys. Rev. Lett. 122, 208101. doi: 10.1103/PhysRevLett.122.208101

PubMed Abstract | CrossRef Full Text | Google Scholar

Fosque, L. J., Williams-Garc-a, R. V., Beggs, J. M., and Ortiz, G. (2021). Evidence for quasicritical brain dynamics. Phys. Rev. Lett. 126:098101. doi: 10.1103/PhysRevLett.126.098101

PubMed Abstract | CrossRef Full Text | Google Scholar

Friedman, N., Ito, S., Brinkman, B. A. W., Shimono, M., DeVille, R. E. L., Dahmen, K. A., et al. (2012). Universal critical dynamics in high resolution neuronal avalanche data. Phys. Rev. Lett. 108:208102. doi: 10.1103/PhysRevLett.108.208102

PubMed Abstract | CrossRef Full Text | Google Scholar

Gautam, S. H., Hoang, T. T., McClanahan, K., Grady, S. K., and Shew, W. L. (2015). Maximizing sensory dynamic range by tuning the cortical state to criticality. PLoS Comput. Biol. 11:e1004576. doi: 10.1371/journal.pcbi.1004576

PubMed Abstract | CrossRef Full Text | Google Scholar

Gerlach, M., and Altmann, E. G. (2019). Testing statistical laws in complex systems. Phys. Rev. Lett. 122:168301. doi: 10.1103/PhysRevLett.122.168301

PubMed Abstract | CrossRef Full Text | Google Scholar

Gireesh, E. D., and Plenz, D. (2008). Neuronal avalanches organize as nested theta- and beta/gamma-oscillations during development of cortical layer 2/3. Proc. Natl. Acad. Sci. U.S.A. 105, 7576–7581. doi: 10.1073/pnas.0800537105

PubMed Abstract | CrossRef Full Text | Google Scholar

Hahn, G., Petermann, T., Havenith, M. N., Yu, S., Singer, W., Plenz, D., et al. (2010). Neuronal avalanches in spontaneous activity in vivo. J. Neurophysiol. 104, 3312–3322. doi: 10.1152/jn.00953.2009

PubMed Abstract | CrossRef Full Text | Google Scholar

Hardstone, R., Poil, S.-S., Schiavone, G., Jansen, R., Nikulin, V., Mansvelder, H., et al. (2012). Detrended fluctuation analysis: a scale-free view on neuronal oscillations. Front. Physiol. 3:450. doi: 10.3389/fphys.2012.00450

PubMed Abstract | CrossRef Full Text | Google Scholar

Hidalgo, J., Grilli, J., Suweis, S., Mu noz, M. A., Banavar, J. R., and Maritan, A. (2014). Information-based fitness and the emergence of criticality in living systems. Proc. Natl. Acad. Sci. U.S.A. 111, 10095–10100. doi: 10.1073/pnas.1319166111

PubMed Abstract | CrossRef Full Text | Google Scholar

Horton, J. C., and Adams, D. L. (2005). The cortical column: a structure without a function. Philos. Trans. R. Soc. B Biol. Sci. 360, 837–862. doi: 10.1098/rstb.2005.1623

PubMed Abstract | CrossRef Full Text | Google Scholar

Kinouchi, O., and Copelli, M. (2006). Optimal dynamical range of excitable networks at criticality. Nat. Phys. 2, 348–351. doi: 10.1038/nphys289

CrossRef Full Text | Google Scholar

Linkenkaer-Hansen, K., Nikouline, V. V., Palva, J. M., and Ilmoniemi, R. J. (2001). Long-range temporal correlations and scaling behavior in human brain oscillations. J. Neurosci. 21, 1370–1377. doi: 10.1523/JNEUROSCI.21-04-01370.2001

PubMed Abstract | CrossRef Full Text | Google Scholar

Mariani, B., Nicoletti, G., Bisio, M., Maschietto, M., Vassanelli, S., and Suweis, S. (2021). On the critical signatures of neural activity. arXiv:2105.05070.

Google Scholar

Mazzoni, A., Broccard, F. D., Garcia-Perez, E., Bonifazi, P., Ruaro, M., and Torre, V. (2007). On the dynamics of the spontaneous activity in neuronal networks. PLoS ONE 2:e439. doi: 10.1371/journal.pone.0000439

PubMed Abstract | CrossRef Full Text | Google Scholar

Meisel, C. (2020). Antiepileptic drugs induce subcritical dynamics in human cortical networks. Proc. Natl. Acad. Sci. U.S.A. 117, 11118–11125. doi: 10.1073/pnas.1911461117

PubMed Abstract | CrossRef Full Text | Google Scholar

Meshulam, L., Gauthier, J. L., Brody, C. D., Tank, D. W., and Bialek, W. (2019). Coarse graining, fixed points, and scaling in a large population of neurons. Phys. Rev. Lett. 123, 178103. doi: 10.1103/PhysRevLett.123.178103

PubMed Abstract | CrossRef Full Text | Google Scholar

Mu noz, M. A. (2018). Colloquium: criticality and dynamical scaling in living systems. Rev. Mod. Phys. 90:031001 doi: 10.1103/RevModPhys.90.031001

CrossRef Full Text | Google Scholar

Neto, J. P., Spitzner, F. P., and Priesemann, V. (2020). A unified picture of neuronal avalanches arises from the understanding of sampling effects. arXiv:1910.09984v09982. doi: 10.1101/759613

CrossRef Full Text | Google Scholar

Nicoletti, G., Suweis, S., and Maritan, A. (2020). Scaling and criticality in a phenomenological renormalization group. Phys. Rev. Res. 2:023144. doi: 10.1103/PhysRevResearch.2.023144

CrossRef Full Text | Google Scholar

Pasquale, V., Massobrio, P., Bologna, L., Chiappalone, M., and Martinoia, S. (2008). Self-organization and neuronal avalanches in networks of dissociated cortical neurons. Neuroscience 153, 1354–1369. doi: 10.1016/j.neuroscience.2008.03.050

PubMed Abstract | CrossRef Full Text | Google Scholar

Petermann, T., Thiagarajan, T. C., Lebedev, M. A., Nicolelis, M. A. L., Chialvo, D. R., and Plenz, D. (2009). Spontaneous cortical activity in awake monkeys composed of neuronal avalanches. Proc. Natl. Acad. Sci. U.S.A. 106, 15921–15926. doi: 10.1073/pnas.0904089106

PubMed Abstract | CrossRef Full Text | Google Scholar

Petersen, C. (2019). Sensorimotor processing in the rodent barrel cortex. Nat. Rev. Neurosci. 20, 533–546. doi: 10.1038/s41583-019-0200-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Poil, S.-S., Hardstone, R., Mansvelder, H. D., and Linkenkaer-Hansen, K. (2012). Critical-state dynamics of avalanches and oscillations jointly emerge from balanced excitation/inhibition in neuronal networks. J. Neurosci. 32, 9817–9823. doi: 10.1523/JNEUROSCI.5990-11.2012

PubMed Abstract | CrossRef Full Text | Google Scholar

Ponce-Alvarez, A., Jouary, A., Privat, M., Deco, G., and Sumbre, G. (2018). Whole-brain neuronal activity displays crackling noise dynamics. Neuron 100, 1446.e6–1459.e6. doi: 10.1016/j.neuron.2018.10.045

PubMed Abstract | CrossRef Full Text | Google Scholar

Priesemann, V., Wibral, M., Valderrama, M., Prpper, R., Le Van Quyen, M., Geisel, T., et al. (2014). Spike avalanches in vivo suggest a driven, slightly subcritical brain state. Front. Syst. Neurosci. 8:108. doi: 10.3389/fnsys.2014.00108

PubMed Abstract | CrossRef Full Text | Google Scholar

Quiroga, R. Q., Nadasdy, Z., and Ben-Shaul, Y. (2004). Unsupervised spike detection and sorting with wavelets and superparamagnetic clustering. Neural Comput. 16, 1661–1687. doi: 10.1162/089976604774201631

PubMed Abstract | CrossRef Full Text | Google Scholar

Ren, J., Aika, Y., Heizmann, C., and Kosaka, T. (1992). Quantitative analysis of neurons and glial cells in the rat somatosensory cortex, with special reference to gabaergic neurons and parvalbumin-containing neurons. Exp. Brain Res. 92, 1–14. doi: 10.1007/BF00230378

PubMed Abstract | CrossRef Full Text | Google Scholar

Ribeiro, T. L., Copelli, M., Caixeta, F., Belchior, H., Chialvo, D. R., Nicolelis, M. A. L., et al. (2010). Spike avalanches exhibit universal dynamics across the sleep-wake cycle. PLoS ONE 5:e14129. doi: 10.1371/journal.pone.0014129

PubMed Abstract | CrossRef Full Text | Google Scholar

Rojas, M. J., Navas, J. A., and Rector, D. M. (2006). Evoked response potential markers for anesthetic and behavioral states. Am. J. Physiol. Regul. Integr. Comp. Physiol. 291, R189–R196. doi: 10.1152/ajpregu.00409.2005

PubMed Abstract | CrossRef Full Text | Google Scholar

Scarpetta, S., and de Candia, A. (2014). Alternation of up and down states at a dynamical phase-transition of a neural network with spatiotemporal attractors. Front. Syst. Neurosci. 8:88. doi: 10.3389/fnsys.2014.00088

PubMed Abstract | CrossRef Full Text | Google Scholar

Schroder, S., Cecchetto, C., Keil, S., Mahmud, M., Brose, E., Dogan, Ö., et al. (2015). “Cmos-compatible purely capacitive interfaces for high-density in-vivo recording from neural tissue,” in 2015 IEEE Biomedical Circuits and Systems Conference (BioCAS) (Atlanta, GA; IEEE), 1–4.

Google Scholar

Scott, G., Fagerholm, E. D., Mutoh, H., Leech, R., Sharp, D. J., Shew, W. L., et al. (2014). Voltage imaging of waking mouse cortex reveals emergence of critical neuronal dynamics. J. Neurosci. 34, 16611–16620. doi: 10.1523/JNEUROSCI.3474-14.2014

PubMed Abstract | CrossRef Full Text | Google Scholar

Sethna, J., Dahmen, K. A., and Myers, C. (2001). Crackling noise. Nature 410, 242–250. doi: 10.1038/35065675

PubMed Abstract | CrossRef Full Text | Google Scholar

Shew, W. L., Clawson, W. P., Pobst, J., Karimipanah, Y., Wright, N. C., and Wessel, R. (2015). Adaptation to sensory input tunes visual cortex to criticality. Nat. Phys. 11, 659–663. doi: 10.1038/nphys3370

CrossRef Full Text | Google Scholar

Shew, W. L., and Plenz, D. (2013). The functional benefits of criticality in the cortex. Neuroscientist 19, 88–100. doi: 10.1177/1073858412445487

PubMed Abstract | CrossRef Full Text | Google Scholar

Shew, W. L., Yang, H., Petermann, T., Roy, R., and Plenz, D. (2009). Neuronal avalanches imply maximum dynamic range in cortical networks at criticality. J. Neurosci. 29, 15595–15600. doi: 10.1523/JNEUROSCI.3864-09.2009

PubMed Abstract | CrossRef Full Text | Google Scholar

Sorrenti, V., Cecchetto, C., Maschietto, M., Fortinguerra, S., Buriani, A., and Vassanelli, S. (2021). Understanding the effects of anesthesia on cortical electrophysiological recordings: a scoping review. Int. J. Mol. Sci. 22:1286. doi: 10.3390/ijms22031286

PubMed Abstract | CrossRef Full Text | Google Scholar

Swanson, L. (2003). Brain Maps: Structure of the Rat Brain, 3rd Edn. London, Academic Press.

Google Scholar

Thewes, R., Bertotti, G., Dodel, N., Keil, S., Schrder, S., Boven, K.-H., et al. (2016). “Neural tissue and brain interfacing cmos devices–an introduction to state-of-the-art, current and future challenges,” in 2016 IEEE International Symposium on Circuits and Systems (ISCAS) (Montreal, QC: IEEE), 1826–1829.

Google Scholar

Touboul, J., and Destexhe, A. (2010). Can power-law scaling and neuronal avalanches arise from stochastic dynamics? PLoS ONE 5::e8982. doi: 10.1371/journal.pone.0008982

PubMed Abstract | CrossRef Full Text | Google Scholar

Touboul, J., and Destexhe, A. (2017). Power-law statistics and universal scaling in the absence of criticality. Phys. Rev. E 95, 012413. doi: 10.1103/PhysRevE.95.012413

PubMed Abstract | CrossRef Full Text | Google Scholar

Villegas, P., di Santo, S., Burioni, R., and Mu noz, M. A. (2019). Time-series thresholding and the definition of avalanche size. Phys. Rev. E 100:012133. doi: 10.1103/PhysRevE.100.012133

PubMed Abstract | CrossRef Full Text | Google Scholar

Voss, L. J., Garc-a, P. S., Hentschke, H., and Banks, M. I. (2019). Understanding the effects of general anesthetics on cortical network activity using ex vivo preparations. Anesthesiology 130, 1049–1063. doi: 10.1097/ALN.0000000000002554

PubMed Abstract | CrossRef Full Text | Google Scholar

Wilting, J., and Priesemann, V. (2019). 25 years of criticality in neuroscience–established results, open controversies, novel concepts. Curr. Opin Neurobiol. 58, 105–111. doi: 10.1016/j.conb.2019.08.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, S., Ribeiro, T. L., Meisel, C., Chou, S., Mitz, A., Saunders, R., et al. (2017). Maintained avalanche dynamics during task-induced changes of neuronal activity in nonhuman primates. Elife 6:e27119. doi: 10.7554/eLife.27119

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, S., Yang, H., Nakahara, H., Santos, G. S., Nikolić, D., and Plenz, D. (2011). Higher-order interactions characterized in cortical activity. J. Neurosci. 31, 17514–17526. doi: 10.1523/JNEUROSCI.3127-11.2011

PubMed Abstract | CrossRef Full Text | Google Scholar