## Abstract

The human brain spontaneously generates neural oscillations with a large variability in frequency, amplitude, duration, and recurrence. Little, however, is known about the long-term spatiotemporal structure of the complex patterns of ongoing activity. A central unresolved issue is whether fluctuations in oscillatory activity reflect a memory of the dynamics of the system for more than a few seconds.

We investigated the temporal correlations of network oscillations in the normal human brain at time scales ranging from a few seconds to several minutes. Ongoing activity during eyes-open and eyes-closed conditions was recorded with simultaneous magnetoencephalography and electroencephalography. Here we show that amplitude fluctuations of 10 and 20 Hz oscillations are correlated over thousands of oscillation cycles. Our analyses also indicated that these amplitude fluctuations obey power-law scaling behavior. The scaling exponents were highly invariant across subjects. We propose that the large variability, the long-range correlations, and the power-law scaling behavior of spontaneous oscillations find a unifying explanation within the theory of self-organized criticality, which offers a general mechanism for the emergence of correlations and complex dynamics in stochastic multiunit systems. The demonstrated scaling laws pose novel quantitative constraints on computational models of network oscillations. We argue that critical-state dynamics of spontaneous oscillations may lend neural networks capable of quick reorganization during processing demands.

- spontaneous oscillations
- large-scale dynamics
- temporal properties
- correlations
- scaling behavior
- self-organized criticality
- complexity

Oscillations at various frequencies are a prominent feature of the spontaneous electroencephalogram (EEG) (Berger, 1929; Connors and Amitai, 1997) and are believed to reflect functional states of the brain (Llinás, 1988; Steriade et al., 1993; Arieli et al., 1996; Herculano-Houzel et al., 1999; Tsodyks et al., 1999). These oscillations arise from correlated activity of a large number of neurons whose interactions are generally nonlinear (Steriade et al., 1990, 1993; Lopez da Silva, 1991). The intrinsic neural properties and intricate patterns of connectivity add further complexity to the behavior of neural systems (Llinás, 1988;Connors and Amitai, 1997; Destexhe et al., 1998). The mechanisms and dynamics of network oscillations have been widely studied with electrophysiological recordings (Destexhe et al., 1998, 1999), as well as with computational models (Destexhe et al., 1998; Stam et al., 1999). Neural oscillations *in vivo* exhibit large variability in both amplitude and frequency. The dynamic nature of these fluctuations, however, has remained unclear. Particularly for the human electroencephalogram, 8–13 Hz oscillations have attracted widespread interest in this context. However, the complexity of the EEG has rendered it impossible to reliably distinguish the waxing and waning of oscillations over epochs longer than 2–15 sec from that of filtered white noise (Paluš, 1996; Cerf et al., 1997; Stam et al., 1999). This suggests that the underlying neural populations are unlikely to obey entirely low-dimensional dynamics.

Recent studies have demonstrated that a large variety of complex processes, including forest fires (Malamud et al., 1998), earthquakes (Bak, 1997), financial markets (Mantegna and Stanley, 1995; Lux and Marchesi, 1999), heartbeats (Peng et al., 1995), and human coordination (Gilden et al., 1995; Chen et al., 1997), exhibit unexpected statistical similarities, most commonly power-law scaling behavior of a particular observable. Scaling behavior (or scale-free behavior) means that no characteristic scales dominate the dynamics of the underlying process. Scale-free behavior can be revealed with scaling analysis, which quantifies the fluctuations of a parameter as a function of the scale at which the parameter is evaluated. Scale-free behavior reflects a tendency of complex systems to develop correlations that decay more slowly and extend over larger distances in time and space than the mechanisms of the underlying process would suggest (Bassingthwaighte et al., 1994; Barabási and Stanley, 1995; Bak, 1997). The long-range correlations build up through local interactions until they extend throughout the entire system. After this stage, the dynamics of the system exhibit power-law scaling behavior, and the underlying process operates in a critical state, a phenomenon often termed self-organized criticality (SOC) (Bak et al., 1987, 1988). Unlike deterministic approaches aimed at finding low-dimensional chaos, the SOC framework allows for a high-dimensional character of the dynamics and for the presence of stochastic effects.

We have investigated whether noninvasively recorded spontaneous oscillations in the human brain show scaling behavior. Here we demonstrate the presence of long-range temporal correlations and power-law scaling behavior of oscillations at ∼10 and 20 Hz.

