Abstract

Although the emerging field of functional connectomics relies increasingly on the analysis of spontaneous fMRI signal covariation to infer the spatial fingerprint of the brain's large-scale functional networks, the nature of the underlying neuro-electrical activity remains incompletely understood. In part, this lack in understanding owes to the invasiveness of electrophysiological acquisition, the difficulty in their simultaneous recording over large cortical areas, and the absence of fully established methods for unbiased extraction of network information from these data. Here, we demonstrate a novel, data-driven approach to analyze spontaneous signal variations in electrocorticographic (ECoG) recordings from nearly entire hemispheres of macaque monkeys. Based on both broadband analysis and analysis of specific frequency bands, the ECoG signals were found to co-vary in patterns that resembled the fMRI networks reported in previous studies. The extracted patterns were robust against changes in consciousness associated with sleep and anesthesia, despite profound changes in intrinsic characteristics of the raw signals, including their spectral signatures. These results suggest that the spatial organization of large-scale brain networks results from neural activity with a broadband spectral feature and is a core aspect of the brain's physiology that does not depend on the state of consciousness.

Introduction

In the absence of sensory input, much of the brain remains remarkably active, including even cortical regions dedicated to the processing of sensory signals. A large body of evidence, drawing primarily from functional magnetic resonance imaging (fMRI) studies, has shown that large-scale brain activity measured in the absence of explicit sensory input (here referred to as “spontaneous” activity) exhibits spatial covariations that resemble the functional networks that co-activate during behavioral tasks (Fox and Raichle 2007). Although the study of this activity may offer opportunities to chart the brain's major functional networks under both healthy (Fox and Raichle 2007) and pathological (Zhang and Raichle 2010) conditions, neither its neural origins nor ultimate role is well understood (Leopold and Maier 2012).

Further understanding will almost certainly require electrophysiological studies, which, in principle, allow a more direct measurement of neuronal signaling and provide information not available with fMRI. For example, electrophysiological recordings better capture the temporal structure of neuronal activity, allowing the study of spectral difference that may exist between functional networks (Siegel et al. 2012). Indeed, a number of electrophysiological studies have investigated spontaneous activity (He et al. 2008; Nir et al. 2008; Breshears et al. 2010; Liu et al. 2010; Schölvinck et al. 2010; Brookes et al. 2011; Hipp et al. 2012; Wang et al. 2012) and found some evidence for a spectral distinction between the spontaneous activity in different brain regions (He et al. 2008; Nir et al. 2008; Breshears et al. 2010; Schölvinck et al. 2010; Hipp et al. 2012; Wang et al. 2012). Nevertheless, the relationship between these spectral features and the functional specificity of the fMRI networks remains unclear. This is in part due to the fact that the electrophysiological studies that best capture spectral information require invasive (intracranial) recordings and often have very limited spatial coverage of the electrophysiological signals. In addition, the region-by-region-based analysis approach commonly applied in these studies does not allow capturing the correlational structure of the activity in an unbiased manner (He et al. 2008; Nir et al. 2008; Breshears et al. 2010; Liu et al. 2010; Schölvinck et al. 2010; Hipp et al. 2012; Wang et al. 2012).

To address these shortcomings, we investigated large-scale covariation of electrophysiological signals with electrocorticography (ECoG) arrays implanted over an extensive portion of the left hemispheric surface of macaque monkeys, whose spontaneous brain activity has been extensively studied with fMRI (Vincent et al. 2007; Moeller et al. 2009; Hutchison et al. 2011; Hutchison and Everling 2012; Wey et al. 2014) and was found to have covariational structures that resemble those found in human (Hutchison and Everling 2012; Wey et al. 2014). The ECoG data, which covered a broad spectral range, were analyzed with a novel data-driven clustering approach, which has previously been used for analyzing resting-state fMRI signals (Liu et al. 2012). This was done across various behavioral states, including eyes-closed wakefulness, sleep, ketamine/medetomidine anesthesia and propofol anesthesia.

Materials and Methods

Subjects and Materials

A chronically implanted customized multichannel ECoG electrode array (Unique Medical, Japan) was employed for neural recording (Nagasaka et al. 2011). Each electrode consisted of a 3-mm diameter platinum disc that was insulated with a layer of silicone except for a small (0.8 mm diameter) disc-shaped area at its center. The array covered a large cortical area with approximate inter-electrode distance of 5 mm. The array was implanted in the subdural space in 4 adult macaque monkeys (monkeys C, G, and K are Macaca fuscata, and monkey S is Macaca mulatta). The position of ECoG electrodes in monkeys C, K, and S was manually remapped onto the brain of monkey G (Figs 1a and 2a). 128-channel ECoG electrodes were implanted in the left hemisphere, covering a contiguous area that included frontal, parietal, temporal, and occipital lobes (Fig. 2a). A reference electrode consisting of a rectangular platinum plate was placed in the subdural space between ECoG array and dura mater. A ground electrode was placed in the epidural space (see [Nagasaka et al. 2011] for the detailed description). ECoG signal was recorded at a sampling rate of 1 kHz by the Cerebus data acquisition system (Blackrock, UT, USA).

Figure 1.

Analysis flow illustrated with data from a representative session. (a) Spatial topography of 128 electrodes in monkey S with 4 exemplary ones (14, 23, 31, and 53) marked in color; (b) raw ECoG signals acquired during a 5-min eyes-closed session with an enlarged 5-s segment for the exemplary electrodes; (c) preprocessed spectrograms of the exemplary electrodes suggest relationship not simply explained by volume conduction effect: distant pairs (14 versus 23, and 31 versus 53) show similar spectrograms whereas the adjacent pair (23 versus 31) displays distinct patterns; (d) BLPs can be derived by averaging the spectrogram within specific frequency ranges; (e) cross-electrode correlations of the spectrogram or BLPs (β-BLP for the example here) form a matrix whose row (and column) values represent correlation values of specific electrode with all others (here called profile); map representation of correlation profiles of the exemplary electrodes (white stars) suggests 2 covariation structures whose boundary lies between the electrodes 23 and 31; (f) clustering electrodes based on the similarity of their correlation profiles and then reordering the matrix accordingly reveal several diagonal blocks that correspond to different covariation structures, whose topographic distribution and correlation strength can be then represented with group-averaged correlation maps.

