Distinct criticality of phase and amplitude dynamics in the resting brain

Converging research suggests that the resting brain operates at the cusp of dynamic instability signified by scale-free temporal correlations. We asked if the scaling properties of these correlations differ between amplitude and phase fluctuations, which may reflect different aspects of cortical functioning. Using source-reconstructed magneto-encephalographic signals, we found power-law scaling for the collective amplitude and for phase synchronization, both capturing whole-brain activity. The temporal changes of the amplitude comprise slow, persistent memory processes, whereas phase synchronization exhibits less temporally structured and more complex correlations, indicating a fast and flexible coding. This distinct temporal scaling supports the idea of different roles of amplitude and phase in cortical functioning.


Introduction
It has been proposed that the brain is in or near a critical state -its dynamics may be positioned at the border between spatiotemporal order and disorder, reminiscent of non-equilibrium phase transitions in thermodynamic systems [1,2,3,4]. The concept of brain criticality is attractive because critical systems display optimal performance on several characteristics such as information transfer [5,2], wide dynamic range [6,7], information capacity [8,9], and longterm stability [10,11]. Criticality relates closely to self-organization [12,10,11], which is considered crucial to cortical functioning [13,14,15,16,17].
A hallmark of critical behavior is the presence of power laws [18,4]. Power laws symbolize scale-free behavior, adopting the same form on all time scales: they are self-similar. Consider the case of a scale-free auto-correlation function AC. The corresponding power law obeys the form AC (s · τ ) = s 2H · AC (τ ), i.e.
if time τ is rescaled to s · τ , then the shape of AC is preserved and only rescaled by a factor s 2H . The scaling exponent H is referred to as the Hurst exponent [19] and qualifies the underlying correlation structure: H = 0.5 corresponds to an uncorrelated, random process whereas H > 0.5 indicates persistent, long-range correlations.
There is accumulating evidence for the presence of power laws in brain activity [20,21]. Neural spikes come in avalanches that display scale-free distributions [2]. Spectral distributions of encephalographic signals have 1/f -structures [22,23,24,25,26] and their auto-correlation structures also show power-law behavior [1,27]. In the spatial domain, scale-free distributions have been observed in functional as well as neuroanatomical connectivity patterns [28,29,30].
Previous work has focused only on spatially local measures of brain activity [22,1,2,27,20,21] or considered pairs of nodes in networks [31,32] rather than analyzing global brain activity. In complex systems, the global activity can be very informative about the generating dynamical structure [8], in particular when studying critical behavior [33]. Here, we adopted these concepts to investigate criticality in the brain. By using source-reconstructed magneto-encephalographic (MEG) signals we sought to disambiguate between the scaling characteristics of amplitude and phase fluctuations because they may resemble different aspects of cortical functioning.

MEG data & outcome variables
Magneto-encephalographic (MEG) signals of ten subjects were recorded and sampled at 1 kHz in eyes-closed resting state for approximately five minutes.
After down sampling to 250 Hz, signals were beamformed onto a ninety-node brain parcellation, yielding ninety time series y k (t) per subject. Data were previously published by Cabral and coworkers [34].
We used two collective variables to capture whole-brain activity per subject.
First, we defined the phase synchronization, and, second, the mean amplitude, We note that R (t, f ) is the modulo of the complex-valued Kuramoto order parameter [35]. We z-scored R (t, f ) and A (t, f ) to reduce between-subject variability such that we could assess subject-averaged behavior by means of a detrended fluctuation analysis (DFA, [36]), as described below.
To relate our study to the previously established results in RSNs, we also examined the amplitude dynamics in more detail (see also [37]). The expression of RSNs is mainly reflected in the low-frequency content of the amplitudes a k (t, f ), whose time scale is comparable to those of the blood-oxygenationlevel-dependent (BOLD) signal [38]. In order to study these slow amplitude dynamics, we evaluated the longer time scales of the a k (t, f ) dynamics by extracting its amplitude a  k (t, f ). Subsequently, we defined the collective variables associated with the amplitude dynamics: Further analysis was identical to that described above; see Figure 1 for illustration. (t,f)