Parts of this work have been published previously in abstract form (Linkenkaer-Hansen et al., 2000).

## MATERIALS AND METHODS

*Recordings and experimental conditions.* Spontaneous brain activity from 10 normal subjects (aged 20–30 years, one female) was recorded simultaneously with magnetoencephalography (MEG) and EEG using a whole-scalp magnetometer array containing 122 planar gradiometers (Knuutila et al., 1993) and a 64-channel EEG cap (Virtanen et al., 1996). The study was approved by the Ethics Committee of the Department of Radiology of the Helsinki University Central Hospital. The subjects were seated in a magnetically shielded room and instructed to relax with eyes open or closed in two 20 min recording sessions. The MEG and EEG data were sampled at 300 Hz and the pass band of 0.3–90 Hz.

*Data analysis.* The amplitude of neural oscillations was estimated with wavelet filtering and subsequently evaluated for the presence of temporal correlations using the autocorrelation function (ACF) and detrended fluctuation analysis (DFA). As control for the neural origin of temporal correlations, we used an MEG reference recording and surrogate data.

*Wavelet filtering.* The signals were filtered with a Morlet wavelet; the modulus of the complex-valued outcome, *W*(*t*,*f*) , represents the amplitude of the signal at a time range centered at *t* and in a frequency band centered at *f* (Torrence and Compo, 1998). For each frequency band, we centered the wavelet at the peak frequency determined individually with amplitude spectra. The widths of the wavelet in the time and frequency domains are expressed as the attenuation by a factor of *e*^{2}and denoted *t _{e}
* and

*f*, respectively: and where

_{e}*m*is the Morlet parameter determining the compromise between time and frequency resolution (here,

*m*= 6). For a typical alpha oscillation at

*f*= 10 Hz, the signal is thus integrated by the wavelet for ∼262 msec in the time domain and 3.4 Hz in the frequency domain (i.e., the effective pass band is 8.3–11.7 Hz).

*Temporal correlations.* Temporal correlations of oscillations were quantified with the ACF and the DFA applied to the modulus of the wavelet filtered signals, i.e., to the amplitude envelope of the oscillations at a given frequency.

The autocorrelation function gives a measure of how a signal is correlated with itself at different time delays. When normalized, the autocorrelation attains its maximum value of one at zero time lag, decays toward zero with increasing time lag for (finite) correlated signals, and fluctuates around and close to zero at time lags free of correlations.

The detrended fluctuation analysis has been developed for quantifying correlation properties in nonstationary signals, e.g., in physiological time series, because long-range correlations—revealed by an ACF analysis—can arise also as an artifact of the “patchiness” of nonstationary data (Peng et al., 1994, 1995). In DFA, the modulus of the wavelet-transformed signal at center frequency *f* is first integrated to produce a vector *y* of the cumulative sum of the signal amplitude around its average value:
where *N* is the number of samples in the signal:
Equation 1The integrated signal is then divided into time windows of size τ. For each window, the least-squares fitted line (the local trend) is computed; the *y* coordinate of this line is denoted*y _{τ}
*(

*t*). The integrated signal,

*y*(

*t*), is detrended by subtracting the local trend,

*y*(

_{τ}*t*), in each window. The average root-mean-square fluctuation, 〈

*F*(τ)〉, of this integrated and detrended time series is computed as: Equation 2This procedure is repeated for all time window sizes and with 50% overlap between windows to estimate how the average fluctuation 〈

*F*(τ)〉 scales with window size. The scaling is often of a power-law form: Equation 3The scaling exponent

*α*, also termed the “self-similarity parameter” (Lux and Marchesi, 1999), is extracted with linear regression in double-logarithmic coordinates using a least-squares algorithm. A self-similarity parameter of

*α*= 0.5 characterizes the ideal case of an uncorrelated signal, whereas 0.5 <

*α*< 1.0 indicates power-law scaling behavior and the presence of temporal correlations over the range of τ, where Equation 3 is valid. Periodic signals have

*α*= 0.0 for time scales larger than the period of repetition. The above procedure is illustrated in Figure 1.

*Reference data.* Broadband environmental noise is often temporally correlated. To verify that intrinsic sensor noise or environmental disturbances did not cause any of the effects reported here, a 20 min MEG recording without a subject in the instrument was performed and subjected to identical analyses as the real data.