Figure 1.

Analysis flow illustrated with data from a representative session. (a) Spatial topography of 128 electrodes in monkey S with 4 exemplary ones (14, 23, 31, and 53) marked in color; (b) raw ECoG signals acquired during a 5-min eyes-closed session with an enlarged 5-s segment for the exemplary electrodes; (c) preprocessed spectrograms of the exemplary electrodes suggest relationship not simply explained by volume conduction effect: distant pairs (14 versus 23, and 31 versus 53) show similar spectrograms whereas the adjacent pair (23 versus 31) displays distinct patterns; (d) BLPs can be derived by averaging the spectrogram within specific frequency ranges; (e) cross-electrode correlations of the spectrogram or BLPs (β-BLP for the example here) form a matrix whose row (and column) values represent correlation values of specific electrode with all others (here called profile); map representation of correlation profiles of the exemplary electrodes (white stars) suggests 2 covariation structures whose boundary lies between the electrodes 23 and 31; (f) clustering electrodes based on the similarity of their correlation profiles and then reordering the matrix accordingly reveal several diagonal blocks that correspond to different covariation structures, whose topographic distribution and correlation strength can be then represented with group-averaged correlation maps.

Figure 2.

Spatial patterns of broadband ECoG power covariation. (a) Spatial topography of electrode coverage in individual monkeys is presented on cortical surface rendering reconstructed from the anatomical MRI of Monkey G. Exemplary ECoG raw signals from 2 electrodes in the Monkey C underwent a dramatic transition from a high-frequency, low-amplitude voltage pattern during eyes-closed wakefulness to a low-frequency, high-amplitude one under ketamine/medetomidine anesthesia (b), which is confirmed by distinct power spectra under these 2 conditions (c). However, correlation patterns of broadband ECoG power during (d) the eyes-closed wakefulness and (e) ketamine/medetomidine anesthesia similarly resemble (f) the fMRI resting-state networks (RSNs) from isoflurane-anesthetized macaques (adapted from [Hutchison et al. 2011]). Corresponding patterns or RSNs are aligned in the same row.

Figure 2.

Spatial patterns of broadband ECoG power covariation. (a) Spatial topography of electrode coverage in individual monkeys is presented on cortical surface rendering reconstructed from the anatomical MRI of Monkey G. Exemplary ECoG raw signals from 2 electrodes in the Monkey C underwent a dramatic transition from a high-frequency, low-amplitude voltage pattern during eyes-closed wakefulness to a low-frequency, high-amplitude one under ketamine/medetomidine anesthesia (b), which is confirmed by distinct power spectra under these 2 conditions (c). However, correlation patterns of broadband ECoG power during (d) the eyes-closed wakefulness and (e) ketamine/medetomidine anesthesia similarly resemble (f) the fMRI resting-state networks (RSNs) from isoflurane-anesthetized macaques (adapted from [Hutchison et al. 2011]). Corresponding patterns or RSNs are aligned in the same row.

All experimental and surgical procedures were performed in accordance with the experimental protocols (No. H24-2-203(4)) approved by the RIKEN ethics committee and the recommendations of the Weatherall report, “The use of non-human primates in research”.

Experimental Procedures

A total of 20 experiments were conducted over different days and monkeys (Supplementary Table 1). Each experiment included the eyes-closed waking condition and 1 of the following 3 conditions: ketamine/medetomidine anesthesia, propofol anesthesia, or natural sleep. Ketamine is an N-methyl-D-aspartate antagonist known to induce a state of anesthetic dissociation. It is commonly combined with the α-2 adrenergic agonist medetomidine to improve muscle relaxation and analgesia. During the ketamine/medetomidine condition, monkeys were seated in a dedicated chair with arms and head movement restricted and eyes covered by an eye mask. The monkey was sitting calmly, and the ECoG signals were recorded for up to 20 min during this eyes-closed condition. After that, we started to monitor heart rate. Ketamine (∼5.0 mg/kg for Macaca fuscata, ∼8.8 mg/kg for Macaca mulatta) and medetomidine (∼0.016 mg/kg for Macaca fuscata, ∼0.053 mg/kg for Macaca mulatta) were then given through intramuscular injection whereas heart rate was kept below 120 beats per minute. The loss of consciousness (LOC) was defined behaviorally: this was done by recording the time point at which the monkey became unresponsive (i.e., absence of somatic movement) to touching its nostril and philtrum by a rolling pin and to opening the monkey's hands. After LOC was established, neural activity, which is characterized by slow-wave activity, was recorded for ∼25 min with heart rate and breathing monitored carefully. In contrast to the ketamine/medetomidine combination, propofol is thought to induce general anesthesia primarily by potentiating GABAA ion channels. The experiments with the propofol anesthesia were conducted in a similar manner with the intravenous injection of propofol (5.2 mg/kg for Macaca fuscata) instead of ketamine/medetomidine.

For the experiments that included the natural sleep condition, electro-oculography (EOG) and electromyography (EMG) signals were also recorded at a sampling rate of 1 kHz by the Cerebus data acquisition system. The EOG signal was recorded from muscles of the right eye in the horizontal direction. Two electrodes (Nihon Kohden, Disposable electrode for ECG Monitoring M-150) were put on the tail and the inner corner of right eye, respectively, and the potential difference between them was used as EOG signal. Another 2 electrodes were put on the chin, and their potential difference was used as EMG signal. The data were first acquired continuously for up to 1.5 h with the eyes-covered monkey sitting in a quiet and dark environment, and the slow-wave oscillations appeared intermittently in ECoG during this time period. After that, the light was turned on and the experimenter made noise to wake up the monkey, and the noise was repeated as long as the slow-wave activity appeared in the real-time ECoG recording. The data acquisition was continued under this eyes-closed condition for ∼10 min. After that, the eye mask was removed and the data were acquired for another 10 min with the monkey's eyes opened.

Definition of Sleep State

In order to distinguish waking and sleeping conditions, we monitored the presence of slow -wave activity. We did not attempt to distinguish different sleep states, as we had difficulty in detecting clear sleep spindles and K-complexes that would be required for this. Nevertheless, we observed strong, intermittent slow-wave oscillations, which were generally associated with decreased EOG and EMG signals. To define sleep, we quantified slow-wave activity as follows.