Data analysis -DFA
DFA is considered robust against non-stationarity rendering it suitable for analyzing the temporal autocorrelation structure of encephalographic activity, in general, and global amplitude and phase synchronization, in particular. We employed a modified form of DFA including (Bayesian) model selection to verify the presence of power-law behavior [39].
In a nutshell, to quantify the autocorrelation structure of (the cumulative sum of) a signal Y (t), one divides it into non-overlapping segments Y i (t), with This definition yields a set of 'realizations' of fluctuations F i that, in the presence of a power law, scale like F i (n · τ ) = n α ·F i (τ ), which is equivalent to log (F i ) = α · log (n) + const. That is, in a log-log representation these fluctuations have a linear relationship with segment size. DFA seeks to identify the scaling exponent α that provides an estimate for the aforementioned Hurst exponent H.
Instead of computing the mean value of F i as in conventional DFA [36], we here determined the probability density function p n (F i ) for every segment length n (see Figure 2). This approach allows for quantifying the appropriateness of a model f θ (ñ) for fitting the fluctuation structureF = log (F i ) as a function ofñ = log (n) by means of the log-likelihood function ln (L) = n ln (p n (f θ )), where the tilde indicates a transformation to logarithmic coordinates. Here, the model f θ (ñ), parametrized by the set θ, may obey any arbitrary form including the linear one, which corresponds to a power law. We tested this linear relationship against a set of alternative models ( Table 1) using both the Bayesian information criterion and the Akaike information criterion (BIC and AIC c , respectively). The model resulting in the least value of the information criterion was selected as the proper model. Whenever this yielded the linear model f 1 θ , we considered power-law behavior to be present and identified the scaling exponent with its slope [39].
is the form a power law would adopt. The alternative models f 2 θ -f 7 θ constitute polynomials up to third order. With f 8 θ (x) and f 9 θ (x) we considered two models resembling a (un)stable linear stochastic dynamics.
Since we were interested in the subject-averaged scaling exponents, we determined p n for every subject individually and averaged over subjects to obtain p n . These averaged probability density functions were used both for model selection and to determine the scaling exponent α. The scaling range was given by [n min , n max ] = 1.875 · 10 2 , 1.875 · 10 4 [0.75, 75] seconds, i.e. it spanned two decades. In this range we used one hundred equally spaced window sizes (on a logarithmic scale).

Statistics -surrogate data
Bootstrapping served to establish statistical significance using three types of surrogate data. Two of them consisted of randomly permuting temporal order   For the first type of surrogates, we randomly permuted the order parameter time series for each subject. With this we evaluated our DFA and fitting procedure, since the surrogate time series lacked any temporal correlation structure and therefore should result in α = 0.5 [40]. With the second type of surrogate we evaluated whether the filtering procedure and Hilbert analysis could have biased the results. For this, we permuted all original time series y k (t).
By permuting y k (t), all temporal structure was destroyed and therefore this constituted a rather weak null. In the third type of surrogates we performed a random cycling of y k (t) by shifting the time indices of y k (t) for each k, but keeping their order intact. In this way we retained the original auto-correlation structure of φ  k (t, f ) but destroyed the cross-correlations. Subsequent analyses for all collective variables were identical to those for the original data.

Results
The presence of power-law scaling was evidenced by clear linear relationships of A (t, f ) and R (t, f ) fluctuations in log-log scale; see   In the alpha band the scaling exponents were 0.88 and 0.62, and in the beta band 0.80 and 0.57, for amplitude A and phase synchronization R, respectively.
That is, both A (t, f ) and R (t, f ) showed persistent behavior but amplitude had increased perseverance compared to phase (Figure 3).
We also found long-range temporal correlations in the variables A (a) (t, f ) and R (a) (t, f ) as shown in Figure 4 and confirmed by the AIC c and BIC values by preferring the linear model, except in case of the AIC c for A (a) (t, f ) in the alpha band (Tables 2-3). The difference in scaling behavior between A (a) (t, f ) and R (a) (t, f ) was similar to that between A (t, f ) and R (t, f ).   Testing against surrogate data confirmed the significance of these correlations. Figure 5   To further highlight the peculiar role of amplitude and phase, we finally contrasted our results with the scaling of fluctuations of the mean MEG activity.
For this, we applied our DFA to Y (t, f ) = 1 90 90 k=1 y k (t, f ), i.e. we considered not the amplitude and phase but the 'raw' MEG signals. This mean activity did not display long-range correlations but rather anti-persistent ones (α = 0.02, Figure 6).