*Surrogate data.* For the EEG data, we used so-called surrogate data as control, which are commonly used as a control when probing a signal for a nonrandom temporal structure (Ivanov et al., 1996). Surrogate signals have identical power spectra with the original signals but do not have temporal correlations; they were generated by first computing the Fourier transforms of the original signals, randomizing the Fourier phases while preserving the moduli, and then performing inverse Fourier transforms.

## RESULTS

### Oscillatory activity in occipital and rolandic regions

Amplitude spectra were computed for the 122 MEG and the 64 EEG channels. For all 10 subjects and both conditions, we detected prominent peaks in the alpha frequency band (8–13 Hz) in MEG and EEG channels over the occipital and parietal regions (Fig.2*A*,*B*). For both MEG and EEG data, we selected the four channels with the largest alpha rhythm amplitude for further analysis for each subject and condition (the same channels were used for the “eyes-open” and the “eyes-closed” conditions). The peak alpha frequency was 10.4 ± 0.6 Hz (mean ± SD).

Mu rhythm (8–13 Hz) was detected in MEG channels over the right somatosensory cortex in nine subjects (Fig. 2*C*). Additionally, in these subjects, one (∼21 Hz in six subjects) or two (∼16 and ∼21 Hz in three subjects) peaks in the beta frequency band (15–25 Hz) were observed over the somatosensory region (Fig.2*C*). The four channels with the largest amplitude of mu rhythm were selected for further analysis for each subject; for the three subjects having two peaks in the beta range, we analyzed the component having the higher frequency. The peak frequencies were 10.7 ± 0.5 Hz (mu rhythm) and 21.3 ± 1.2 Hz (beta rhythm).

In terms of scaling analysis, the pronounced peaks in the amplitude spectra at 10 Hz show that the dynamics of broadband spontaneous activity is not scale-free; rather, it is dominated by a characteristic time scale of ∼100 msec. In the following sections, we address whether also the amplitude fluctuations of spontaneous oscillations have characteristic scales, which would imply a typical duration of oscillatory bursts.

### Fractal appearance of spontaneous alpha oscillations

Wavelet analysis was used to estimate the amplitude of the signals in narrow frequency bands (Fig. 3) (see Materials and Methods). The wavelet was centered at the peak frequency of a given frequency band determined from the amplitude spectra of individual subjects. Highly irregular amplitude fluctuations were observed in both conditions for the occipital MEG alpha rhythm (Fig.3*A*,*B*). To visualize the trend of the amplitude fluctuations at different scales of temporal resolution, the wavelet-filtered original and surrogate signals were first down-sampled from 300 to 15 Hz. Both the original and the surrogate signals were highly irregular at time scales <12 sec (Fig.3*C*, *top*). Down-sampling the signals to 1.5 Hz and enlarging the displayed interval to 120 sec reveals larger variability in the alpha activity than for the surrogate data (Fig. 3*C*,*middle*). Finally, the display of the entire 1200 sec at a resolution of ∼10 sec still shows large amplitude changes for the alpha but only minor ripples for the surrogate data (Fig.3*C*, *bottom*).

The appearance of large variability at many scales (as in Fig.3*C*) is epitomical of fractal objects and is increasingly being acknowledged to hint about the presence of spatial and temporal correlations at many scales (Bassingthwaighte et al., 1994;Barabási and Stanley, 1995; Bak, 1997). This is in contrast to the variability of signals from uncorrelated or strongly noise-dominated processes that appear even when measured at coarse scales.

### Temporal correlations of spontaneous alpha oscillations