First, we estimated δ-band (1 − 4 Hz) ECoG power within 1-s time bins (shifted every 200 ms) using the multi-taper method (Thomson 1982) implemented in Chronux, a Matlab package for the analysis of neural data (Mitra and Bokil 2008). This δ-band power was normalized by the median value of δ-band power under the eyes-opened condition. The normalized δ-band power was binarized according to an arbitrary threshold of 3, which we regarded as a level of significant slow-wave activity. We further estimate the spatial extent to which the slow-wave oscillations spread using a spatial synchronization index (SSI), defined as the cross-electrode average of the binarized δ-band power. A high SSI is associated with strong slow-wave activity (Supplementary Fig 1) as well as low levels of EOG and EMG (Supplementary Figs 1 and 2).

We defined the sleeping state as time periods when the SSI exceeded 0.3, which means that at this stage more than 30% electrodes have δ-band activity at least 3 times higher than the median level of δ-band activity under the eye-opened condition.

Temporal Variation of ECoG Power

The raw ECoG signals were first visually inspected. Serious line-noise contaminations were observed for 1 or 2 channels in some experiments, and these channels were excluded from further analysis. For the remaining channels, the primary frequency (50 Hz) and harmonics (multiples of 50 Hz) of the line noise were removed using a sine-wave-fitting method provided by Chronux. The ECoG signals were re-referenced to the common mean of all the channels, with the data from 1 experiment of Monkey S also being re-referenced using Laplacian method separately (Nunez and Pilgreen 1991). Examination of the raw voltage traces from the electrode array (Fig. 1b) revealed that ECoG activity was coordinated across recording sites over a large range of spatial and temporal scales, consistent with findings from other studies (Nir et al. 2008; Liu et al. 2010; Hipp et al. 2012). Based on previous work, we focused analysis on slow fluctuations of band-limited power signals (BLPs), which may have particular relevance to fMRI hemodynamic signal (Logothetis et al. 2001). Multi-taper time-frequency transformation (Thomson 1982) (parameters: window length, 1 s; step, 0.2 s; and the number of tapers, 5) was applied to each channel using Chronux functions, resulting in spectrograms representing temporal evolution of ECoG powers at any specific frequency between 1 and 100 Hz with a spectral resolution of ∼1 Hz. The power was converted into decibel (dB) units with a logarithmic operation. To balance the degrees of freedom for the subsequent correlation analysis, the spectrograms were divided into 5-min sessions. For the sleep condition, the time periods not defined as the sleep (the time periods with the SSI <0.3) were removed before dividing the whole spectrogram into multiple sessions.

To compensate for the unequal spectral power of the ECoG signal at different frequencies (often following an ∼1/f dependency, in part resulting from mechanistical (nonbiological) origins inherent to the recording technique), the ECoG power at each frequency bin was centered by removing the mean and then normalized by the temporal standard deviation; as a result, the normalized spectrograms showed the same fluctuation amplitude for each frequency. To reduce contributions from global power changes of non-neuronal origin, the cross-electrode average was subtracted from the spectrogram of each electrode. It is worth noting that this procedure may remove some global signals with neuronal origin as well (Schölvinck et al. 2010); however, this was not a major concern as the main goal of the analysis was the characterization of region-specific information. To derive BLP signals, each normalized spectrogram was averaged within the following frequency bands: δ, 1 − 4 Hz; θ, 5 − 8 Hz; α, 9 − 12 Hz; β, 13 − 30 Hz; and γ, 31 − 100 Hz. Covariation of the resulting signals provided a measure of synchrony over timescales matching those of hemodynamic changes and allowed for an assessment of the spatial organization of slow, spontaneous neural modulations in the hemisphere. Note that this analysis method differs the more common practice of computing neural covariation based on the moment-to-moment changes of the raw traces (e.g., magnitude-squared coherence of the voltage), which depends upon precise and consistent phase relationship over millisecond timescale.

Clustering-Based Correlation Analysis

The correlation analysis was adapted from a clustering method used previously for the analysis of resting-state fMRI signals (Liu et al. 2012). First, an inter-electrode cross-correlation matrix was generated by determining the correlation between signals from each pair-wise combination of electrodes based on either the normalized spectrogram (correlations were calculated based on the whole spectrogram, referred to as the broadband power later) or BLP, and its rows/columns then represent correlation profiles (spatial correlation maps) of different electrodes (Fig. 1e). If region-specific correlations exist, the electrodes in the same correlation structures are expected to exhibit similar correlation profiles: strong correlations with electrodes inside the structure and weak correlations with electrodes outside (Fig. 1e). Based on this notion, the electrodes were classified into multiple groups based on the similarity of their correlation profiles using k-means clustering method (Forgy 1965). Clustering is a procedure for classifying a set of objects into different groups such that within-group differences are smaller than across-group differences. The k-means clustering partitioned the rows/columns of the correlation matrix {p1, p2, …, pn} into k clusters R = {R1, R2, , Rk} such that the sum of within-cluster distances J is minimized

$J=∑i=1k∑pj∈Rid(pj,μi)$
where μi is the mean of correlation profiles in Ri, and d(•) represents the distance between 2 profiles, defined as 1 minus their Pearson's correlation coefficient.

After clustering, the electrode correlation profiles (the rows/columns of the correlation matrix) were averaged within groups while excluding the diagonal elements (to reduce possible biases caused by auto-correlations), and the resulting maps show how the correlation strength of this group of electrodes varied across space (Fig. 1f). For the group-level analysis, the correlation matrices of individual monkeys were first expanded by spatially interpolating the correlation profiles (rows/columns) to high-dimensional space of the electrodes from all the monkeys (or from Monkey C and G for the propofol anesthesia and sleep conditions). The expanded correlation matrices were averaged across experiments and individuals, and the clustering analysis was then applied to the final average. After reshuffling the rows/columns of the correlation matrix according to the clustering results, block structures along the diagonal appeared, reflecting a structural organization of the inter-region correlations (Fig. 1f).

To determine the number of clusters k, both the Bayesian information criterion (BIC) (Schwarz 1978) and silhouette index (Rousseeuw 1987) were calculated at k values ranging from 2 to 50. While the BIC and silhouette index in general monotonically change as k increases, the changes decelerated when k exceeded values in the range of 6–10. The k-means clustering was then repeated using k values within this narrow range, and the results were compared. It was found that the results for k equal to 9 or 10 tend to include correlation structures with similar correlation profiles/maps. Therefore, k was fixed at 8 for the final analysis.