Discussion
We report power-law scaling in both amplitude and phase of collective neural activity on long time scales, which is consistent with the hypothesis that the brain operates in a critical state. Operating in a critical state is not the only way a system can generate power-law scaling. Systems in subcritical states [32] or merely stochastic systems [41,42] may also display power laws. Biological systems display sub-and supercritical dynamics but they can be tuned into criticality [43,11,44]. We favor the interpretation of critical states, also because it is consistent with scale-free auto-correlation structures of single channel EEG activity [1,27], size and duration of neural avalanches [2,16,45] and, in the spatial domain, degree distributions of neuroanatomical and functional connectivity networks [28,29]. Such scaling laws in neuronal dynamics are also correlated with those found in behavior [46]. been reported using pair-wise synchronization measures like PLI and ∆ 2 (t, ∆t) by, e.g., Kitzbichler and coworkers [31,47] and Farmer [48]. In contrast, we evaluated scaling behavior on a global brain scale by using overall amplitude and phase synchrony as collective variables. Analyzing the scaling behavior of these variables enabled us to directly compare phase and amplitude behavior in brain activity.
What does the difference in power-law scaling of amplitudes and phases tell us about information processing in the brain? Despite the fact that amplitude and phase differ in their perseverance, both show long-range correlations over a scale of hundreds of seconds, suggesting that they reflect 'memory' of cortical states. A higher scaling exponent reflects a slower decay of auto-correlations and hence a more predictable signal with decreased entropy [49]. Therefore, we speculate that the difference in scaling exponents and their associated complexity reflect that amplitude and phase play different roles in information process- Volume conduction can be a confounder in analyzing encephalographic recordings [50]. Several methods to mitigate its effects have been proposed, generally relying on removing the instantaneous interactions that signify volume conduction [51,52,53]. In consequence these methods can only be applied in a pairwise fashion, such that they are not applicable when considering the variables R (·) (t, f ) and A (·) (t, f ). We note, however, that volume conduction does not significantly influence the auto-correlation structure of the signals under study.
This finds support by Shriki and coworkers [45] who showed that mere linear mixing cannot 'transform' uncorrelated activity to power-law scaling. In fact, if activity displays a power law, linear mixing does not alter this apart from slightly lowering the scaling exponent.
The occurrence of power laws is not only consistent with the 'criticality hypothesis' [3]. The macroscopic behavior of self-organizing processes canin general -be cast into a low-dimensional system when critical [54]; this in fact motivated looking at collective amplitude and phase synchrony. Several modeling studies support the seminal role of self-organization in neural dynamics and often highlight self-organized criticality. For example, synaptic plasticity under the influence of a simple learning rule leads to scale-free networks [13,14,17] and power-law distributions of avalanche dynamics [16]. In resting state, the mechanism of local feedback mediated inhibition increases model performance [55] and may also be interpreted as a form of self organization. Furthermore, the slow evolution of the order parameters in self-organizing systems is consistent with the time scale on which RSNs evolve [38,56,57,37].
While previous work has shown that RSNs fluctuate at slow (> 1sec) time scales [38,37], it has been recently shown that RSNs are also expressed in MEG activity on faster (< 1sec) time scales [58]. Alongside the results in this paper, this represents mounting evidence that RSNs are expressed across a range of time scales, i.e. they are time scale invariant.
In summary we have shown the presence of persistent long-range correlations in the evolution of global brain dynamics, i.e. the auto-correlation function obeys a power-law with scaling exponents exceeding those corresponding to random processes without memory. This adds further support to the hypothesis that the brain is in a (permanently) critical state. The here-reported scaling exponents clearly discriminate amplitude and phase dynamics, suggesting their differential role in whole-brain information processing.