To quantify the temporal structure of the alpha rhythm amplitude fluctuations, we used power spectrum and autocorrelation analyses. Power spectrum analysis measures the contribution of different frequencies to the total power of a signal. In the amplitude envelope of alpha oscillations, the presence of preferred modulation frequencies of oscillations would thus produce peaks in the power spectrum*P*(*f*). We, however, observed a linear decay of power with increasing frequency in double-logarithmic coordinates in the range of 0.005–0.5 Hz; i.e., a 1/*f*^{β} type of a power spectrum: *P*(*f*) ≈*f*^{−β} (Fig.4). For the MEG data, power-law exponents were *β*_{closed} = 0.44 ± 0.09 (mean ± SD; *r*^{2} = 0.94) and *β*_{open} = 0.52 ± 0.12 (*r*^{2} = 0.89). The reference recording gave rise to a white-noise spectrum with*β*_{ref} = 0.03 (*r*^{2} = 0.02), thus ruling out 1/*f*^{β} type of modulation of background 10 Hz noise. The EEG data yielded exponents*β*_{closed} = 0.36 ± 0.17 and*β*_{open} = 0.51 ± 0.12. The differences in exponents between the two conditions (eyes-closed vs eyes-open) or between recording modalities were not significant (two-tailed *t*-tests; all nonsignificant differences in this paper have a *p* > 0.1). As a further control, we used surrogate data (see Materials and Methods); this also resulted in a white-noise spectrum: *β*_{sur} = 0.05 (*r*^{2} = 0.05). The 1/*f*^{β} power spectra indicate a lack of a characteristic time scale for the duration and recurrence of oscillations and are characteristic for fractal signals.

We then computed autocorrelation functions for the wavelet-filtered MEG and EEG data. The autocorrelation analysis indicated the presence of statistically significant correlations up to time lags of >100 sec (Fig. 5). The decay of the autocorrelation function was slow over two decades and well fitted by a power law: ACF(*t*) ≈*t*^{−γ} (Fig. 5). The MEG data yielded γ_{closed} = 0.58 ± 0.23 (*r*^{2} = 0.99) and γ_{open} = 0.73 ± 0.31 (*r*^{2} = 0.97), whereas the EEG data yielded exponents γ_{closed} = 0.52 ± 0.35 (*r*^{2} = 0.97) and γ_{open} = 0.81 ± 0.32 (*r*^{2} = 0.98). The behavior of the autocorrelation functions is in congruence with 1/*f*^{β} type of power spectra. The differences between the two conditions and between the exponents derived from MEG versus EEG data were not significant.

These results indicate that the irregular patterns of amplitude fluctuations of alpha oscillations evident from Figure 3 are embedded with correlations at many time scales. The decrease in correlation with temporal distance appears to be governed by a power-law.

### Spontaneous alpha oscillations exhibit robust power-law scaling behavior

The power spectrum analysis and autocorrelation analyses used in the previous section are not optimal for the quantification of correlations in potentially nonstationary data, because long-range correlations (revealed by an autocorrelation analysis) can arise also as an artifact of the “patchiness” of nonstationary data. Thus, to further consolidate the presence of long-range correlations, we implemented the detrended fluctuation analysis (see Materials and Methods).

DFA was applied to the same amplitude time series of alpha oscillations as analyzed in the previous section. The self-similarity parameter *α* of the DFA is the power-law exponent characterizing the temporal correlations; uncorrelated signals yield a self-similarity parameter *α* = 0.5. This was confirmed using identically wavelet-filtered reference recordings and surrogate data, which yielded *α*_{ref} = 0.508 and *α*_{sur} = 0.496, respectively (Fig. 6*A*,*B*). For the α oscillations, on the other hand, we discovered robust power-law scaling behavior across conditions and recording modalities in 10 of 10 subjects (Fig.6*A*,*B*). The onset of the log-log linear increase of the DFA-fluctuation parameter, *F*, was at a window size of ∼5 sec (this is the lower limit in the DFA method because of the integration by the wavelet in the time domain), and the power-law scaling persisted until at least 300 sec. To obtain reliable scaling statistics for time scales larger than ∼300 sec, longer data series would be needed because of the generally large variability of the root-mean-square fluctuation from one window to the other. The MEG data yielded*α*_{closed} = 0.71 ± 0.06 and*α*_{open} = 0.71 ± 0.05 for conditions eyes-closed and eyes-open, respectively (Fig.6*A*). The EEG data yielded*α*_{closed} = 0.68 ± 0.07 and*α*_{open} = 0.70 ± 0.04 (Fig.6*B*). The difference in self-similarity parameters between the two conditions was not statistically significant, and very similar self-similarity parameters were obtained also for the two recording modalities, despite the different sensitivity of MEG and EEG to the underlying currents (Hari and Ilmoniemi, 1986).

The self-similarity parameter *α* may be viewed as an index of the dynamics of the neural oscillations, whereas the mean amplitude relates to the strength of oscillatory activity. That these two measures convey complementary information about neural activity was indicated by the analysis of their correlation (Fig. 6*C*). The mean amplitude and *α*_{open} were not correlated in either MEG or EEG data, whereas*α*_{close} was weakly, albeit significantly, negatively correlated with the mean amplitude for both MEG and EEG (*p* < 0.04;*r*^{2} < 0.51; Spearman correlation). These correlations are surprising because a decrease in signal-to-noise ratio biases the estimated self-similarity parameter toward that of the reference recording (*α*_{ref} ≅ 0.5). This thus indicates that noise (either environmental or from the sensors) has negligible contribution to the self-similarity parameters estimated for alpha oscillations.

To quantify the apparent stability of the self-similarity parameters across subjects, conditions, and recording techniques relative to the variability of the mean amplitudes, we compared for each subject the ratio*α*_{closed}/*α*_{open}with the corresponding ratio of the mean oscillation amplitude. This normalization eliminates amplitude effects caused by intersubject variability in head size, position in the instrument, etc. The amplitude ratio varied considerably across subjects but was always larger than unity (MEG, 1.98 ± 0.87; EEG, 1.81 ± 0.98) (Fig. 6*D*), reflecting the well known high level of alpha rhythm activity when eyes are closed (Fig. 2*A*). The ratios of scaling exponents, on the other hand, were near unity (MEG, 1.04 ± 0.12; EEG, 1.00 ± 0.13) (Fig.6*D*). Linear correlations between amplitude and scaling exponent ratios were nonsignificant in both MEG and EEG data, and the SD of the amplitude ratios was significantly larger than for the exponent ratios (*p* < 0.0001; Fisher's test).

The DFA results indicate that, for the 10 Hz oscillations from the occipital region, spontaneous activity is robustly characterized by long-range temporal correlations that decay as power-law functions and with remarkably invariant scaling exponents. It has been pointed out recently that the scaling parameters of genuine long-range correlated processes obey the following relation: *α* = (2 −*γ*)/2 = (1 + *β*)/2 (Rangarajan and Ding 2000). Using *γ* and *β* from the previous section, good agreement is found for the predicted and the measured*α*. Thus, together, DFA, autocorrelation, and power spectrum analyses provide robust evidence in support for power-law scaling behavior, as well as for the lack of characteristic time scales for the modulation of the alpha rhythm amplitude.

### Generality of long-range temporal correlations and power-law scaling behavior of spontaneous oscillations

To test whether scaling behavior was unique to alpha rhythmicity or a more general property of large-scale network oscillations, we applied the DFA and autocorrelation analyses to oscillations detected with MEG at the mu and beta frequency bands over the right somatosensory region in the eyes-closed condition (Fig.2*C*).

Robust power-law scaling was indeed evident for both mu and beta oscillations over a range of approximately two decades. The self-similarity parameters obtained for mu and beta oscillations were significantly different: *α*_{mu} = 0.73 ± 0.09 and *α*_{beta} = 0.68 ± 0.07 (*p* < 0.005) (Fig.7*A*); however, comparing*α*_{mu} with*α*_{closed} (the occipital alpha) indicated that 10 Hz oscillations from the rolandic and occipital regions had similar scaling properties (significance level of the difference, *p* > 0.3). The power-law decays of the autocorrelation functions were characterized by exponents γ_{mu} = 0.46 ± 0.35 and γ_{beta} = 0.70 ± 0.36 (Fig. 7*B*); the difference in these exponents, as well as the differences in mean amplitudes of the mu and beta oscillations, were significant (*p* < 0.04). Nevertheless, correlations were not found between the magnitude of the self-similarity parameters and of the mean amplitudes for either the mu or beta rhythms, nor did the ratios of the exponents and of the amplitudes correlate linearly (Fig.7*C*). The lack of correlation between the self-similarity parameters and amplitudes makes it unlikely that the difference in scaling exponents results from the lower signal-to-noise ratio of beta oscillations. Differential scaling parameters of mu and beta oscillations suggest that the neural mechanisms and/or networks underlying these two rhythms are distinct. This interpretation is in agreement with reports on differential reactivity and anatomical origin of somatosensory mu and beta oscillations (Hari and Salmelin, 1997). In line with the results on alpha oscillations, we also found for the somatosensory oscillations that the ratios of mu and beta rhythm scaling exponents were more stable than the corresponding mean amplitude measures (exponent ratio, 1.09 ± 0.08; amplitude ratio, 1.41 ± 0.18; the difference in SD of the ratios,*p* < 0.002; Fisher's test).

The presence of power-law scaling behavior in amplitude fluctuations of mu and beta frequency bands in the somatosensory region indicates that these statistical characteristics are not unique to the occipitoparietal alpha rhythm.

## DISCUSSION

We have investigated the large-scale dynamics of network oscillations in the normal human brain. To the best of our knowledge, this is the first characterization of the temporal correlations in spontaneous oscillations at time scales ranging from a few seconds to several minutes. Our results indicate that spontaneous alpha, mu, and beta oscillations have significant temporal correlations for at least a couple of hundred seconds during resting conditions (eyes-open and eyes-closed). The decay of correlation was characterized by power-law scaling. The self-similarity parameters obtained with the detrended fluctuation analysis were highly invariant across subjects. The simultaneously recorded MEG and EEG agreed quantitatively in their estimates of the scaling exponents characterizing the occipital alpha rhythm dynamics. Oscillations at 8–13 and 15–25 Hz had different scaling properties, which suggests that distinct neural networks and/or mechanisms underlie these oscillations.

### The correlated nature of spontaneous oscillations

It has often been noted that spontaneously occurring synchrony in cell populations appears in an irregular manner both temporally and spatiotemporally (Traub et al., 1989; Destexhe et al., 1999; Tsodyks et al., 1999). Previous studies have reported that 10 Hz oscillations are generated with great variability every 5–20 sec and last for ∼0.5–10 sec (Lopez da Silva, 1991; Destexhe et al., 1998). It has remained unknown, however, to what extent oscillatory activity beyond these time scales is statistically dependent. The present scaling analyses indicate that successive oscillations indeed are correlated, even over thousands of oscillation cycles (Figs. 4-7).

Scaling analysis is used increasingly in many fields of science to characterize complex phenomena. It can be used to test a model for its ability to generate scale-free behavior (Ivanov et al., 1998). Alternatively, transitions in scaling behavior from one parameter range to another can reveal scales at which different mechanisms influence the system dynamics (Barabási and Stanley, 1995; Peng et al., 1995; Ivanov et al., 1996). The stability of scaling exponents obtained here (Figs. 6*C*,*D*, 7*C*) suggests that scaling exponents may indeed be useful quantitative hallmarks of also the dynamic processes underlying spontaneous brain oscillations. The good quantitative agreement of the scaling exponents derived from MEG and EEG data reflects that scaling exponents are quantitative indices of relative fluctuations and do not depend on the unit of choice or the method used to measure the underlying dynamic process.

Moreover, the results of the power spectrum analysis indicated that bursts of oscillations are not modulated at any characteristic or dominant time scale (Fig. 4). The correlated nature of these oscillations suggests that “a burst” is only a part of a series of connected events and that the fractal structure of the signal reflects a hierarchy of bursts within bursts rather than a succession of individual or independent bursts. This we shall address further in the next section.

### Local interactions as a mechanism underlying long-range temporal correlations and scaling behavior

One of the defining aspects of population oscillations is the ability of neural networks to establish spatiotemporal correlations with millisecond range precision and over large distances mainly through local synaptic connections (Traub et al., 1989). Here we describe a general framework of how local interactions create large-scale dynamics, which could account for the long-range temporal correlations and the power-law scaling behavior observed for spontaneous oscillations.

Since the first reports on self-organized criticality (Bak et al., 1987, 1988), ample evidence has indicated that several complex systems self-organize through local interactions to a critical state with long-range spatiotemporal correlations (Bak et al., 1989; Mantegna and Stanley, 1995; Boettcher and Paczuski, 1996; Paczuski et al., 1996;Bak, 1997; Malamud et al., 1998; Lux and Marchesi, 1999). This state is termed “critical” because similar scaling behavior can be observed in equilibrium systems when fine-tuning a parameter to the point at which the system undergoes a phase transition. In nonequilibrium systems, however, this complex state can be “self-organized” and emerge purely under the dynamics of the system. In this case, the local rules of interaction sculpt the dynamics across many scales, and no characteristic scale can be identified.

Neural networks host the common features of SOC systems, such as a large number of units (neurons), local and nonlinear interactions between neurons, externally imposed perturbations, a certain amount of stochastic variation of internal parameters, and ability to store information in spatial patterns. In analogy with the scale-free behavior of SOC systems, we propose that the power-law form of the amplitude dynamics of spontaneous oscillations may not be highly dependent on the specific mechanisms underlying the generation of the population oscillations. A fundamental prerequisite for the emergence of a critical state, however, is that the network oscillations are associated with synaptic plasticity. Structural memory affecting the recruitment of neurons into future population oscillations is critical to ensure a degree of correlation. The exact values of the power-law exponents, on the other hand, may be related to both the biophysical mechanisms and the neural architecture underlying the oscillations. In line with this, our results showed that mu and alpha oscillations scaled similarly, but beta oscillations had a significantly smaller scaling exponent than mu and alpha.

A self-organized critical process, as the source of the temporal power laws, would further suggest that similar power laws exist also for parameters in the spatial domain. Based on the analogy with other SOC systems, one prediction is power-law statistics for the probability that a number of neurons, *n*, is recruited into an oscillatory event. Quantification of spatial correlations may, however, require invasive studies with greater spatial resolution than the present MEG and EEG measurements.

### Possible functional significance of self-organized scale-free dynamics

The functional significance of the scale-free behavior of oscillations may be diverse. Temporal correlations of spontaneous network oscillations, as we have described here, may be the physiological correlate of behavioral results such as the 1/*f*^{β} noise in the human judgement of temporal intervals (Gilden et al., 1995) and the long-range correlations observed for synchronization errors in human coordination (Chen et al., 1997). Thus, in terms of mechanisms, it may be the dynamic structural memory of the neural networks (see previous section) that constrain perception and behavior to power-law statistical patterns, even in situations in which humans attempt to avoid such correlations.

From a theoretical point of view and based on simulations, it has been argued that a state of criticality would be optimal for a network to swiftly adapt to new situations (Alstrøm and Stassinopoulos 1995;Stassinopoulos and Bak, 1995; Bonabeau, 1997; Chialvo and Bak, 1999). In the critical state, the spatiotemporal correlations are highly susceptible to perturbations; the dynamics may be viewed as balancing between a predictable stable pattern of activity and uncorrelated random behavior. Thus, if the “fractal” structure of neural oscillations indeed arises from self-organized neural network dynamics poised at criticality, one would expect the ongoing activity to be effectively disrupted by externally imposed perturbations. This, in fact, has been observed. During event-related desynchronization (ERD), spontaneous 8–13 Hz oscillations are suppressed rapidly (approximately within one cycle) by sensory stimulation (Hari and Salmelin, 1997;Nikouline et al., 2000), memory search (Kaufman et al., 1991), or motor activity (Pfurtscheller, 1989; Crone et al., 1998). Furthermore, the mapping of ERD on the cortical surface has revealed transitions from spatially diffuse to focused and somatotopically specific patterns of alpha suppression (Crone et al., 1998), consistent with the picture of spontaneous cortical states being driven into stimulus specific configurations of correlated neural activity (Tsodyks et al., 1999). We suggest that the widespread and rapid onset of ERD reflects long-range spatial correlations in the neural networks. Because all oscillations studied here showed surprisingly robust scaling behavior, we tentatively propose that, under normal physiological conditions, neural networks in general may operate in a critical state, thereby lending them capable of quick reorganization during processing demands.

Further studies are needed to determine how the power-law scaling exponents are affected by various experimental, pharmacological, or pathological conditions and whether current computational models of spontaneous network oscillations agree qualitatively and quantitatively with the present findings.

## Footnotes

This work was supported by the Danish Natural Science Research Council, The Danish Research Agency, Center for International Mobility (Helsinki), and Helsinki University Central Hospital Research funds. J.M.P. was supported by the Academy of Finland and by the Juselius Foundation. We thank T. L. van Zuijen, J. Sinkkonen, and M. Kesäniemi for discussions. Wavelet software provided by C. Torrence and G. Compo (http://paos.colorado.edu/research/wavelets) was modified for this study.

Correspondence should be addressed to Klaus Linkenkaer-Hansen, BioMag Laboratory, P.O. Box 442, FIN-00029 HUS, Finland. E-mail:klaus{at}oliivi.huch.fi.