Power Spectral Analysis

To investigate spectral differences across clusters and conditions, the spectrograms were first averaged over time to obtain power spectra for each electrode. This was done prior to frequency normalization that counteracts the effect of ∼1/f power dependency (see above). The mean of each spectrum was then subtracted to eliminate inter-electrode variations in absolute power level, which could be caused by differences in electrode sensitivity or in brain activity level (Manning et al. 2009). For each experiment, the spectra were then averaged, according to network parcellation results of the clustering analysis, to represent the power spectra for different functional clusters. The results from individual experiments were further averaged to generate the group-level results. The above-mentioned analysis was also repeated without mean subtraction step (Supplementary Fig. 3) or with averaging according to a fixed network parcellation (Supplementary Fig. 4). This fixed network parcellation was an intersection of the parcellations under the eyes-closed condition and ketamine/medetomidine anesthesia, with the exclusion of 1 network that has only a few electrodes after intersection operation.

The power spectra of ECoG BLPs were also investigated. After the removal of the line-noise and re-referencing to the common mean, the continuous ECoG signals acquired under different conditions were band-pass-filtered into the frequency bands specified earlier. The amplitude of the band-pass-filtered signals was extracted using the Hilbert transform, and their power spectra were then estimated using the multi-taper method with a window length of 300 s. The BLP was not extracted directly from the spectrograms because the time-frequency analysis used to generate them introduces an effective low-pass filtering, and as a result, the high-frequency fluctuations in BLPs may be underestimated. The BLP power spectra of different electrodes were then averaged for each experiment according to the network parcellations derived from the clustering analysis and further averaged across experiments to generate the group-level results.

Statistical Analysis of Network-Specific Power Spectra

Based on the division of k electrode groups (clusters), the power spectrum of a particular electrode can be indicated as Sij(f), where i is the index of groups and j is the within-group index for electrodes. For each experiment, we quantified the between-network differences at a specific frequency f0 using the mean between-group-sum-of-square (MBGSS), a statistic commonly used in one-way ANOVA

$MBGSS(f0)=∑i=1kni(S¯i(f0)−S¯(f0))2k−1=∑i=1k(∑i=1kSij(f0))2ni(k−1)−(∑i=1k∑j=1niSij(f0))2n(k−1)$
where ni is the number of electrodes in group i and n is the total number of electrodes. This statistic follows a chi-square distribution with the degrees of freedom equal to k 1. A MBGSS line can be calculated to represent MBGSS values as a function of frequency f. The MBGSS lines can be averaged across experiments but the degrees of freedom needs to be adjusted accordingly by multiplying the number of experiments. The ratio of averaged MBGSS lines under 2 different conditions, whose values follow the F distributions under the null hypothesis, can then be used to test whether the between-network difference is significantly higher under 1 of the conditions.

Similarly, we used the mean residual-sum-of-square (MRSS) to quantify the within-network variations in spectral power.

$MRSS(f0)=∑i=1k∑j=1ni(Sij(f0)−S¯i(f0))2n−k=1n−k(∑i=1k∑j=iniSij(f0)2−∑i=1k(∑j=iniSij(f0))ni2)$

The smaller this quantity is, the better the model fits the data, and the more similar the spectra within the same group are. For a single experiment, this statistic follows a chi-square distribution with the degree of freedom of nk, which also needs to be adjusted when averaging across experiments. The MRSS lines were computed for 2 models (parcellations) R and R, and their ratio was then used to test whether the within-group variation in power spectra is significantly smaller if we divide the electrodes following 1 of the models (Supplementary Fig. 5).

Results

Slow Fluctuations in ECoG Power Replicate fMRI Resting-State Networks

Analysis of the spatio-temporal structure of the ECoG signals based on broadband power (see Materials and Methods) revealed 6–10 clusters of electrodes with similar power variation during eye-closed wakefulness, of which the 8-cluster result was the most robust across animals (see Materials and Methods). These clusters (Fig. 2d) were localized in the prefrontal, frontal opercular, dorsal precentral, superior temporal, supramarginal, parietal, occipito-temporal, and occipital cortex. These correlational patterns could be derived from just a few minutes of data collection in individual animals. In addition, their spatial distribution was virtually independent of signal preprocessing steps such as local electrode voltage re-referencing (to minimize volume conduction effects) or subtraction of the mean power signal (to remove the contribution of global electrophysiological fluctuations) (Supplementary Fig. 6). Visual comparison of the clustering results with a previous fMRI study of the anesthetized macaque (Fig. 2f) revealed a striking correspondence of the cluster patterns with 8 of 11 patterns of fMRI signal covariation (Hutchison et al. 2011). The similarity of these covariational patterns, derived in one case from hemodynamic responses and in the other from the power fluctuations in the ECoG signal, suggests that the fMRI resting-state networks (RSNs) may reflect regional electrophysiological power fluctuations in the cerebral cortex.

To investigate whether the observed correlational patterns originated from electrical activity in specific spectral frequency bands, a possibility indicated by earlier work (Mantini et al. 2007; Hipp et al. 2012; Siegel et al. 2012), we restricted the power signals to conventional spectral bands, commonly known as δ, θ, α, β, or γ (see Materials and Methods) (Fig. 1df). Applying our clustering method to each frequency band individually resulted in maps that were, by and large, similar to one another (Fig. 3 and Supplementary Fig. 7). The most obvious difference was the stronger long-range correlation in the γ-BLP compared with the other frequency ranges, which was observed in coupling between frontal, temporal, and parietal areas in individual monkeys (Fig. 4a), and was also visible in certain γ-BLP-derived clusters at the group level (Fig. 3 and Supplementary Fig. 7: note the uppermost cluster in the γ column). The general similarity of the maps derived from different frequency bands suggests that the neurophysiological mechanisms underlying the power correlations are broadband in nature.

Figure 3.

Spatial patterns of ECoG power covariation at different frequency bands during (a) the eyes-closed wakefulness and (b) ketamine/medetomidine anesthesia. Gray dots indicate electrode locations and white lines show the brain contour.

