{'fraction_non_alphanumeric': 0.031962129830097975, 'fraction_numerical': 0.01700855014494881, 'mean_word_length': 5.217545060460871, 'pattern_counts': {'":': 1, '<': 8, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Intrinsic Functional Connectivity Resembles Cortical Architecture at Various Levels of Isoflurane Anesthesia

Abstract Cortical single neuron activity and local field potential patterns change at different depths of general anesthesia. Here, we investigate the associated network level changes of functional connectivity. We recorded ongoing electrocorticographic (ECoG) activity from temporo-parieto-occipital cortex of 6 ferrets at various levels of isoflurane/nitrous oxide anesthesia and determined functional connectivity by computing amplitude envelope correlations. Through hierarchical clustering, we derived typical connectivity patterns corresponding to light, intermediate and deep anesthesia. Generally, amplitude correlation strength increased strongly with depth of anesthesia across all cortical areas and frequency bands. This was accompanied, at the deepest level, by the emergence of burst-suppression activity in the ECoG signal and a change of the spectrum of the amplitude envelope. Normalization of functional connectivity to the distribution of correlation coefficients showed that the topographical patterns remained similar across depths of anesthesia, reflecting the functional association of the underlying cortical areas. Thus, while strength and temporal properties of amplitude co-modulation vary depending on the activity of local neural circuits, their network-level interaction pattern is presumably most strongly determined by the underlying structural connectivity.


Introduction
In the absence of specific tasks or stimuli the brain's spontaneous, ongoing activity is structured on different spatiotemporal scales. Various measures of functional connectivity reflecting different modes of intrinsically generated coupling between different brain areas have been developed (Deco et al. 2011;Siegel et al. 2012;Engel et al. 2013). Among these intrinsic coupling modes (ICMs), amplitude envelope correlations (in the following termed "envelope ICMs") represent the level of simultaneity of activity fluctuations between neural populations in 2 cortical areas at a particular carrier frequency (Engel et al. 2013). Amplitude envelopes are typically modulated at frequencies below 0.1 Hz and show low dependence on behavioral state (Leopold et al. 2003). In the awake state, high envelope ICMs are found between bilateral homologous areas (Nir et al. 2008;Hipp et al. 2012) as well as within functional areas of one hemisphere (He et al. 2008). The spatial specificity of the observed patterns can be increased by excluding effects of volume conduction through orthogonalization (see Methods; Hipp et al. 2012). Convolution of the amplitude envelope with the hemodynamic response function yields a close estimate of the functional magnetic resonance imaging (fMRI) blood-oxygen level dependent (BOLD) signal Martin 2014). Accordingly, envelope ICMs across 2 brain areas observed with electrophysiological methods are closely linked to the correlations of their BOLD signals (He et al. 2008;Schölvinck et al. 2010;Keller et al. 2013). Thus, for instance, networks defined by envelope ICMs recorded with magnetoencephalography (MEG) show a high similarity to fMRI-resting networks (de Pasquale et al. 2010;Hipp et al. 2012).
Ongoing activity is severely altered during general anesthesia, that is, the reversible loss of consciousness induced through the administration of hypnotic agents. Patterns of cortical activity under general anesthesia change with the level of the hypnotic agent. After loss of consciousness, the EEG first shows continuous patterns with increased delta-and alphaand decreased beta power (Brown et al. 2010). In deeper anesthesia, brain activity changes to burst suppression, which is characterized by an isoelectric line interrupted by bursts of activity (Llinás and Steriade 2006;Brown et al. 2010). The frequency of bursts decreases with deepening anesthesia, until complete suppression is reached (Brown et al. 2010). Networklevel changes caused by anesthesia have been examined using various recording techniques and analysis approaches. Studies focusing on the awake-anesthetized transition have suggested that loss of consciousness is associated with global changes in functional connectivity (Cimenser et al. 2011) and distinct alterations of phase-amplitude coupling (Lewis et al. 2012;Mukamel et al. 2014). Large-scale envelope ICMs under different anesthesia depths have been studied using fMRI. While some of these studies observe an increase of connectivity with deepening anesthesia (Liu et al. 2013), others report a breakdown of coupling (Bettinardi et al. 2015). Micro-scale networks on the cellular level have been characterized using extra-and intracellular microelectrode recordings under different concentrations of anesthetic agents (Llinás and Steriade 2006). Nevertheless, little is known about the corresponding connectivity patterns in mesoscale cortical networks consisting of spatially and functionally related cortical areas.
Here, we characterize envelope ICMs in visual, auditory and parietal areas of the ferret at different levels of isoflurane anesthesia and a constant nitrous oxide administration. We recorded local field potentials from the cortical surface using a custom-built 64-contact electrocorticographic (ECoG) array. Unlike EEG, these recordings are not affected by the distortion of the electrical field in space and frequency by the surrounding tissues (Buzsáki et al. 2012). The isoflurane concentration was varied between 0.4% and 1.6% in steps of 0.2%. For each anesthesia level, envelope ICMs were quantified in different frequency bands using orthogonalized amplitude correlation (Hipp et al. 2012). Similar coupling patterns were visible at different isoflurane concentrations across animals, reflecting inter-individual differences in drug sensitivity. We therefore applied hierarchical clustering to all connectivity matrices of each frequency band and identified similar coupling-defined network states across animals. The resulting 3 clusters corresponded to light, medium, and deep anesthesia. Overall, the strength of envelope ICMs increased across anesthesia depths.
Normalized connectivity patterns, however, remained similar across anesthesia depths, resembling in their topology known patterns of structural connectivity in the cortex.

Methods
Data presented in this study were collected from 6 adult female ferrets (Mustela putorius). All experiments were approved by the independent Hamburg state authority for animal welfare (BUG-Hamburg) and were performed in accordance with the guidelines of the German Animal Protection Law.

ECoG-Array
We employed a micro-ECoG-array co-developed with the University of Freiburg (ECoG-array; IMTEK, Freiburg) covering a large portion of the posterior, parietal, and temporal surface of the left ferret brain hemisphere (Rubehn et al. 2009). The probe was designed to record from posterior early and higher visual areas, auditory, somatosensory, and parietal areas, located in the medial, suprasylvian, and lateral gyrus, respectively (Manger et al. 2002;Bizley et al. 2007). Sixty-four platinum contacts (ø: 250 μm; 0.05 mm 2 , 7-25 kΩ @ 1 kHz) were arranged equidistantly (1.5 mm) in a hexagonal manner on a 10 μm thin polyimide foil. The approximately 12 × 10.5 mm sized foil was subdivided in 3 "fingers", allowing for a smooth adaptation to the slightly curved surface of the ferret cortex (Fig. 1A).

Surgery
Animals were initially anesthetized with an injection of ketamine (15 mg/kg) and medetomidine (0.08 mg/kg). After tracheotomy, a Y-shaped glass tube was placed in the trachea to allow artificial ventilation and the administration of isoflurane (0.4-1.6%, 1:2 N 2 O-O 2 mix; see in more detail below). N 2 O was added for its anesthetic and analgetic effect. Before fixing the animal's head in the sterotaxic frame, a local anesthetic (Lidocain, 10%) was The ECoG-array was designed to cover the temporo-parieto-occipital cortex. (C) During the course of the experiment, isoflurane concentration was varied in 13 steps across 7 levels. Each level was kept for 20 min, of which the last 4 min was chosen for analysis. Data acquisition started at least 5 h after the initial anesthesia with ketamine/meditomidine. (D) Average power (normalized to the maximum per frequency) during the course of the experiment for all 6 animals.
While frequencies below around 30 Hz showed a relative increase with higher isoflurane concentrations, higher frequencies displayed the opposite effect. applied to the ear cannula. To prevent dehydration throughout the experiment, a cannula was inserted into the femoral vein to deliver a continuous infusion of 0.9% NaCl, 0.5% NaHCO 3 , and pancuronium bromide (6 μg/kg/h). Body temperature was monitored rectally and automatically maintained at 38°C with a homeothermic blanket. After the initial sedation, 1.0% isoflurane was added to the ventilated air to maintain anesthesia throughout the following surgical procedures. Stereotaxic ear-and eyebars were removed after surgery before the level of anesthesia was decreased. An anesthetic monitoring device (Datex-Ohmeda S/5, Datex-Ohmeda Medical Instrumentation Inc., Tewksbury, USA) was used to continuously measure inspiratory and expiratory concentrations of isoflurane (Supplementary Table 2), MAC, N 2 O, O 2 , and CO 2 throughout the course of the experiment. All data were digitized and logged to disk via the electrophysiological recording systems.
To get access to the left parieto-occipital cortex, first, the temporalis muscle was reflected laterally. A saline-cooled and soft-tissue-safe micro-saw (Mectron, Piezosurgery) was used to create a craniotomy (~13 × 11 mm) over the left posterior cortex. After careful removal of the dura, the ECoG-array was placed on the cortex. The ECoG-array was covered with artificial dura (Lyoplant/Braun) and the previously extracted bone flap was returned into place and fixed with KWIK-CAST (World Precision Instruments, WPI), a tissue-safe silicon-based twocomponent glue. Finally, the temporal muscle and the skin were flapped back in anatomical position and provisionally fixed with a suture. This ensured that the cortex would reside in physiological cortico-spinal fluid and avoided cooling of the cortex to nonphysiological temperatures throughout the course of the experiment (Kalmbach and Waters 2012).
After the end of the whole experimental session (~24-36 h), the animal was deeply anesthetized with 5% isoflurane and sacrificed.

Data Acquisition
The continuous stream of neuronal data (~4.5 h) was processed with an AlphaLab "SnR" electrophysiological recording system. Signals were amplified 200-fold, hardware bandpass filtered (0.5 Hz; 2-pole-8 kHz; 3-pole) and digitized at 44.6 kHz sampling rate. It was then immediately software lowpass filtered (347 Hz), downsampled by a factor of 32 (1.4 kHz) and saved to disk. Concurrently, vital signs and anesthesia parameters as described above were sampled with the same system (2.8 kHz) and also stored to disk.
Data acquisition started between 5 h and 9.5 h after the initial ketamine anesthesia following at least 4 h of isoflurane ventilation. Before the recordings started, the isoflurane level was kept at 0.4% for at least 1 h. During the recording period the isoflurane concentration was changed every 20 min in steps of 0.2%, first increasing up to a level of 1.6% and then decreasing again down to 0.4%. This resulted in 13 different conditions of "increasing" and "decreasing" levels of isoflurane. After each dosage change, an equilibration of isoflurane concentration between the inspiratory and expiratory airflow was reached after about 2-4 min. The administered isoflurane concentrations differed on average less than 0.001% from the preset values (STD: 0.018%; see Supplementary Table 2).

Data Analysis
All data analysis was performed with Matlab (The MathWorks Inc.), partially using the Fieldtrip (Oostenveld et al. 2011) and the Chronux toolbox (http://chronux.org/; Mitra and Bokil 2008). Offline data were bandpass filtered in a phase preserving manner between 2 Hz and 300 Hz and further downsampled to 700 Hz. Notch filters were applied at 50 Hz, 100 Hz, and 150 Hz (1 Hz width, 4-pole butterworth). Broken contacts, as indicated by low signal amplitudes or high-amplitude noise, were excluded from further analysis. All data were visually inspected before further analysis to exclude ECG-or other electrical artifacts.
We computed the power spectra with a multitaper (4/1) analysis of all channels and conditions between 2 Hz and 256 Hz at 6 frequency bands per octave. Power in the classical EEG frequency bands-theta (4-8 Hz), alpha (8-13 Hz), beta (13-30 Hz), low gamma (30-64 Hz), and high gamma (64-128 Hz) -was calculated as the average of the respective frequency ranges.
Envelope ICMs were quantified by computing orthogonalized amplitude correlations. As shown previously, these are not affected by zero-phase lagged components of the channels which might be a result of volume conduction or other, external common noise (for a detailed description of the method, see Hipp et al. 2012). Based on a time-frequency representation of the 2 signals, orthogonalization of signal B to A for time t and frequency f was calculated as Orthogonalized amplitude correlations were then calculated as The underlying time-frequency representation was calculated based on Morlet wavelets. The temporal resolution was set to 1/2 of the maximum frequency's cycle duration (7.8 ms), while the initial frequency resolution was at 8 logarithmically spaced bands per octave, resulting in 56 Morlet wavelets of 7 cycles. For further analysis, we calculated the mean correlation coefficients of each channel pair in the theta, alpha, beta, lowgamma, and high-gamma band as specified above.
We analyzed the spectral properties of the amplitude modulation by calculating a multitaper spectrum of the amplitude envelope of each channel. To estimate the amplitude correlation in a given modulation frequency band, we orthogonalized one signal with respect to the other as described above and applied a bandpass filter to both amplitude time courses before calculating the correlation coefficient.

Assignment of Contacts to Cortical Areas
To account for inter-animal differences in the position of the ECoG-array and thus make cross-animal comparisons possible, area-level correlation matrices were calculated for crossanimal analysis. Using the sulci as landmarks, a schematic representation of the array position for each experiment was generated relative to a functional map of the ferret cortex (Bizley et al. 2007). We then assigned each contact to a cortical area by a "winner-gets-all" method and calculated the average connectivity values between areas. See Supplementary Fig. 1 and Supplementary Table 1 for more precise description of the method.

Hierarchical Clustering
For each frequency band, we applied agglomerative hierarchical clustering as implemented in Matlab to the area-level coupling matrices of all animals and conditions to identify common network states. To ensure the robustness of the result, we compared the emerging clusters for 2 different distance measures ("pdist": "euclidian" and "cityblock") and using 2 linkages ("ward" and "average"; Ward 1963). An example of a clustering result is shown in Supplementary Fig. 2.
We chose to partition the dataset of 78 conditions (6 subjects × 13 isoflurane levels) into 3 clusters of 44, 19, and 15 members, respectively. At this level, the consistency of the defined clusters across algorithms was high, and the dendrograms showed a clear separation of the clusters from the next level. Further branching would have led to stronger, inconsistent parcellation with groups of very few members (<6), such that the comparability would suffer from statistically noisy results. The results presented here are based on "average" clustering by "euclidean distance".

Normalization of Correlation Values
To examine correlation differences which could not be explained by the distance between contacts, each correlation coefficient in the contact-level connectivity matrix was re-expressed relative to the correlations of all contact pairs with the same intercontact distance. Contact distances were rounded to the first decimal place. We only considered distances which were represented by at least 25 contact pairs in all animals. All distances above the first one not fulfilling this criterion were excluded. This resulted in 18 distances between 1.5 mm and 9.8 mm ( Supplementary Fig. 3). The z-score z = (x − μ)/σ was calculated for each of the correlation values.

Burst Detection and Avalanche Definition
We performed an avalanche-style analysis of epochs of burst suppression (Beggs and Plenz 2003;Shew et al. 2009). The timedomain data of the individual channels were segmented in periods of "burst" and "suppression" using a modified version of the method proposed by Lewis et al. (2013). We analyzed data between 1.0% and 1.6% isoflurane, since under these concentrations all animals displayed burst-suppression patterns. All data were z-transformed. For each animal, periods of unambiguous suppression were labeled in the first 120 s of the recording from one contact under 1.0% isoflurane. Data were filtered above 8 Hz and the instantaneous amplitude was calculated via a Hilbert transform and smoothed with a moving average filter of a width of 200 ms. The threshold for bursts epochs was set 4 standard deviations above the mean of manually labeled suppression epochs. Suprathreshold epochs of least 200 ms were defined as bursts. Suppression phases were defined as epochs in which the value fell below the threshold for at least 200 ms. The result was cross-checked by visual inspection. Avalanches of neural activity were defined as epochs of burst activity propagating across the grid. They were by definition separated by periods of global suppression, in which none of the channels was in "burst mode". Channels which displayed bursting activity for more than 95% of the recording time were discarded before avalanche definition.

Results
We recorded ongoing cortical activity from the left parietotemporo-occipital cortex of 6 anesthetized female ferrets using a custom-made 64-channel subdural electrode array (Fig. 1A). Animals were ventilated with a mixture of isoflurane and N 2 O. We determined the changes of the brain network connectivity at different isoflurane levels.

Isoflurane Increases Low-Frequency and Decreases High-Frequency Power
We examined the effect of deepening anesthesia on power change relative to the first data segment in the typical frequency bands. For this purpose, we calculated the modulation index (MI) (A − B)/(A + B) within each contact-with "A" being the power of the segment of interest, "B" of the first data segment-scaling the values to a range from −1 to 1. With increasing isoflurane concentration all animals showed a relative increase in theta, alpha, and beta power and a relative decrease in gamma power (Fig. 1D).
We next assessed whether there were any topographical differences in the power change occurring from light to deep anesthesia. We pooled the power values across 5 functional groups as follows (cf. Supplementary Fig. 1; Bizley et al. 2007): early visual (areas 17, 18, 19), higher visual (area 21, SSY), auditory (areas A1, AAF, ADF, PPF, PSF), posterior parietal (areas PPr, PPc), and somatosensory (areas SII, SIII). For these areal groups, we compared the MI at 1.6% isoflurane. For all frequency bands, the power change relative to 0.4% isoflurane was significant for all functional groups (Komolgorov-Smirnov test (KS-test), Bonferroni-Holm corrected, P < 0.05). However, there were only a few significant differences between the functional groups (Fig. 2). Overall, the Figure 2. Modulation index (MI) of the power within different groups of cortical areas at 1.6% isoflurane relative to the starting condition at 0.4% isoflurane. In all groups, the power showed a significant difference to the starting condition (MI different from 0, KS-test, P < 0.05, Bonferroni-Holm-corrected). While power in all groups increased in the theta-, alpha-, and beta band, it decreased in the gamma-and high-gamma band. Abbreviations: EVis, early visual, HVis, higher visual, PP, posterior parietal, SoSe, somatosensory, Aud, auditory. observed relative change in power occurred relatively homogeneously across the whole temporo-parieto-occipital cortex. This result is in accordance with other studies showing an occipitofrontal power shift under anesthesia rather than local differences between the areas of the occipital cortex (Cimenser et al. 2011;Supp et al. 2011;Sellers et al. 2013).

Cluster Analysis of Envelope ICMs Reveals 3 Anesthesia-Dependent Network States
For each 4-min data segment, we calculated orthogonalized amplitude correlations between all channels. This measure excludes the zero-phase-lag component of the 2 signals and is thus not susceptible to potential effects of volume conduction (Hipp et al. 2012). Correlation values were calculated on a wavelet-transformed time-frequency representation and then pooled to the following frequency bands: theta (4-8 Hz), alpha (8-13 Hz), beta (13-30 Hz), low gamma (30-64 Hz), and high gamma (64-128 Hz).
For these frequency bands, we then computed coupling matrices reflecting the envelope ICMs at different anesthesia depths. It has long been known from clinical settings that different anesthetic concentrations are needed to elicit the same anesthetic effect in different individuals. We hypothesized that due to inter-individual differences in responsiveness to the anesthetic agent, similar coupling patterns might emerge in the individual animals under different concentrations of isoflurane. Therefore, we used hierarchical clustering to identify similar coupling patterns for each frequency band across animals and anesthesia depths.
For each frequency band, we calculated the envelope ICM coupling matrix for each data segment and animal. Amplitude correlation values were averaged within 14 cortical areas covered by the ECoG-array in all animals. This approach eliminated inter-individual differences of the contact positions and the number of recording sites per cortical area, allowing crossanimal comparison. For each frequency band, the matrices were vectorized and subjected to hierarchical clustering ( Supplementary Fig. 2). Based on the dendrograms, we chose to partition the dataset into 3 clusters corresponding to network states defined by functional connectivity under light, medium, and deep anesthesia (Fig. 3A). The deep-anesthesia connectivity state was consistent across all frequency bands but appeared in only 4 of the 6 animals. In contrast, the light and medium anesthesia connectivity pattern were observed in all animals. They were nearly identical for the theta, alpha, and beta bands. In the low-and high-gamma band, the medium anesthesia connectivity state occurred at lower isoflurane concentrations than in the lower frequency bands. This indicates a frequency band related difference in network-level effects of isoflurane.
For each anesthesia network state, we calculated the arealevel coupling matrices averaged across animals (Fig. 3B). Within each frequency band, global correlation strength increased from the light to the deep-anesthesia network state for all fields of the coupling matrix (light: 0.08, medium: 0.23, deep: 0.41; HedgesG: light-medium: −1.57, medium-deep: −1.66, light-deep: −3.25). A KS-test between the clusters confirmed that this difference was significant for all connections in the matrix (Bonferroni-Holm corrected, P < 0.05).
The structure of the coupling matrices suggested that particularly high correlations occurred among functionally closely related areas. To test this hypothesis, we pooled the singlecontact correlation values into 5 functional groups (early visual, higher visual, posterior parietal, auditory, and somatosensory) and compared correlation values within these functional groups with those across functional groups (Fig. 3C). Within functional group correlations were consistently higher than across functional group correlations for all network states and frequency bands (KS-test, Bonferroni-Holm corrected, P < 0.05).

Network States are Associated with Changes of Burst-Suppression Characteristics
Inspection of the raw ECoG signal for the 3 connectivity states revealed patterns known from single-cell and EEG recordings (Steriade et al. 1994). The data segments in the light anesthesia state mostly revealed continuous activity patterns, which in some animals were overlaid by bursts of higher-frequency activity. For the medium anesthesia state frequent bursts with short periods of suppression were observed and for deepanesthesia groups of bursts interrupted by long suppression intervals ( Supplementary Fig. 4). At the lowest isoflurane dose of 0.4%, 3 animals already showed burst-suppression patterns. At the highest concentration of 1.6%, all animals had entered burst suppression, 4 displaying the deepest burst-suppression state. Spectral analysis of the anesthesia states showed an increase of low-frequency power and a decrease of highfrequency power in the medium and deep-anesthesia state (Fig. 1D, see also Supplementary Fig. 5).
We further characterized the burst-suppression segments of the dataset using the approach of Lewis et al. (2013). Data segments recorded under isoflurane concentrations of 1.0% and higher were included since under these conditions all animals displayed burst-suppression patterns.
We calculated the average burst frequency, mean burst length, and mean suppression length for each data segment. Burst frequency and suppression length showed a distribution similar to the medium and deep-anesthesia state as obtained through our connectivity-based analysis ( Supplementary Fig. 6A, B). Burst lengths showed no clear dependence on isoflurane concentration ( Supplementary Fig. 6C).
We quantified the high visual similarity of connectivitybased anesthesia states, burst frequency, and suppression length through linear correlation. For all frequency bands, the matrices were significantly correlated at P < 0.001 (Bonferroni-Holm corrected, Supplementary Table 3). Thus, the anesthesiadependent changes of burst-suppression characteristics of local activity were closely linked to the states observed on the network level. Deep-anesthesia states were characterized by a decrease of burst frequency and an increase of the length of suppression episodes.

Overall Strength of Envelope ICMs Shows Dependence on Distance and Frequency Range
To capture general properties of envelope ICMs within the 3 network states, we examined their dependence on frequency and spatial distance of the contacts. For this purpose, we calculated the mean correlation values across all electrode combinations and within 6 contact distance categories for the 3 clusters (Fig. 4). For all 3 anesthesia-dependent network states, correlation showed a drop-off as a function of distance and frequency. In the light anesthesia network state, envelope ICMs showed a peak in the alpha-frequency range and continuously diminished with increasing frequency. During medium depth anesthesia, the highest correlation values were observed in the beta band. During deep anesthesia, mean strength of envelope ICMs was similar in the theta, alpha, and beta bands, and only dropped within the gamma and high-gamma band.

Spatial Structure of Envelope ICMs Remains Robust to Anesthesia Depth and Distance Effects
As described above, the most obvious effect of anesthesia was an overall increase of envelope ICMs between and within all areas. The relative spatial pattern within the connectivity matrix, however, appeared similar within each frequency band across anesthesia-dependent network states. To further examine this aspect, we re-expressed each value within each area-level coupling matrix as a z-score relative to all correlation values of the particular coupling matrix. The resulting mean z-scored connectivity matrix (Fig. 5A), representing only the spatial structure of the envelope ICMs separated from the absolute coupling  strength, remained similar across anesthesia depths. The topography of the envelope ICMs was similar for all frequency bands.
Since the functional groups comprised neighboring cortical areas, the general decrease of amplitude correlation with distance would induce a bias to higher intra-compared with cross-group correlations. Therefore, we tested whether the spatial patterns visible in the coupling matrices could be explained exclusively by the lower distance of the contacts comprised in areas and functional groups. For this purpose, we re-expressed each correlation coefficient of the contact-level coupling matrix Corresponding within-group ICMs (ig) were higher compared with crossgroup coupling (cg) in all 3 states. (C) Each correlation value was re-expressed as a z-score to the correlations of the contact pairs with a similar inter-contact distance (cf. Supplementary Fig. 3). Connectivity patterns were preserved in the resulting distance-corrected matrices. (D) Coupling within groups of cortical areas was still higher than cross-group coupling. as a z-score relative to all correlations between contact pairs with the same distance ( Supplementary Fig. 3). On the resulting values, we recomputed the coupling matrices for the 3 network states defined above. As shown in Fig. 5C, the topographical patterns observed previously were still largely preserved. Furthermore, within-group z-scores were still higher than across-group z-scores (KS-test, Bonferroni-Holm corrected, P < 0.05). This suggests that higher envelope ICMs between functionally related areas cannot exclusively be attributed to their spatial proximity.

Envelope Modulation Frequency Changes with Anesthesia Depth
To further investigate the temporal dynamics of amplitude envelopes in the 3 network states, we calculated the mean power spectrum of the amplitude envelope for each frequency band and each state (Fig. 6A). The spectrum of the signal amplitude was relatively flat for all frequency bands in the light anesthesia network state. During medium anesthesia, a peak around 0.7 Hz was visible, which was most pronounced in the gamma band. For the deep-anesthesia state, the peak frequencies of the amplitude modulation were below 0.1 Hz.
This change in modulation frequency raised the question whether amplitude correlations differed for specific modulation frequencies, and whether this dependency changed with anesthesia. We therefore bandpass filtered the time course of the amplitude envelope in eleven frequency bands between 0.04 Hz and 2 Hz and calculated the mean correlation within each of them for each frequency and each cluster (Fig. 6B). Correlation at all modulation frequencies increased with anesthesia. During light anesthesia, the resulting spectrum was flat. For the medium anesthesia state, it showed a peak around 0.5 Hz, which was more pronounced in the beta-and gamma-band than in the theta-and alpha-band. During the deep-anesthesia network state, correlation showed a constant decay with modulation frequency.

Deep-Anesthesia Network State is Characterized by Global Spread of Avalanches
To link our analysis of burst-suppression epochs to the global network dynamics we used an avalanche-based approach, established in the analysis of activity of neuronal and slice cultures (Beggs and Plenz 2003;Shew et al. 2009). Avalanches are cascades of activation propagating through a network. Once one node of the network is activated (in our case, when the area of cortical tissue under one ECoG contact switched into "burst mode"), it can activate other nodes of the network, thus initiating a cascade of activity. This cascade can reach from the activation of only one neighboring node to the activation of all nodes of the network. By definition, individual avalanches are separated by epochs of global "suppression" across all contacts of the ECoG.
Data were segmented into burst-and suppression epochs as described above. We included data recorded at isoflurane concentrations of 1.0% and above, under which all the time-series data showed burst-suppression patterns. Avalanches were defined as sequences of noninterrupted, overlapping burst activity in at least one of the channels of the grid. They were by definition separated by epochs of global suppression, in which none of the channels was in "burst-mode". The number of avalanches per data segment of 4 min was 87.2 ± 37 std. As a parameter relatively independent on the small number of avalanches in each data segment, we chose the ratio of avalanche activity time to suppression time to further classify the data. This ratio showed a pattern similar to the connectivity states (Fig. 7A). Hierarchical clustering of the values showed 2 states almost identical to the medium and deep-anesthesia connectivity state (Fig. 7B). Linear correlation between the state matrices was significant for all frequency bands at P < 0.001 (Supplementary Table 3).
Avalanches differed in their spatial extent, ranging from activity in at least one channel without further spread of activity to activation of all channels during one avalanche. We analyzed the distribution of the number of channels active during anesthesia for the different frequency bands. In low anesthesia, the spectrum is relatively flat. In medium anesthesia a peak emerges around 0.5 Hz which is most pronounced in the gamma band, and in deep anesthesia, power shifts to slower frequencies below 0.1 Hz. (B) Envelope ICMs computed at different modulation frequencies of the amplitude envelope. The shaded area represents the values between the 25th and the 75th percentile. Note the correlation peak around 0.5 Hz in medium anesthesia and the increased correlation at slower frequencies in deep anesthesia. avalanches for the 2 anesthesia states. For both medium-and deep-anesthesia state, the shape of the distribution was similar. It had a small peak at bursts including only few channels and a larger peak at global bursts including almost all channels of the grid. In deep anesthesia, the peak for global bursts was steeper and higher than in medium anesthesia (Fig. 7C). Significance was tested using a bootstrap statistic. We randomly shuffled the data segments into clusters of similar size to the original ones 10 000 times. For each randomization, we calculated the area between the channel number distributions of the 2 random clusters. We then compared the original area between the curves (see Fig. 7C) with the resulting distribution of areas from the randomizations; they differed significantly (P = 0.0018; Supplementary Fig. 7C). Distribution of avalanche length showed a decay from smaller to larger values in both medium and deep anesthesia ( Supplementary Fig. 7A) without a significant difference between the states (P = 0.063; Supplementary Fig. 7B).
Overall, avalanche-based analysis of the burst-suppression data segments showed the same anesthesia states as the connectivity-based approach. Under deep-anesthesia avalanches showed a higher tendency to spread globally across the whole cortex, which coincides with the high global connectivity observed in this state.

Discussion
This study shows, using the ferret cortex as a model, that quantification of envelope ICMs allows to characterize network states emerging across different concentrations of isoflurane anesthesia under nitrous oxide/oxygen-ventilation. Based on hierarchical clustering of the coupling matrices, we could differentiate 3 network states under light, medium, and deep anesthesia. The particular anesthetic dosage at which a state was reached varied across animals, reflecting individual differences in the response to the anesthetic. Overall, our results suggest that envelope ICMs in ongoing activity reflect the functional organization of cortical networks in their topographical structure. Although the average coupling increased with anesthesia depth, the topography of functional networks was conserved across different network states. The dominant modulation frequencies of the amplitude envelopes shifted across anesthesia depths, and so did the frequency bands at which the strongest envelope ICM was observed. Anesthesia states classified by envelope ICMs were associated with different dynamics of local neural population activity visible in the electrocorticogram. The light anesthesia state with low global correlations was characterized by continuous activity patterns, the high correlation levels of the deep-anesthesia state were related to low-frequency burst suppression.

Methodological Considerations
It has to be emphasized that activity we observe is not exclusively a result of isoflurane anesthesia, but as well of the coadministered nitrous oxide. The network we examined consists of early as well as higher auditory, visual, and somatosensory areas, thus comprising a relatively large fraction of the cortex. However, since frontal and occipital areas show different reactions to anesthetic agents in humans (Cimenser et al. 2011;Supp et al. 2011) as well as in the ferret (Sellers et al. 2013), our results cannot be generalized to the entire neocortex without further experimental evidence.
We partition the signal into frequency bands that can be differentiated in ferrets during the awake state and sleep (Stitt et al. 2017) and are as well observed in primates, cats, and rodents. This provides a common framework for comparison across states, but at the cost of neglecting possible state-specific boundaries. The approach for quantification of the envelope ICMs used here includes the orthogonalization of the phases at each frequency which eliminates effects of volume conduction and common input, but comes at the cost of also removing additional physiological in-phase activity, which plays an important role in neural processing (Singer 1999;Engel et al. 2001). As shown before, this also implies that our approach to most underestimates functional connectivity (Hipp et al. 2012). To evaluate the effects of orthogonalization on connectivity, we also calculated the envelope ICMs using non-orthogonalized amplitude correlations. The clustering result was similar and the ICMs displayed similar coupling patterns with a slight overall increase of connectivity ( Supplementary Fig. 8). This finding points against the possibility that orthogonalization eliminates crucial physiological characteristics of network activity under anesthesia. Furthermore, it shows that the ECoG signal is only minimally influenced by volume conduction. Our state classification by hierarchical clustering interprets deepening anesthesia as switches between clearly distinguishable states. Due to the relatively large logical distance between the 3 clusters visible in the dendrograms, we believe this approach to the data is justified, even though it ignores minor within-state differences of the coupling patterns. Furthermore, our approach does not take into account spontaneous fluctuations in network connectivity strength over time (

Changes of Local Dynamics and Spectral Characteristics Under Anesthesia
The ECoG signal can predominantly be attributed to the synchronous population activity of cortical pyramidal neurons (Buzsáki et al. 2012). We observed a change of the ECoG signal from continuous activity under light anesthesia to burstsuppression patterns under higher isoflurane concentrations, an effect long known from surface EEG recordings in clinical settings (Brown et al. 2010). Spectral characterization of local activity revealed changes known from previous studies that have employed isoflurane anesthesia (Sellers et al. 2013). With increasing isoflurane concentration we observed a relative increase in theta, alpha, and beta power and a relative decrease in gamma power. These changes in the local dynamics relate to potentiation of GABAergic inhibition in the thalamo-cortical system (Lukatch et al. 2005). At the single-cell level, bursts correspond to epochs of phasic depolarization, while periods of suppression correspond to electrical silence resulting from anesthetic-mediated suppression of glutamatergic transmission (Steriade et al. 1994;Lukatch et al. 2005).
The spectrum of the amplitude envelopes changes across anesthesia depth, likely representing the drastic alteration of activity pattern visible in the surface ECoG signal. The continuous spontaneous activity of the low anesthesia state was reflected in "envelopes", not raw or band-passed activity, predominantly modulated at frequencies below 0.1 Hz, similar to those described for the awake brain (Leopold et al. 2003;Nir et al. 2008). In medium anesthesia, the "envelope spectrum" showed an overall power increase with a peak around 0.6 Hz especially for the alpha-, beta-, and gamma-band, likely a spectral fingerprint of the emerging burst-suppression patterns. The further increase of modulation frequency, especially below 0.1 Hz in deep anesthesia, might be the correlate of the sparser occurrence and longer duration of bursts, however, a more systematic examination of burst properties would be necessary to verify this assumption. The correlation of the envelopes is highest at the dominant modulation frequencies of the particular state.

Envelope ICMs and Cortical Network State
Complementary to the states of local neuronal population activity visible in the patterns of the raw ECoG signal, we defined network states through amplitude envelope correlations. With deepening anesthesia we observed a global increase of coupling strength with constant relative spatial patterns. The spatial specificity of envelope ICMs we observed across anesthesia depths has been described before in the awake state and sleep. Homologous areas display high interhemispheric correlations in envelope ICMs calculated from MEG (Hipp et al. 2012). Within the same hemisphere, amplitude correlations between ECoG contacts within the sensorimotor cortex have been observed to be higher than those to areas of other modalities in the awake state, slow-wave-, and REM-sleep (He et al. 2008). On a finer spatial scale, the resemblance of tonotopy by envelope ICMs has been shown in the macaque (Fukushima et al. 2012). To our knowledge, our study is the first to examine high-resolution envelope ICMs at different depths of isoflurane anesthesia.
In the medium-and deep-anesthesia network states, which were characterized by burst-suppression patterns, cortical activity was not equally coupled across the whole set of recorded areas, but amplitude correlation was still higher for neighboring cortical sites than for distant ones. Correspondingly, synchrony and similarity of bursts have been shown to depend on cortical distance. Bursts occur nearly simultaneously in adjacent cortical areas (Contreras et al. 1997). However, large-scale ECoG recordings in humans have shown that burst suppression is not completely homogeneous across the whole cortex (Lewis et al. 2013). First, not all areas enter and exit burst suppression simultaneously. Second, under burst suppression, not all bursts are global, and the global ones are less synchronized for more distant cortical areas. These results are well compatible with our observation that even under deep burst suppression, the cortex is not homogeneously synchronized, but the spatial and functional structure of the network is preserved.
In the awake brain, amplitude envelopes and their correlations have been suggested to remain stable across behavioral states (Leopold et al. 2003). Our results show that anesthesia modulates some of their basic characteristics, while others remain constant: While mean correlation increases with anesthesia depth, it shows a negative dependence on distance between recording sites in all 3 states. Similarly, a drop-off of envelope ICMs for all frequency bands as a function of distance has been described in the awake rhesus monkey (Leopold et al. 2003). Furthermore, we observed a drop-off of correlation as a function of frequency, which suggests a more local amplitude co-modulation of higher frequencies. In contrast, in the awake brain particularly high coherence values have been reported for amplitude correlations of the gamma band (Leopold et al. 2003). Whether this is due to the different scale of our recordings, which cover a larger area of much smaller ferret brain, caused by the orthogonalization, or is a property of the anesthetized brain remains open to further investigation.

Relation to fMRI Studies on Functional Connectivity Under Anesthesia
To our current understanding, connectivity revealed by amplitude correlation of electrophysiological signals is comparable to that defined by BOLD correlation in resting-state fMRI (Leopold and Maier 2012;Engel et al. 2013). The low-frequency modulation of neural activity represented in the amplitude envelope is tightly linked to the temporal dynamics of local cerebral blood flow and oxygen extraction mirrored in the BOLD signal via the hemodynamic response function (Martin 2014). ECoG amplitude-and BOLD-derived correlations of ongoing activity are highly correlated (He et al. 2008;Schölvinck et al. 2010;Keller et al. 2013), and resting networks defined by independent component analysis (ICA) of MEG power envelopes resemble those known from fMRI studies (de Pasquale et al. 2010(de Pasquale et al. , 2012. The close association of ECoG amplitude envelopes and cerebral blood flow via the hemodynamic response function has been described in isoflurane-anesthetized rats and is preserved in the anesthetized state ).
However, it should also be noted that the comparability of fMRI studies to electrophysiological results is limited by some methodological differences. Most resting-fMRI studies segment activity patterns into distinct networks using spatial ICA or examine networks of few, sometimes distant seed regions. Due to the limited temporal resolution of fMRI, they typically do not consider frequencies above 0.1 Hz, which we show to contribute to amplitude correlations under anesthesia. For the same reason, these studies usually do not examine the temporal dynamics of the signals.
Since different anesthetics potentially differ in their effects on brain connectivity (Brown et al. 2010;Williams et al. 2010;Jonckers et al. 2013;Grandjean et al. 2014), we confine our discussion to studies using isoflurane. It has been shown that the spatial specificity of fMRI-coupling is preserved under isoflurane anesthesia (Vincent et al. 2007;Hutchison et al. 2010). An increase of functional connectivity with deepening anesthesia was observed in a study examining BOLD patterns in the sensorimotor cortex of the rat under 3 isoflurane concentrations (1.0%, 1.5%, 1.8%) (Liu et al. 2013). Distinct sensorimotor subnetworks were distinguishable under the lowest isoflurane concentration and BOLD correlations decreased as a function of distance. With increasing dosage, spatial specificity of correlation patterns decreased and high brain-wide correlations emerged. The same study described a similar effect on power correlations from 3 recording sites, again in accordance with our own findings (Liu et al. 2013).
These strong long-range BOLD correlations observed under deep anesthesia have been shown to be associated with electrophysiological burst suppression and to breakdown under isoflurane concentrations leading to an isoelectric local field potential . High global correlations have been described under 1.5% isoflurane, which is likely comparable to our medium-or deep-anesthesia state (Kalthoff et al. 2013). While there are no previous electrophysiological data on envelope ICMs under isoflurane anesthesia, increased interhemispheric amplitude correlations have been observed during sleep (Nir et al. 2008). Although these findings are well compatible with our results, other studies report a breakdown of functional connectivity with increasing isoflurane concentration (Wang et al. 2011;Hutchison et al. 2014), and similar observations have been made under other anesthetics (Lu et al. 2007;Bettinardi et al. 2015).
Taken together, it is difficult to derive a consistent picture of resting-fMRI networks under anesthesia from the available literature. These contradictory findings might be due to variations in experimental and analytical methodology. Experimental protocols differ in species, regions between which coupling is calculated, isoflurane concentrations, and co-medication. This makes an assignment of similar brain states across studies very difficult. In addition, neurovascular coupling and thus the BOLD signal are altered by blood CO 2 -and O 2 -level (Nasrallah et al. 2015), a mechanism by which different ventilation protocols might influence the observed coupling patterns. Additionally, different fMRI preprocessing strategies seriously influence the observed functional connectivity. Since high correlations in deep anesthesia are based on synchronous activity across large brain areas, they disappear if the global brain signal is removed in the preprocessing of the fMRI data, giving the opposed impression of sparsely correlated networks (Kalthoff et al. 2013;Liu et al. 2013).

Avalanche Analysis of Burst-Suppression
The high connectivity levels observed in the deep-anesthesia state were associated with a higher tendency of avalanches to spread across the whole network. Avalanches of multi-unit spikes have been characterized in recordings of spontaneous activity from neuronal cultures and cerebral slices cultures (Beggs and Plenz 2003;Shew et al. 2009). The authors report that avalanches in these networks follow power-law distributions in the number of active electrodes and duration. Their further findings suggest that the network operates in a critical state with a branching coefficient close to 1, which would maximize the processing capacities of the network (Haldeman and Beggs 2005). Notably, reduction of inhibition by a ɣ-aminobutyric acid (GABA) A -antagonist causes a peak a shift from a power-lawdistribution to a bimodal distribution (Beggs and Plenz 2003;Shew et al. 2009). We observed similar effects under burst-suppression caused by the GABA A -agonist isoflurane. This can likely be explained by the fact that isoflurane also acts on a variety of other ion channels.

Envelope ICMs and Cortical Network Function
At first glance, it might seem counterintuitive that deepening anesthesia leads to increased functional connectivity. However, increased envelope ICMs in deep anesthesia are associated with the rather stereotypic local patterns of burst suppression, which are highly synchronized throughout the network, as already suggested by Liu et al. (2011). Both the more stereotyped local dynamics and the overall increase in functional connectivity are in line with previous hypotheses on the loss of consciousness under deep anesthesia (Alkire et al. 2008;Supp et al. 2011). It has been suggested that anesthesia has 2 effects on the cortex: a breakdown of integration, and the entrainment to stereotypic activity patterns leading to a loss of information (Alkire et al. 2008). Our data agree with this notion and support the following interpretation. While the topographic patterns of the envelope ICMs remain stable, the signals propagating across the network change drastically with the transition into deeper anesthesia. In the awake state and under light anesthesia, dynamics of ongoing activity is complex and external stimuli elicit differentiated responses representing stimulus features. Under burst suppression, specific local responses are replaced by stereotypic bursts similar to those occurring spontaneously (Hudetz and Imas 2007). These evoked bursts spread to areas otherwise unresponsive to stimuli of the particular modality (Hudetz and Imas 2007;Land et al. 2012). Thus, on the one hand, the stereotypic patterns prevailing throughout the network hinder the emergence and propagation of "informative" intrinsic and stimulus-related neural activity. On the other hand, global connectivity increases, possibly leading to a reduction of the functional segregation of areas. Both aspects lead to a reduction in information processing capacity which may account for the loss of consciousness under anesthesia (Alkire et al. 2008).
Our results reinforce the hypothesis that envelope ICMs closely reflect patterns of structural connectivity between the underlying areas (Engel et al. 2013). We show that this effect is robust even to a severe perturbation such as that induced by anesthesia. For resting-fMRI, the high similarity of functional to structural connectivity has been shown by tracer injection (Vincent et al. 2007). Nonetheless, as suggested previously, ongoing activity has a role in shaping processing capacities of functional networks (Engel et al. 2013). Contrary to classical views, ongoing activity is likely to carry information and to be endowed with meaningful spatiotemporal structure, which reflects previous learning and can bias the processing of stimuli (Engel et al. 2001;Deco et al. 2011). While some aspects of this information are carried by envelope ICMs, others are represented in complementary connectivity measures such as phase ICMs (Engel et al. 2013), the s estimator (Carmeli et al. 2005(Carmeli et al. , 2007 or transfer entropy measurements (see Wibral et al. 2013, Wollstadt et al. 2014). These aspects of ongoing connectivity are modulated on faster time scales than those reflected in envelope ICMs and thus may be responsible for adapting stimulus processing to rapidly changing contextual settings. Envelope ICMs seem to represent coherent excitability fluctuations that may lead to coordinated changes in the activation of brain areas. As we have proposed, this might regulate the availability of neuronal populations or regions for participation in an upcoming task and gate the emergence of more specific interactions for processing of stimulus-or response-related signals (Engel et al. 2013). However, these aspects clearly await further investigation.

Supplementary Material
Supplementary material is available at Cerebral Cortex online.

Notes
We thank Dorrit Bystron for excellent technical support. Conflict of Interest: None declared.