Figure 3.

Spatial patterns of ECoG power covariation at different frequency bands during (a) the eyes-closed wakefulness and (b) ketamine/medetomidine anesthesia. Gray dots indicate electrode locations and white lines show the brain contour.

Figure 4.

Long-range γ-BLP correlation at individual level. (a) 3 examples of single-session patterns from 3 monkeys consistently show long-range correlations among the frontal, superior temporal, and parietal regions; and (b) the long-range correlation between a few electrodes at the frontal eyes field (FEF, white arrows) and the visual association areas are also consistently observed across individuals.

Figure 4.

Long-range γ-BLP correlation at individual level. (a) 3 examples of single-session patterns from 3 monkeys consistently show long-range correlations among the frontal, superior temporal, and parietal regions; and (b) the long-range correlation between a few electrodes at the frontal eyes field (FEF, white arrows) and the visual association areas are also consistently observed across individuals.

The Effects of Sleep and Anesthesia on ECoG-Derived Networks

Repeating the above-mentioned analysis under conditions of sleep and anesthesia demonstrated that functional clusters were remarkably robust to both natural and artificial changes to the conscious state. For example, the ECoG-derived networks obtained during wakefulness (Fig. 2d) showed a close correspondence to those obtained during ketamine/medetomidine anesthesia (Fig. 2e), as well as propofol anesthesia and natural sleep (Supplementary Fig. 8b). Moreover, and despite significant changes in the overall spectral characteristics of the raw ECoG signals (Fig. 2b,c), it was still possible to derive the same networks from each of the frequency bands (Fig. 3b, Supplementary Fig. 7). This similarity extended to the long-range coupling of γ-range power observed during wakefulness (uppermost clusters in the γ column), which is surprising given the known diminution of γ-range power during diminished states of consciousness (John et al. 2001).

We then examined to what extent the basic characteristics of spontaneous activity are similar or different across different networks. We found that during wakefulness, there was a substantial variation in the spectral composition between the different functional clusters (Fig. 5a, see also Supplementary Fig. 5), which is in agreement with a previous observation of regional variation in the power spectrum of the neural signal over the cortical surface (Rosanova et al. 2009). This variation was most pronounced in the ∼5 − 25 Hz frequency range, corresponding roughly to the α and β bands, which showed highest relative power in the occipital and lateral cortex and the lowest in the frontal and medial areas, consistent with previous findings (Rosanova et al. 2009). Importantly, these spectral signatures were unique to the waking state. During natural sleep and under each of the 2 anesthetic regimes, the entire cortex was characterized by a distinct, state-specific spectrum, in each case with higher amplitude in low-frequency components, which was similar across the different networks (Fig. 5b, Supplementary Figs. 4 and 9). This reduction of network-specific spectral signatures was highly significant and strongest for the propofol anesthesia regime (Fig. 5c, see also Supplementary Fig. 3).

Figure 5.

Dependence of spectral difference between networks on brain state. Electrodes are assigned to 8 networks based on the ECoG broadband power covariations and displayed using different colors (blue to red from anterior to posterior). The demeaned power spectra are averaged across electrodes belonging to the same network and displayed together for (a) the eyes-closed wakefulness, (b) ketamine/medetomidine anesthesia, propofol anesthesia, and sleep conditions. (c) The upper panel shows the MBGSS lines that quantify cross-network difference in power spectra, and the lower panel presents the ratio of 2 MBGSS lines, that is, the F-lines. The dashed lines are F-values (wakefulness versus ketamine/medetomidine: F140,77 versus propofol: F140,28 and versus sleep: F140,35) corresponding to a significance level of P < 0.01. Shadows represent areas within 1 S.E.M. across experiments (n = 20, 11, 4, and 5 for the wakefulness, ketamine/medetomidine, propofol, and sleep conditions, respectively).

Figure 5.

Dependence of spectral difference between networks on brain state. Electrodes are assigned to 8 networks based on the ECoG broadband power covariations and displayed using different colors (blue to red from anterior to posterior). The demeaned power spectra are averaged across electrodes belonging to the same network and displayed together for (a) the eyes-closed wakefulness, (b) ketamine/medetomidine anesthesia, propofol anesthesia, and sleep conditions. (c) The upper panel shows the MBGSS lines that quantify cross-network difference in power spectra, and the lower panel presents the ratio of 2 MBGSS lines, that is, the F-lines. The dashed lines are F-values (wakefulness versus ketamine/medetomidine: F140,77 versus propofol: F140,28 and versus sleep: F140,35) corresponding to a significance level of P < 0.01. Shadows represent areas within 1 S.E.M. across experiments (n = 20, 11, 4, and 5 for the wakefulness, ketamine/medetomidine, propofol, and sleep conditions, respectively).

Finally, we performed spectral analysis on ECoG BLPs, in addition to those on raw voltage traces (Fig. 5). No any network-specific patterns were found except for the β- and γ-BLPs under ketamine/medetomidine anesthesia (Supplementary Fig. 10). In general, the BLP modulations were dominated by a very low-frequency component with temporal characteristics comparable to those of fMRI signals (Biswal et al. 1995; Fox and Raichle 2007). In contrast to the robustness of maps derived from the BLP covariation, the spectrum of BLP fluctuations differed markedly between the conditions and also depended strongly on the frequency band of the underlying BLP.

Discussion

In the current study, we investigated the nature of spontaneous neuro-electrical activity by adapting an analysis method previously used for the extraction of patterns of correlated activity in fMRI data (Liu et al. 2012). Applied to large-scale ECoG recordings from macaque, slow (seconds-scale), correlated fluctuations were observed in spatial clusters, many of which resembled those found with fMRI. Similar clusters were found when pre-filtering the data in any of the classical spectral frequency bands. In the following, each of these findings will be discussed in more detail.

Functional ECoG Clusters and the Similarity to fMRI Networks

Although the nature of the fMRI and ECoG data (they were derived from different studies) precluded a detailed quantitative comparison, the apparent similarity between the covariational patterns between these different modalities is consistent with previous observations using ECoG recordings over more limited regions of the human cortex (He et al. 2008; Nir et al. 2008; Wang et al. 2012). These observations suggest that spontaneous variations in the ECoG and fMRI may originate from the same or similar neurophysiological phenomenon. Judged from the ECoG data, the timescale of this phenomenon appears to be in the order a few seconds, much slower than the 1–100 millisecond timescale at which much of the cortico-cortical communication occurs. Possible candidates that may occur on this slow timescale are neuromodulatory processes such as fluctuations in local arousal and awareness (Marrocco et al. 1994; Jones 2003) and purely homeostatic process associated with so-called slow waves (Massimini et al. 2004).

Rather unexpectedly, similar spatial clusters were found for each of the spectral frequency bands investigated (Fig. 3, Supplementary Fig. 7), with the exception of a modest increase in long-range connections for the γ-band. This contrasts to previous findings of band-specific correlations in the αβ range (Brookes et al. 2011; Hipp et al. 2012; Wang et al. 2012) or γ range (He et al. 2008; Nir et al. 2008). The discrepancy could have multiple reasons. Poor sensitivity of skull recordings to high-frequency γ-band, which could be further aggravated with typical orthogonalization procedure, may explain the lack of correlations in this particular range when using magnetoencephalography (Brookes et al. 2011; Hipp et al. 2012). Meanwhile, quantitative differences in spatial specificity (Nir et al. 2008) and correlation strength (Wang et al. 2012) risk over-interpretation as exclusive contributions from dominant frequency bands. Owing to the broad coverage of the electrode array and thus the removal of global nonspecific correlation, we were able to focus on the spatial organization of region-specific correlations, which are clearly independent of spectral frequency.

A possible explanation for the apparent broadband nature of neuro-electrical process underlying the observed network activity is the presence of neuromodulatory and homeostatic process mentioned above, each of which can lead to an overall increase of firing rate (Cash et al. 2009; Lewis et al. 2012), independent of spectral frequency. These and other processes with broadband character could have an event-like neurophysiological signature including the well-known phenomena of slow waves and K-complexes, both of which may have a similar cortical origin (Steriade et al. 1993; Cash et al. 2009). Another possibility is the occurrence of so-called avalanches, which have been studied extensively in a range of electrophysiology preparations (Plenz and Thiagarajan 2007). These were found to occur during both waking and anesthesia and to modulate spectral power in a wide frequency range (Hahn et al. 2010; Thiagarajan et al. 2010). In addition, their spatial extent is thought to reflect the underlying neuronal connectivity and cover local as well as remote cortical regions (Plenz and Thiagarajan 2007). Preliminary fMRI evidence for avalanche-like spontaneous activity has been reported recently in human brain (Tagliazucchi et al. 2012; Liu et al. 2013; Liu and Duyn 2013; Petridou et al. 2013), raising the possibility that it may be a common source for the ECoG and fMRI findings in macaque. Dedicated analysis to identify the presence of such activity in the ECoG data will be required to confirm the contribution, if any, of event-like activity (e.g. slow waves, K-complexes, or avalanches) to the spatial patterns found in the ECoG data.

(In)dependence of Correlational Structure on Conscious State

An intriguing finding is the similarity in cluster distribution between the various states of consciousness. This parallels findings from fMRI, whose network structures also appear largely immune to changes in consciousness (Vincent et al. 2007; Horovitz et al. 2008; Boly et al. 2009; Wey et al. 2014). This is difficult to understand in light of the profound influences that sleep and anesthesia have on nearly all other aspects of brain physiology, including the raw ECoG signals that serve as the basis of the covariation analysis. Based on previous electrophysiological findings, there is every reason to believe that sleep and anesthesia disrupt some of the spatial structure in cortical activity. How to explain this apparent contradiction?

One possibility is that state of reduced consciousness during sleep and anesthesia indeed disfacilitates some of the cortico-cortical communication typical of the normal waking condition, but other neuro-electrical phenomena may persist (or get facilitated) and use the same neuro-anatomical substrate. For example, given the close relationship between the brain's anatomical connectivity and functional network architecture (Honey et al. 2009; Adachi et al. 2012), it is possible that different types of brief electrophysiological events are propagated along the same anatomical connections during different states of consciousness, in each case resulting in a broadband pattern of covariation with a spatial distribution reflecting the anatomical connections. Clearly more research is needed in this area.

Spectral Changes with Changes in Consciousness

In contrast to invariant cluster distribution across the various conditions of consciousness, notable spectral changes were observed. During wakefulness, individual clusters exhibited unique spectral shapes, with the most prominent differences in the αβ frequency ranges (Fig. 5). During sleep and anesthesia, these differences largely disappeared, and the different clusters adopting a common spectral shape that depends on specific condition. It is therefore tempting to speculate that one hallmark of the conscious state is a differentiation in the spectral composition of activity in the various functional networks. Our observation that the signature of individual networks most prominently involves the 5–25 Hz frequency range may owe to differing levels of top-down input to each region, which previously has been associated with modulation in that frequency range (von Stein et al. 2000; Buschman and Miller 2007; Klimesch et al. 2007; Palva and Palva 2007). A reduction of such top-down input with reductions in consciousness may explain the observed spectral changes. However, in absence of further electrophysiological and behavioral evidence, these suggested explanations remain entirely speculative.

Implications for fMRI Measurement of “functional connectivity”

The invariance of ECoG clusters and fMRI networks across different brain conditions puts into question their relevance to behavior. However, it should be realized that the current findings are limited to the spatial distribution of clusters, rather than the underlying level of activity. It is possible that quantitative level of synchronization in each network may depend on condition, and there is indeed evidence for this from several fMRI studies (Boly et al. 2008; Horovitz et al. 2009; Zhang and Raichle 2010). Therefore, signal correlation derived from fMRI (or ECoG for that matter) may provide a behaviorally relevant measures of “functional connectivity”. At the same time, an exploration of resting-state fMRI signals in aspects other than correlations is probably needed in order to obtain more state-sensitive information. A candidate measure that may be useful in this regard is the dynamic aspect of inter-areal synchronization, which has been shown to provide extra network information that is not available with simple static correlations (Chang and Glover 2010; Rack-Gomer and Liu 2012; Liu et al. 2013; Liu and Duyn 2013).

Notes

The authors appreciate Dr Picchioni for helpful discussion on sleep electrophysiology and Dr Hutchison for his help on fMRI results presentation. Conflict of Interest: None declared.

References

Y
T
Sporns
O
Watanabe
T
Matsui
T
Miyamoto
K
Miyashita
Y
2012
.
Functional connectivity between anatomically unconnected areas is shaped by collective network-level effects in the macaque cortex
.
Cereb Cortex
.
22
:
1586
1592
.
Biswal
B
Yetkin
FZ
Haughton
VM
Hyde
JS
1995
.
Functional connectivity in the motor cortex of resting human brain using echo-planar MRI
.
Magn Reson Med
.
34
:
537
541
.
Boly
M
Phillips
C
Tshibanda
L
Vanhaudenhuyse
A
Schabus
M
Dang-Vu
TT
Moonen
G
Hustinx
R
Maquet
P
Laureys
S
2008
.
Intrinsic brain activity in altered states of consciousness: how conscious is the default mode of brain function?
.
1129
:
119
129
.
Boly
M
Tshibanda
L
Vanhaudenhuyse
A
Noirhomme
Q
Schnakers
C
Ledoux
D
Boveroux
P
Garweg
C
Lambermont
B
Phillips
C
et al
2009
.
Functional connectivity in the default network during resting state is preserved in a vegetative but not in a brain dead patient
.
Hum Brain Mapp
.
30
:
2393
2400
.
Breshears
JD
Roland
JL
Sharma
M
Gaona
CM
Freudenburg
ZV
Tempelhoff
R
Avidan
MS
Leuthardt
EC
2010
.
Stable and dynamic cortical electrophysiology of induction and emergence with propofol anesthesia
.
.
107
:
21170
21175
.
Brookes
MJ
Woolrich
M
Luckhoo
H
Price
D
Hale
JR
Stephenson
MC
Barnes
GR
Smith
SM
Morris
PG
2011
.
Investigating the electrophysiological basis of resting state networks using magnetoencephalography
.
.
108
:
16783
16788
.
Buschman
TJ
Miller
EK
2007
.
Top-down versus bottom-up control of attention in the prefrontal and posterior parietal cortices
.
Science
.
315
:
1860
1862
.
Cash
SS
Halgren
E
Dehghani
N
Rossetti
AO
Thesen
T
Wang
C
Devinsky
O
Kuzniecky
R
Doyle
W
JR
et al
2009
.
The human K-complex represents an isolated cortical down-state
.
Science
.
324
:
1084
1087
.
Chang
C
Glover
GH
2010
.
Time-frequency dynamics of resting-state brain connectivity measured with fMRI
.
Neuroimage
.
50
:
81
98
.
Forgy
EW
1965
.
Cluster analysis of multivariate data: efficiency versus interpretability of classifications
.
Biometrics
.
21
:
768
769
.
Fox
MD
Raichle
ME
2007
.
Spontaneous fluctuations in brain activity observed with functional magnetic resonance imaging
.
Nat Rev Neurosci
.
8
:
700
711
.
Hahn
G
Petermann
T
Havenith
MN
Yu
S
Singer
W
Plenz
D
Nikolic
D
2010
.
Neuronal avalanches in spontaneous activity in vivo
.
J Neurophysiol
.
104
:
3312
3322
.
He
BJ
Snyder
AZ
Zempel
JM
Smyth
MD
Raichle
ME
2008
.
Electrophysiological correlates of the brain's intrinsic large-scale functional architecture
.
.
105
:
16039
16044
.
Hipp
JF
Hawellek
DJ
Corbetta
M
Siegel
M
Engel
AK
2012
.
Large-scale cortical correlation structure of spontaneous oscillatory activity
.
Nat Neurosci
.
15
:
884
890
.
Honey
CJ
Sporns
O
Cammoun
L
Gigandet
X
Thiran
JP
Meuli
R
Hagmann
P
2009
.
Predicting human resting-state functional connectivity from structural connectivity
.
.
106
:
2035
2040
.
Horovitz
SG
Braun
AR
Carr
WS
Picchioni
D
Balkin
TJ
Fukunaga
M
Duyn
JH
2009
.
Decoupling of the brain's default mode network during deep sleep
.
.
106
:
11376
11381
.
Horovitz
SG
Fukunaga
M
de Zwart
JA
van Gelderen
P
Fulton
SC
Balkin
TJ
Duyn
JH
2008
.
Low frequency BOLD fluctuations during resting wakefulness and light sleep: a simultaneous EEG-fMRI study
.
Hum Brain Mapp
.
29
:
671
682
.
Hutchison
RM
Everling
S
2012
.
Monkey in the middle: why non-human primates are needed to bridge the gap in resting-state investigations
.
Front Neuroanat
.
6
:
29
.
Hutchison
RM
Leung
LS
Mirsattari
SM
Gati
JS
Menon
RS
Everling
S
2011
.
Resting-state networks in the macaque at 7T
.
Neuroimage
.
56
:
1546
1555
.
John
ER
Prichep
LS
Kox
W
Valdes-Sosa
P
Bosch-Bayard
J
Aubert
E
Tom
M
di Michele
F
Gugino
LD
2001
.
Invariant reversible QEEG effects of anesthetics
.
Conscious Cogn
.
10
:
165
183
.
Jones
BE
2003
.
Arousal systems
.
Front Biosci
.
8
:
s438
s451
.
Klimesch
W
Sauseng
P
Hanslmayr
S
2007
.
EEG alpha oscillations: the inhibition-timing hypothesis
.
Brain Res Rev
.
53
:
63
88
.
Leopold
DA
Maier
A
2012
.
Ongoing physiological processes in the cerebral cortex
.
Neuroimage
.
62
:
2190
2200
.
Lewis
LD
Weiner
VS
Mukamel
EA
Donoghue
JA
Eskandar
EN
JR
Anderson
WS
Hochberg
LR
Cash
SS
Brown
EN
et al
2012
.
Rapid fragmentation of neuronal networks at the onset of propofol-induced unconsciousness
.
.
109
:
E3377
E3386
.
Liu
X
Chang
C
Duyn
JH
2013
.
Decomposition of spontaneous brain activity into distinct fMRI co-activation patterns
.
Front Syst Neurosci
.
7
:
101
.
Liu
X
Duyn
JH
2013
.
Time-varying functional network information extracted from brief instances of spontaneous brain activity
.
.
110
:
4392
4397
.
Liu
X
Zhu
XH
Qiu
P
Chen
W
2012
.
A correlation-matrix-based hierarchical clustering method for functional connectivity analysis
.
J Neurosci Methods
.
211
:
94
102
.
Liu
Z
Fukunaga
M
de Zwart
JA
Duyn
JH
2010
.
Large-scale spontaneous fluctuations and correlations in brain electrical activity observed with magnetoencephalography
.
Neuroimage
.
51
:
102
111
.
Logothetis
NK
Pauls
J
Augath
M
Trinath
T
Oeltermann
A
2001
.
Neurophysiological investigation of the basis of the fMRI signal
.
Nature
.
412
:
150
157
.
Manning
JR
Jacobs
J
Fried
I
Kahana
MJ
2009
.
Broadband shifts in local field potential power spectra are correlated with single-neuron spiking in humans
.
J Neurosci
.
29
:
13613
13620
.
Mantini
D
Perrucci
MG
Del Gratta
C
Romani
GL
Corbetta
M
2007
.
Electrophysiological signatures of resting state networks in the human brain
.
.
104
:
13170
13175
.
Marrocco
RT
Witte
EA
Davidson
MC
1994
.
Arousal systems
.
Curr Opin Neurobiol
.
4
:
166
170
.
Massimini
M
Huber
R
Ferrarelli
F
Hill
S
Tononi
G
2004
.
The sleep slow oscillation as a traveling wave
.
J Neurosci
.
24
:
6862
6870
.
Mitra
P
Bokil
H
2008
.
Observed Brain Dynamics
.
New York
:
Oxford University Press
.
Moeller
S
Nallasamy
N
Tsao
DY
Freiwald
WA
2009
.
Functional connectivity of the macaque brain across stimulus and arousal states
.
J Neurosci
.
29
:
5897
5909
.
Nagasaka
Y
Shimoda
K
Fujii
N
2011
.
Multidimensional recording (MDR) and data sharing: an ecological open research and educational platform for neuroscience
.
PLoS One
.
6
:
e22561
.
Nir
Y
Mukamel
R
Dinstein
I
Privman
E
Harel
M
Fisch
L
Gelbard-Sagiv
H
Kipervasser
S
Andelman
F
Neufeld
MY
et al
2008
.
Interhemispheric correlations of slow spontaneous neuronal fluctuations revealed in human sensory cortex
.
Nat Neurosci
.
11
:
1100
1108
.
Nunez
PL
Pilgreen
KL
1991
.
The spline-Laplacian in clinical neurophysiology: a method to improve EEG spatial resolution
.
J Clin Neurophysiol
.
8
:
397
413
.
Palva
S
Palva
JM
2007
.
New vistas for alpha-frequency band oscillations
.
Trends Neurosci
.
30
:
150
158
.
Petridou
N
Gaudes
CC
Dryden
IL
Francis
ST
Gowland
PA
2013
.
Periods of rest in fMRI contain individual spontaneous events which are related to slowly fluctuating spontaneous activity
.
Hum Brain Mapp
.
34
:
1319
1329
.
Plenz
D
Thiagarajan
TC
2007
.
The organizing principles of neuronal avalanches: cell assemblies in the cortex?
Trends Neurosci
.
30
:
101
110
.
Rack-Gomer
AL
Liu
TT
2012
.
Caffeine increases the temporal variability of resting-state BOLD connectivity in the motor cortex
.
Neuroimage
.
59
:
2994
3002
.
Rosanova
M
Casali
A
Bellina
V
Resta
F
Mariotti
M
Massimini
M
2009
.
Natural frequencies of human corticothalamic circuits
.
J Neurosci
.
29
:
7679
7685
.
Rousseeuw
PJ
1987
.
Silhouettes: a graphical aid to the interpretation and validation of cluster analysis
.
Comput Appl Math
.
20
:
53
65
.
Schölvinck
ML
Maier
A
Ye
FQ
Duyn
JH
Leopold
DA
2010
.
Neural basis of global resting-state fMRI activity
.
.
107
:
10238
10243
.
Schwarz
GE
1978
.
Estimating the dimension of a model
.
Ann Stat
.
6
:
461
464
.
Siegel
M
Donner
TH
Engel
AK
2012
.
Spectral fingerprints of large-scale neuronal interactions
.
Nat Rev Neurosci
.
13
:
121
134
.
M
Nunez
A
Amzica
F
1993
.
A novel slow (<1 Hz) oscillation of neocortical neurons in vivo: depolarizing and hyperpolarizing components
.
J Neurosci
.
13
:
3252
3265
.
Tagliazucchi
E
Balenzuela
P
Fraiman
D
Chialvo
DR
2012
.
Criticality in large-scale brain FMRI dynamics unveiled by a novel point process analysis
.
Front Physiol
.
3
:
15
.
Thiagarajan
TC
Lebedev
MA
Nicolelis
MA
Plenz
D
2010
.
Coherence potentials: loss-less, all-or-none network events in the cortex
.
PLoS Biol
.
8
:
e1000278
.
Thomson
DJ
1982
.
Spectrum estimation and harmonic analysis
.
ProcIEEE
.
70
:
1055
1096
.
Vincent
JL
Patel
GH
Fox
MD
Snyder
AZ
Baker
JT
Van Essen
DC
Zempel
JM
Snyder
LH
Corbetta
M
Raichle
ME
2007
.
Intrinsic functional architecture in the anaesthetized monkey brain
.
Nature
.
447
:
83
86
.
von Stein
A
Chiang
C
Konig
P
2000
.
Top-down processing mediated by interareal synchronization
.
.
97
:
14748
14753
.
Wang
L
Saalmann
YB
Pinsk
MA
Arcaro
MJ
Kastner
S
2012
.
Electrophysiological low-frequency coherence and cross-frequency coupling contribute to BOLD connectivity
.
Neuron
.
76
:
1010
1020
.
Wey
HY
Phillips
KA
McKay
DR
Laird
AR
Kochunov
P
Davis
MD
Glahn
DC
Duong
TQ
Fox
PT
2014
.
Multi-region hemispheric specialization differentiates human from nonhuman primate brain function
.
Brain Struct Funct
.
219
:
2187
2194
.
Zhang
D
Raichle
ME
2010
.
Disease and the brain's dark energy
.
Nat Rev Neurol
.
6
:
15
28
.