- 1Department of Anesthesiology, Huashan Hospital, Fudan University, Shanghai, China
- 2Institute of Neuroscience, State Key Laboratory of Neuroscience, Key Laboratory of Primate Neurobiology, CAS Center for Excellence in Brain Science and Intelligence Technology, Shanghai Institutes for Biological Sciences, Chinese Academy of Sciences, Shanghai, China
- 3University of Chinese Academy of Sciences, Beijing, China
- 4Kunming Institute of Zoology, Chinese Academy of Sciences, Kunming, China
Animal functional magnetic resonance imaging (fMRI) has provided key insights into the physiological mechanisms underlying healthy and diseased brain states. In non-human primates, resting-state fMRI studies are commonly conducted under isoflurane anesthesia, where anesthetic concentration is used to roughly infer anesthesia depth. However, within the recommended isoflurane concentration range (1.00–1.50%), the brain state can switch from moderate anesthesia characterized by stable slow wave (SW) electroencephalogram (EEG) signals to deep anesthesia characterized by burst suppression (BS), which is electrophysiologically distinct from the resting state. To confirm the occurrence rate of BS activity in common setting of animal fMRI study, we conducted simultaneous resting-state EEG and fMRI experiments on 16 monkeys anesthetized using 0.80–1.30% isoflurane, and detected BS activity in two of them. Datasets either featured with BS or SW activity from these two monkeys were analyzed to investigate the intrinsic functional connectivity (FC) patterns during BS. In datasets with BS activity, we observed robust coupling between the BS pattern (the binary alternation between burst and suppression activity in EEG signal) and filtered BOLD signals in most brain areas, which was associated with a non-specific enhancement in whole brain connectivity. After eliminating the BS coupling effect by regressing out the BS pattern, we detected an overall increase in FC with a few decreased connectivity compared to datasets with SW activity. These affected connections were preferentially distributed within orbitofrontal cortex, between orbitofrontal and prefrontal/cingulate/occipital cortex, and between temporal and parietal cortex. Persistence of the default mode network and recovery of thalamocortical connections were also detected under deep anesthesia with BS activity. Taken together, the observed spatially specific alterations in BS activity induced by isoflurane not only highlight the necessity of EEG monitoring and careful data preprocessing in fMRI studies on anesthetized animals, but also advance our understanding of the underlying multi-phased mechanisms of anesthesia.
Introduction
Functional magnetic resonance imaging (fMRI) has become a commonly used noninvasive technique for brain activity research in both animals and humans. Resting-state fMRI (rs-fMRI) measures functional connectivity (FC) between different brain regions in the non-stimulus state and is a valuable technique for exploring complex brain networks in different brain states or disease models (Greicius, 2008; Van Den Heuvel and Hulshoff Pol, 2010). Anesthesia is an effective method to prevent head motion and physiological stress when collecting blood oxygenation level dependent (BOLD) fMRI signals, particularly in animal studies. On the other hand, it is also a critical confounding factor because anesthetic selection and dosage can affect or even inverse the intrinsic BOLD fluctuations readout (Hutchison and Everling, 2012; Grandjean et al., 2014; Lv et al., 2016; Paasonen et al., 2017, 2018; Wu et al., 2017; Bukhari et al., 2018). For instance, isoflurane, the most prevalent inhaled anesthetic for monkey fMRI study, has a significant, concentration-dependent influence on FC networks (Vincent et al., 2007; Nallasamy and Tsao, 2011; Grandjean et al., 2014; Hutchison et al., 2014; Bukhari et al., 2018). A limited range of isoflurane concentration (1.00–1.50%) is therefore recommended to obtain stable FC in anesthetized nonhuman primate resting-state experiments (Hutchison et al., 2014).
However, monkey fMRI studies with simultaneously recorded electroencephalogram (EEG) signals have revealed that relatively deep anesthesia characterized by burst suppression (BS) occasionally occurs at the recommended dosage of isoflurane (end-tidal concentration 1.25–1.50%) (Vincent et al., 2007). BS is a stereotypic electrical activity pattern in the brain that presents as continuous alternation between two states, high-voltage waves (burst) and isoelectric epoch (suppression), and is fundamentally different from the slow wave (SW) activity most commonly observed in the resting state (Steriade et al., 1994). Electrophysiological studies have demonstrated wide synchronization and hyper-excitation across the whole brain during BS activity (Steriade et al., 1994). A robust coupling effect between BS pattern (the binary alternation between burst and suppression activity in EEG signal) and BOLD fluctuations has been revealed in sevoflurane-anesthetized humans (Golkowski et al., 2017). Enhanced BOLD connectivity in the somatosensory network has also been reported in rats during isoflurane-induced BS activity (Liu et al., 2011). Further investigation of the same network showed that connectivity became less spatially specific as anesthesia depth increased to BS status, evidencing functional reorganization of the brain (Liu et al., 2013). These findings indicate multiphasic progression from light to deep anesthesia, suggesting a non-linear relationship between brain activity strength and anesthesia depth (Amzica, 2015).
Here, we aim to explore the signature of the whole-brain functional network connectivity during the isoflurane-induced BS activity in adult macaque monkeys. We are also concerned with the influence of preprocessing strategy on the spatially non-specific increase in FC reported by Liu et al. (2013). We are particularly interested in alterations of the default mode network (DMN) and thalamocortical connections, both of which are sensitive to anesthesia (Alkire et al., 1999; Fiset et al., 1999; Vogt and Laureys, 2005). Despite observed decreases in both DMN (Greicius et al., 2008; Huang et al., 2014) and thalamocortical connectivity (White and Alkire, 2003; Huang et al., 2014) under anesthesia compared to the awake state in humans, it is unknown whether or not these two networks are present during the EEG-defined BS activity in nonhuman primates.
Materials and Methods
Participants
All experimental procedures for nonhuman primate research in this study were approved by the Institutional Animal Care and Use Committee at the Institute of Neuroscience and the Biomedical Research Ethics Committee, Shanghai Institutes for Biological Sciences, Chinese Academy of Sciences, and conformed to National Institutes of Health guidelines for the humane care and use of laboratory animals.
The data included in this work were selected from a series of simultaneous resting-state EEG and fMRI experiments of 16 monkeys (Macaca fascicularis, six males, 10 females, and ages 4–5) anesthetized using a recommended concentration range of isoflurane. After careful scrutiny of the EEG traces, BS activity was detected in two monkeys (monkey O and monkey T) at 0.8–1.3% isoflurane during four separated experiments, out of which three experiments were performed on monkey O with a minimum separation of one and half month. The occurrence of BS was more frequent in monkey O than monkey T (Table 1). Further analysis on FC was based on all the data from monkey O and monkey T, including 20 runs with stable BS activity (BS group) and 8 runs with stable slow wave activity (SW group) from monkey O, 7 and 3 runs, respectively, from monkey T. Each run lasted 7 min with an inter-run interval of approximately 1 min within each experiment.
Animal Preparation
The animal preparation procedure was conducted in a similar manner to our previous work (Wang et al., 2013; Lv et al., 2016). Induction of anesthesia was achieved by intramuscular injection with ketamine (10 mg/kg, Gutian Pharma Co., Ltd., China) before data collection. After intubation, each monkey was ventilated with a mixture of isoflurane (2–2.5%, Lunan Pharma Co., Ltd., China) and pure oxygen via either a standard ventilator (CWE, Inc., Ardmore, PA, United States) in the preparation room or an MRI-compatible ventilator (CWE Inc., Weston, WI, United States) inside the scanner room. The monkey was maintained with intermittent positive-pressure ventilation to ensure a constant respiration rate (25–35 breaths/min). Vital signs including blood oxygenation, electrocardiogram (ECG), rectal temperature (Small Animal Instruments, Inc., Stony Brook, New York), respiration rate and end-tidal CO2 (Smiths Medical ASD Inc., Dublin, Ohio) were continuously monitored. Oxygen saturation was kept at over 95% and body temperature was kept constant using a heated water blanket (Gaymar Industries Inc., Orchard Park, New York). Lactated Ringer’s solution was given with a maximum rate of 10 ml/kg/hour during the anesthesia process (Logothetis et al., 1999).
To ensure adequate EEG data quality, we prepared the animal scalp by shaving and thorough cleaning with abrasive gel and alcohol swabs. A custom EEG cap made of stretchable material was fitted over the scalp and then fastened by a chinstrap. An ECG electrode was attached close to the heart to facilitate off-line removal of cardio-ballistic artifacts. To increase the signal-to-noise ratio, we injected conductive gel and ensured a low impedance (<5 kilo-ohms) at each electrode. The electrode filling holes were subsequently covered with medical tape to prevent the gel from drying out during recording. After setting up the EEG cap, the monkey was restrained within a water blanket in a sphinx-like position with the head protruding and facing forward. The animal’s head was secured using a custom-built MRI-compatible stereotaxic frame after local anesthetic (5% lidocaine cream) was applied around the ears to block peripheral nerve stimulation.
The isoflurane concentration was initially set to 1.2–1.3% within the recommended range in rs-fMRI research on isoflurane-anesthetized monkeys (Vincent et al., 2007; Hutchison et al., 2014), and adjusted by 0.05% increase or decrease to keep the continuously monitored physiological parameters within normal ranges (oxygen saturation: >95%; heart rate: 85–160 beats/min; temperature: 36–38°C; respiration rate: 22–35 breaths/min; end-tidal CO2: 24–30 mmHg) throughout the experiment. Datasets were collected at least 15 min after the adjustment of isoflurane level. Details of isoflurane concentration, BS ratio and physiological parameters for all included datasets were listed in Table 1.
Simultaneous EEG-MRI Data Acquisition
MRI images were acquired at the Institute of Neuroscience on a 3T whole-body scanner (Trio; Siemens Healthcare, Erlangen, Germany) running with an enhanced gradient coil insert (AC88; 80 mT/m maximum gradient strength, 800 mT/m/s maximum slew rate) and a custom-built bird-cage volume coil with 8-channel array receiver. Whole-brain resting-state fMRI data were collected using a gradient echo planar sequence (TR = 2000 ms; TE = 29 ms; flip angle = 77°; slices = 32; matrix = 64 × 64; field of view = 96 mm × 96 mm; 1.5 mm × 1.5 mm in plane resolution; slice thickness = 2.5 mm; GRAPPA factor = 2). For each session, 5–10 runs were acquired and each run consisted of 200 functional volumes. A pair of gradient echo images (echo time: 4.22 and 6.68 ms) with the same orientation and resolution as EPI images were acquired to generate a field map for distortion correction of EPI images. High-resolution T1-weighted anatomical images were acquired using a MPRAGE sequence (TR = 2500 ms; TE = 3.12 ms; inversion time = 1100 ms; flip angle = 9°; acquisition voxel size = 0.5 mm × 0.5 mm × 0.5 mm; 144 sagittal slices). Six whole-brain anatomical volumes were acquired and further averaged for better brain segmentation and 3D cortical reconstruction.
Simultaneous EEG scalp recordings were acquired with BrainVision Recorder software using a BrainAmp MR amplifier and a 28-channel EEG cap customized for macaques (Brain Products GmbH, Gilching, Germany) with sintered Ag/AgCl ring electrodes. EEG signals from 21 active channels were sampled at 5000 Hz with a resolution of 0.5 μV per bit and measuring range of ±16 mV. The sampling clocks of the MR and EEG systems were synchronized using the SyncBox (Brain Products GmbH, Gilching, Germany), thus providing the time of fMRI volume acquisition for later gradient artifact removal in EEG signals.
Data Analysis
MRI Data Preprocessing
Preprocessing of functional MRI images was conducted using the SPM8 toolbox1 and the FMRIB Software Library toolbox (FSL2). The first 10 volumes were discarded before preprocessing. The field map images of each session were applied to compensate for the geometric distortion of EPI images caused by magnetic field inhomogeneity using FSL FUGUE. After slice timing correction and motion correction, the corrected images were normalized to standard space of the monkey F99 atlas3 using an optimum 12-parameter affine transformation and nonlinear deformations, and then resampled to 1.5 mm cubic voxels and spatially smoothed with a 3 mm full-width at half-maximum (FWHM) isotropic Gaussian kernel. Linear drift of the volumes was removed, followed by regression of nuisance covariates (six head motion parameters, white matter and ventricle signals) and temporal filtering (0.01–0.1 Hz). To eliminate the BS coupling effect, we conducted another regression using BS pattern as an additional regressor.
EEG-MRI Coupling Analysis
We first assessed the coupling effect of BS pattern on filtered BOLD fluctuations. For each dataset with BS activity in EEG signal, the suppression episodes were automatically defined by thresholding the normalized amplitude summed across all available EEG channels. All labels of suppression onset or offset were visually inspected by an experienced anesthesiologist. The BS pattern was defined as a binary signal with zero indicating suppression epochs and one indicating non-suppression epochs. The time course of the binary BS pattern was then convoluted with the canonical hemodynamic response function (HRF) implemented in SPM8 as a regressor of interest. The significance of BS coupling effect on each voxel was evaluated via one-sample t-test. To account for the effect of repeated measure from the same subject, binary variables indicating whether the datasets were collected from a specific animal were included as covariates in the SPM model. Family-wise error (FWE) correction was applied to account for multiple comparisons with voxel-wise P < 0.001.
Whole-Brain Functional Connectivity Analysis
For the whole-brain FC analysis, we first parcellated the monkey brain into 82 cortical areas based on the Regional Map template (Kotter and Wanke, 2005; Bezgin et al., 2012) and 12 subcortical areas based on the INIA19 (Rohlfing et al., 2012) (see Supplementary Table S1 for a complete list of anatomical labels). Pearson’s correlation coefficients between the mean time courses of any pair of regions were calculated to represent their FC, resulting in a 94 × 94 connectivity network matrix for each dataset. Statistical significance of group differences were evaluated using linear regression analysis with binary variables of individuals as covariates in the GRETNA toolbox4. The edge-wise threshold of the significance level was set at P = 0.001, and cluster-level correction of P < 0.05 was applied to adjust the multiple comparison using the network-based statistic (Zalesky et al., 2010). Effect sizes measured via Hedges’ g values were calculated to evaluate the extent of group difference. In addition, we assessed the normalized distribution of disrupted edges within and between lobes. The bias caused by the unbalanced number of brain regions within different lobes was adjusted using the standardized residuals (Z scores), defined as the raw residuals (the difference between the observed edge count and expected edge count) divided by the square root of the expected edge count (Sheskin, 2003). The significance of the Z scores was estimated via non-parametric permutation (5,000 times). Specifically, the disrupted edges were randomly assigned across the whole brain network and the standardized residuals were re-calculated to generate the null distribution. The percentage of the assignments with a larger or equal Z score was defined as the P value. Bonferroni correction was applied for the correction of multiple comparisons.
Seed-Based Functional Connectivity Analysis
As the functional connections within the DMN and thalamocortical network are of particular interest during the anesthesia process, we further conducted seed-based functional network analysis. Bilateral posterior cingulate cortex (PCC) and bilateral thalamus defined in Regional Map parcellation were selected as seed regions, respectively. The averaged time course of fMRI signals within the seed region was then treated as a regressor of interest in the generalized regression model in the SPM toolbox. The connectivity strength between each voxel and the seed region was estimated via the beta coefficient in the regression model and statistically tested in each group using a one-sample t-test with covariates indicating individual animals. The group difference in voxel-wise connectivity was tested using a one-way repeated measure ANOVA. FWE correction was applied to account for multiple comparisons with voxel-wise P < 0.001.
Results
Brain-Wide BOLD Fluctuations Are Coupled With BS Activity
Figure 1 shows robust coupling between the BS pattern in EEG signals (Figure 1A) and the spontaneous fluctuations in BOLD signals after temporal filtering (Figure 1B, Z scored) in one typical dataset. Group analysis results indicated that the BOLD fluctuations of most neocortex were positively correlated with the EEG BS pattern (Figure 1C, voxel-wise P < 0.001, FWE correction), including the prefrontal cortex, temporal cortex, parietal cortex, somatosensory cortex, PCC, primary motor and premotor cortex, visural cortex, as well as thalamus.
Figure 1. Coupling between BOLD fluctuations and burst suppression patterns in EEG signals. (A) EEG signals with burst suppression pattern simultaneously recorded during fMRI scanning from the central frontal electrode. The bottom and top panel show 380 and 60 s EEG, respectively. The time onset corresponds to the acquisition of the first MRI volume. (B) Burst suppression pattern with HRF convolution and averaged BOLD signals in gray matter (GM), white matter (WM), and cerebrospinal fluid (CSF). The dashed vertical lines in gray are aligned to the trough of the GM signal. (C) Group statistics (n = 27) shown as voxel-wise t-values of a generalized linear model modeling burst suppression pattern as a regressor of interest were displayed in bspmview (voxel-wise P < 0.001, FWE correction).
Altered Whole-Brain Functional Connectivity During BS Activity Compared to SW Activity
Results of whole-brain comparison based on datasets with or without removal of BS coupling effect were presented in Figure 2 and Supplementary Figure S3, respectively. Group-averaged connectivity networks from datasets with BS activity and stable SW activity are presented in Figure 2A. Group difference in each connection is represented using effect size (Hedges’ g value, bottom-left triangle in Figure 2B). A total of 371 functional connections showed significant differences between BS and SW groups (cluster-level corrected P < 0.00, NBS correction with edge-wise P < 0.001, top-right triangle in Figure 2B), including 317 increased and 54 decreased connections in BS group. As presented in Figure 2C, the BS group showed a substantial increase in cortical connectivity, accompanied by significant decrease mainly in connections between parietal and temporal/prefrontal cortex (PFC), and between insula and occipital/parietal cortex. The normalized distribution of altered connections across different brain lobes is summarized in Figure 2D (bottom-left triangle). The number of disrupted edges between certain lobes is significantly larger than the random distribution (P < 0.05, Bonferroni correction). The majority of within-lobe connections are concentrated in the orbitofrontal cortex (OFC). The majority of across-lobe connections are concentrated between OFC and prefrontal/cingulate/occipital cortex, and between parietal and temporal cortex.
Figure 2. Altered functional connections during burst suppression activity compared to stable slow wave activity after the coupling effect with burst suppression pattern was corrected via regression. (A) Averaged functional connectivity (FC) matrices during burst suppression (BS, bottom-left) and slow wave (SW, top-right) activity. Covariates were regressed out before temporal filtering. (B) Effect sizes (ES, Hedges’ g value) of BS versus SW (bottom-left) and corresponding P values (top-right, P < 0.05, NBS correction with edge-wise P < 0.001). Brain nodes are organized according to the regions/lobes as listed in Supplementary Table S1. (C) Altered functional connections represented with node information. See Supplementary Table S1 for details of brain area abbreviations. (D) Normalized spatial distribution of disrupted connections across the brain (bottom-left) and the corresponding significance (top-right, P < 0.05, Bonferroni correction).
Increased Connectivity in Default Mode Network During BS Activity
The DMN was defined as the voxel-wise FC with bilateral PCC. The statistical results of the connectivity strength in BS and SW groups are presented in Figures 3A,B, respectively (voxel-wise P < 0.001, FWE correction). The BOLD activity in PCC was positively correlated with bilateral parietal cortex, medial and centrolateral prefrontal cortex, superior and ventral temporal cortex, primary motor and dorsolateral premotor cortex, anterior cingulate cortex and visual areas, and anticorrelated with secondary somatosensory and inferior parietal cortex during SW activity (Figure 3B). In comparison, there is an overall increase in DMN connectivity during BS activity (Figure 3A), particularly in the superior temporal cortex, prefrontal cortex, secondary somatosensory cortex and auditory cortex (Figure 3C, voxel-wise P < 0.001, FWE correction).
Figure 3. Altered functional connections across anesthesia levels in default mode network. (A,B) Voxel-wise connectivity with bilateral posterior cingulate cortex (indicated in black circle) in burst suppression (BS) and slow wave (SW) states. Volume metrics were mapped to cortical surface in CARET5. (C) Altered functional connections during BS activity compared to SW activity. All results were corrected for multiple comparisons at a significance level of P < 0.05 via FWE with voxel-wise P < 0.001.
Increased Thalamocortical Connections During BS Activity
Thalamocortical connections were assessed by seed-based functional analysis with the bilateral thalamus as the seed region. Significant connections with the thalamus were observed in the bilateral striatum, medial prefrontal cortex, inferior parietal cortex, superior temporal cortex, auditory cortex, somatosensory cortex, posterior insula, visual cortex, anterior and posterior cingulate cortex during BS activity (Figure 4A, voxel-wise P < 0.001, FWE correction). Few thalamocortical connections to striatum, inferior parietal cortex, superior temporal cortex, retrosplenial and posterior cingulate cortex, somatosensory cortex and anterior visual area were observed during SW activity (Figure 4B, voxel-wise P < 0.001, FWE correction). The group comparison results further confirmed the overall increase during BS activity particularly in bilateral putamen, globus pallidus, right visual area, and left secondary somatosensory cortex (Figure 4B, voxel-wise P < 0.001, FWE correction).
Figure 4. Altered functional connections across anesthesia levels in thalamocortical network. (A,B) Voxel-wise connectivity with bilateral thalamus (indicated in black circle) in burst suppression (BS) and slow wave (SW) states. Volume metrics were mapped to cortical surface in CARET5. (C) Altered functional connections during BS activity compared to SW activity. All results were corrected for multiple comparisons at a significance level of P < 0.05 via FWE with voxel-wise P < 0.001.
Discussion
Robust BS Coupling Effect in the Raw BOLD Fluctuations
In the current study, we show dramatic fluctuations in BOLD signal co-occurring with EEG BS activity in the isoflurane-anesthetized monkey brain. The coupling between simultaneously recorded BOLD and EEG signals can be explained by the highly synchronized excitability of the entire cerebral cortex during BS activity. Convergent research has shown that BS is associated with a hyper-excited cortical state, where bursts can be evoked by sub-threshold stimulus or occur spontaneously (Steriade et al., 1994; Hartikainen et al., 1995; Hudetz and Imas, 2007; Kroeger and Amzica, 2007; Ferron et al., 2009). Electrophysiological studies indicate that neurons in almost the entire cerebral cortex exhibit a similar “on-off” pattern during this state, featuring stereotyped alteration between depolarizing events and electrical silence of the neuronal membrane (Steriade et al., 1994). The hyper-excited burst activity was positively correlated with the cerebral blood flow (CBF) fluctuations, which is the base of BOLD fluctuations, as indicated in a recent isoflurane-anesthetized rat research (Liu et al., 2011). Therefore, the BS coupling effect on BOLD fluctuations is supported by the physiological basis that synchronization of spontaneous high-voltage burst activity across a large population of neurons can result in the synchronization of hemodynamic fluctuations across different brain regions.
Separating Specific and Non-specific Functional Connectivity During BS Activity
Although several lines of evidence suggest a coherent neural activity basis of resting-state FC, whether there is a specific FC signature that underlies the BS anesthesia depth remains unclear. As previously mentioned, neural activity that changes between two electrical phases with enormous disparity can lead to radical changes in BOLD fluctuations at a large time scale, which in turn results in an overwhelming increase in the correlation coefficients between whole-brain BOLD signals during BS activity (Supplementary Figure S3). This non-specific coupling effect of BS pattern is superimposed on the potential specific functional changes associated with anesthesia depth or consciousness loss (Figure 2), and obscures the latter effect.
In removing regular nuisance effect, we demonstrated that regression of WM and CSF signals before temporal filtering is more effective (Supplementary Figure S1). This is consistent with the previous findings in human rs-fMRI studies that nuisance-related variability such as head motion and cardiac artifacts can be better controlled when covariate regression is applied before temporal filtering rather than after filtering in the preprocessing procedure (Hallquist et al., 2013). Further evaluation of the BS coupling effect on preprocessed images indicated that BS related global effect as indicated in Figure 1C cannot be fully eliminated by using either preprocessing procedure (Supplementary Figure S2). Therefore, we added the HRF convoluted BS pattern as an additional regressor in the linear regression model to eliminate the global coupling effect in fMRI datasets with BS activity.
Overall Increase in Functional Connectivity During BS Activity
An overall increase in functional connections during BS activity was observed compared to SW state in both whole-brain analysis (Figure 2) and seed-based analysis (Figures 3, 4). As non-specific increase in the correlation matrix due to BS pattern was well controlled in preprocessing, the remaining alterations across different anesthesia levels may largely reflect specific changes in a subgroup of networks or brain areas. Interestingly, most of the connectivity alterations took place within OFC, and between OFC and other cortical areas. The role of the OFC in BS activity or anesthesia progress merits further investigation. Reduced connectivity strength was observed in a subset of connections between insula and occipital-parietal cortex, suggesting the segregation of the limbic system from the neocortex along with the deepening of anesthesia level. Few alterations were observed within subcortical areas, although a significant increase from subcortical to temporal cortex, parietal cortex, OFC and PFC was detected, which was also verified in the seed-based analysis in thalamus.
Preservation of DMN During BS Activity
The DMN is the most robust intrinsic network of the resting-state brain, which is associated with intrinsic processes such as mind wandering and self-reference (Buckner et al., 2008; Vanhaudenhuyse et al., 2011). Studies of alterations in DMN activity corresponding with different sleep stages indicate a potential correlation between the DMN and consciousness in humans (Raichle et al., 2001). Disconnection between the prefrontal cortex and the PCC has been observed during deep sleep (Horovitz et al., 2009; Samann et al., 2011). Subjects in sedation also exhibit a significant decrease in the default mode activity (Greicius et al., 2008; Martuzzi et al., 2010). Taken together, these findings suggest that the DMN may disappear at a certain stage during loss of consciousness. Nevertheless, preservation of the DMN connectivity under deep anesthesia was recently reported in both anesthetized monkeys and humans (Vincent et al., 2007; Boveroux et al., 2010; Hutchison and Everling, 2012). Our results displayed similar DMN patterns during both BS and SW activity as previous monkey studies (Vincent et al., 2007; Teichert et al., 2010). The robust findings of DMN under deep anesthesia indicate that the DMN may not be tightly associated with consciousness level (Buckner et al., 2008).
Recovery of Thalamocortical Connectivity and the Role of Thalamus During BS Activity
Numerous studies have indicated the influence of anesthesia depth on thalamocortical connections. The presence of thalamocortical connectivity has been reported during light anesthesia where subjects just reached unconsciousness (Martuzzi et al., 2010; Mhuircheartaigh et al., 2010). Recent reports show diminished but detectable thalamocortical connectivity in anesthetized subjects associated with unconsciousness (White and Alkire, 2003; Boveroux et al., 2010; Liang et al., 2012). In the current study, few connections with striatum, cingulate cortex, primary sensory areas, inferior parietal cortex, and superior temporal cortex were observed during SW activity. In contrast, a substantial increase of connections to striatum, somatosensory, and visual cortex were detected during BS activity.
The thalamus is believed to be a crucial region for the synchronized burst activity of neocortex during the BS episodes of EEG (Brown et al., 2010). The first study of the neuronal mechanisms of BS (Steriade et al., 1994) confirmed a close correspondence between neocortical and thalamic activities through intracellular and multisite extracellular simultaneous recordings from cortex, thalamus, and upper brainstem. Despite the overwhelming predominance of cortical electrical silence during EEG suppression epochs, there were still rhythmic oscillations discharged to the cortex from a part of the thalamus. This thalamocortical oscillation was able to facilitate neuron firing or subthreshold depolarizing potentials as well as the revival of EEG activity. However, other studies indicated that cortical hyper-excitability during BS activity is favored by diminished cortical inhibition rather than thalamus induced excitation (Ferron et al., 2009). In our results, high synchronization between thalamic activity and BS pattern (Figure 1C) was revealed. Additionally, after the removal of BS coupling effect, there were significant connections from thalamus to primary sensory networks, which were weaker during SW activity. These enhanced thalamocortical connections might underlie the phenomenon that facilitated global burst activity could be induced by various kinds of subthreshold sensory stimuli under deep anesthesia with BS activity (Kroeger and Amzica, 2007). As for the role of thalamus in BS activity, further studies are desiderated.
Anesthetic Sensitivity and the Importance of Simultaneous EEG Monitoring
For monitoring and defining anesthesia level, most fMRI studies generally employ two kinds of anesthesia depth monitoring indexes. One is behavioral endpoint, such as “loss of responsiveness” for human or “immobility” for animals. The other is end-tidal concentration or minimum alveolar concentration (MAC) for volatile anesthetics, both of which represent the anesthetic concentration level in subjects for inhibiting noxious stimulus (Eger et al., 1965). However, these indexes may not reflect the actual brain state (Antognini and Carstens, 2002; Ryu et al., 2018). The improved comprehension of EEG signals under anesthesia makes it a promising method suitable for more accurate control of anesthesia level (Freye and Levy, 2005; Purdon et al., 2013, 2015; Marchant et al., 2014).
Based on post-hoc analysis of the coherence of rs-fMRI correlation at various anesthesia levels, a limited range of isoflurane concentration (1.00–1.50%) is recommended to obtain stable FC (Hutchison et al., 2014). However, as mentioned previously, monkey fMRI studies with simultaneous EEG monitoring have recorded BS activity at 1.25–1.5% isoflurane (Vincent et al., 2007). In current work, we recorded BS activity in two out of 16 monkeys even at 0.8–1.3% isoflurane. The occurrence of BS activity at a commonly used concentration of isoflurane in fMRI studies suggests individual variance in anesthetic sensitivity at the brain level. Anesthetic sensitivity can be affected by both pharmacokinetic and pharmacodynamics factors and varies across different subjects (Chemali et al., 2015). Even in the same subject, there are many other confounding factors, including physiological status on the experiment day, and sleep quality before the experiment. The etiology and mechanism of anesthetic sensitivity in the brain merits further investigation.
Although BS activity induced by an isoflurane level as low as 0.8% rarely happens, it makes a strong argument for the potential occurrence of BS activity when delivering a widely-used range of concentrations in the field. Simultaneous EEG is necessary to monitor the actual brain state in regular rs-fMRI studies on anesthetized animals. More importantly, BS activity-induced changes in FC is robust. In current work, the occurrence of BS activity was spontaneous, resulting in inconsistent burst suppression ratio (BSR, range 0.56–38.9%). BSR is a commonly used parameter to describe the intensity of BS activity with larger BSR indicating higher probability of BS occurrence and deeper anesthesia level. In current study, robust increase of FC was observed across BS runs with different BSR. Greater enhancement was observed under a stable BS state with frequent and long-lasting suppression episodes than those with occasional occurrence of BS event. Further studies with deliberately controlled BSR are needed to investigate the relationship between BSR and BS effect on FC.
Limitations
There are several limitations worth noting in this research. First, in consideration of the difficulty in controlling the occurrence of BS event and the individual variation in anesthetics sensitivity, present finding of the notably weak thalamocortical connectivity during SW activity could be affected by the relative small sample size. Second, another preliminary analysis on the current datasets suggested that the global synchronization strength during BS activity is potentially associated with BSR, which is not deliberately controlled in current study. Future studies with more subjects and precise control of BSR are needed to probe the biological underpinnings of the BOLD fluctuations at this specific anesthesia level. Third, to explore the influence of preprocessing procedures, we compared two conventional strategies used in monkey fMRI data analysis (Supplementary Figure S1). Other widely debated preprocessing steps such as global signal regression were not included. However, Scholvinck et al. (2010) showed that the often discarded global component of resting-state BOLD fluctuations is tightly coupled with underlying neural activity and may affect connectivity analysis. Fourth, there is a great chance that the global BS coupling effect is not fully excluded by regressing out the HRF convoluted BS pattern, although regression is the most commonly used approach to deal with nuisance variables with global influence. Lastly, pure oxygen was used in current experiments, which is considered as a probable confounding factor leading to more widespread connectivity in rs-fMRI studies (Nasrallah et al., 2015). However, we expect that the potential effect of hyperoxia is comparable in both BS and SW runs, as should impose very limited influence on the group difference in FC. Current results on specific alterations in FC during BS activity compared to SW activity should be interpreted with caution. However, the trend of overall increase in intrinsic FC during BS activity is robust and worth attention from researchers in the field.
Conclusion
In summary, we present evidence showing a dramatic increase in inter-regional connectivity in fMRI network under deep anesthesia with BS activity. The enhanced connectivity can be explained by the coupling effect of BS pattern on whole-brain BOLD signals. This non-specific coupling effect can be well controlled by covariates regression conducted before temporal filtering with the HRF convoluted BS pattern as an additional regressor. After the coupling effect was separated from the connectivity matrix, we detected an overall increase in the functional connections under deep anesthesia with BS activity comparing to a lighter level with SW activity, including connections in DMN and thalamocortical networks. The orbitofrontal lobe are the most affected brain areas during BS activity. The non-linear changes in FC from light to deep anesthesia levels highlight the importance of future investigation on the physiological basis underlying BS activity. This is essential to clarify the mechanisms of anesthesia and coma states with BS events.
Data Availability
The data that support the findings of this study are available on request from the corresponding author (ZW).
Author Contributions
ZZ, ZheW, and YW designed research. ZZ, D-CC, and ZhiW performed research. D-CC and ZhiW analyzed the data. ZZ, D-CC, KZ, and ZheW wrote the manuscript.
Funding
This work was supported by the National Key R&D Program of China (No. 2016YFC1201802 to YW, No. 2017YFC1310400 to ZheW), the Strategic Priority Research Program of Chinese Academy of Sciences (No. XDB32000000), grants from National Natural Science Foundation (81671058 and 81730031 to YW, 81571300 and 31771174 to ZheW), Shanghai Municipal Science and Technology Major Project (No. 2018SHZDZX05) and Natural Science Foundation and Major Basic Research Program of Shanghai (No. 16JC1420100).
Conflict of Interest Statement
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Acknowledgments
We thank Hu Zhang and Qinying Jiang for their assistance in monkey data collection.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fnins.2019.00296/full#supplementary-material
Footnotes
- ^http://www.fil.ion.ucl.ac.uk/spm
- ^http://www.fmrib.ox.ac.uk
- ^http://sumsdb.wustl.edu/sums/macaquemore.do
- ^https://www.nitrc.org/projects/gretna
References
Alkire, M. T., Pomfrett, C. J., Haier, R. J., Gianzero, M. V., Chan, C. M., Jacobsen, B. P., et al. (1999). Functional brain imaging during anesthesia in humans: effects of halothane on global and regional cerebral glucose metabolism. Anesthesiology 90, 701–709. doi: 10.1097/00000542-199903000-00011
Amzica, F. (2015). What does burst suppression really mean? Epilepsy Behav. 49, 234–237. doi: 10.1016/j.yebeh.2015.06.012
Antognini, J. F., and Carstens, E. (2002). In vivo characterization of clinical anaesthesia and its components. Br. J. Anaesth. 89, 156–166. doi: 10.1093/bja/aef156
Bezgin, G., Vakorin, V. A., Van Opstal, A. J., Mcintosh, A. R., and Bakker, R. (2012). Hundreds of brain maps in one atlas: registering coordinate-independent primate neuro-anatomical data to a standard brain. Neuroimage 62, 67–76. doi: 10.1016/j.neuroimage.2012.04.013
Boveroux, P., Vanhaudenhuyse, A., Bruno, M. A., Noirhomme, Q., Lauwick, S., Luxen, A., et al. (2010). Breakdown of within- and between-network resting state functional magnetic resonance imaging connectivity during propofol-induced loss of consciousness. Anesthesiology 113, 1038–1053. doi: 10.1097/ALN.0b013e3181f697f5
Brown, E. N., Lydic, R., and Schiff, N. D. (2010). General anesthesia, sleep, and coma. N. Engl. J. Med. 363, 2638–2650. doi: 10.1056/NEJMra0808281
Buckner, R. L., Andrews-Hanna, J. R., and Schacter, D. L. (2008). The brain’s default network: anatomy, function, and relevance to disease. Ann. N. Y. Acad. Sci. 1124, 1–38. doi: 10.1196/annals.1440.011
Bukhari, Q., Schroeter, A., and Rudin, M. (2018). Increasing isoflurane dose reduces homotopic correlation and functional segregation of brain networks in mice as revealed by resting-state fMRI. Sci. Rep. 8:10591. doi: 10.1038/s41598-018-28766-3
Chemali, J. J., Kenny, J. D., Olutola, O., Taylor, N. E., Kimchi, E. Y., Purdon, P. L., et al. (2015). Ageing delays emergence from general anaesthesia in rats by increasing anaesthetic sensitivity in the brain. Br. J. Anaesth. 115(Suppl. 1), i58–i65. doi: 10.1093/bja/aev112
Eger, E. I. II, Saidman, L. J., and Brandstater, B. (1965). Minimum alveolar anesthetic concentration: a standard of anesthetic potency. Anesthesiology 26, 756–763. doi: 10.1097/00000542-196511000-00010
Ferron, J. F., Kroeger, D., Chever, O., and Amzica, F. (2009). Cortical inhibition during burst suppression induced with isoflurane anesthesia. J. Neurosci. 29, 9850–9860. doi: 10.1523/JNEUROSCI.5176-08.2009
Fiset, P., Paus, T., Daloze, T., Plourde, G., Meuret, P., Bonhomme, V., et al. (1999). Brain mechanisms of propofol-induced loss of consciousness in humans: a positron emission tomographic study. J. Neurosci. 19, 5506–5513. doi: 10.1523/JNEUROSCI.19-13-05506.1999
Freye, E., and Levy, J. V. (2005). Cerebral monitoring in the operating room and the intensive care unit: an introductory for the clinician and a guide for the novice wanting to open a window to the brain. Part I: the electroencephalogram. J. Clin. Monit. Comput. 19, 1–76. doi: 10.1007/s10877-005-0712-z
Golkowski, D., Ranft, A., Kiel, T., Riedl, V., Kohl, P., Rohrer, G., et al. (2017). Coherence of BOLD signal and electrical activity in the human brain during deep sevoflurane anesthesia. Brain Behav. 7:e00679. doi: 10.1002/brb3.679
Grandjean, J., Schroeter, A., Batata, I., and Rudin, M. (2014). Optimization of anesthesia protocol for resting-state fMRI in mice based on differential effects of anesthetics on functional connectivity patterns. Neuroimage 102(Pt 2), 838–847. doi: 10.1016/j.neuroimage.2014.08.043
Greicius, M. (2008). Resting-state functional connectivity in neuropsychiatric disorders. Curr. Opin. Neurol. 21, 424–430. doi: 10.1097/WCO.0b013e328306f2c5
Greicius, M. D., Kiviniemi, V., Tervonen, O., Vainionpaa, V., Alahuhta, S., Reiss, A. L., et al. (2008). Persistent default-mode network connectivity during light sedation. Hum. Brain Mapp. 29, 839–847. doi: 10.1002/hbm.20537
Hallquist, M. N., Hwang, K., and Luna, B. (2013). The nuisance of nuisance regression: Spectral misspecification in a common approach to resting-state fMRI preprocessing reintroduces noise and obscures functional connectivity. Neuroimage 82C, 208–225. doi: 10.1016/j.neuroimage.2013.05.116
Hartikainen, K. M., Rorarius, M., Perakyla, J. J., Laippala, P. J., and Jantti, V. (1995). Cortical reactivity during isoflurane burst-suppression anesthesia. Anesth. Analg. 81, 1223–1228.
Horovitz, S. G., Braun, A. R., Carr, W. S., Picchioni, D., Balkin, T. J., Fukunaga, M., et al. (2009). Decoupling of the brain’s default mode network during deep sleep. Proc. Natl. Acad. Sci. U.S.A. 106, 11376–11381. doi: 10.1073/pnas.0901435106
Huang, Z., Wang, Z., Zhang, J., Dai, R., Wu, J., Li, Y., et al. (2014). Altered temporal variance and neural synchronization of spontaneous brain activity in anesthesia. Hum. Brain Mapp. 35, 5368–5378. doi: 10.1002/hbm.22556
Hudetz, A. G., and Imas, O. A. (2007). Burst activation of the cerebral cortex by flash stimuli during isoflurane anesthesia in rats. Anesthesiology 107, 983–991. doi: 10.1097/01.anes.0000291471.80659.55
Hutchison, R. M., and 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. doi: 10.3389/fnana.2012.00029
Hutchison, R. M., Hutchison, M., Manning, K. Y., Menon, R. S., and Everling, S. (2014). Isoflurane induces dose-dependent alterations in the cortical connectivity profiles and dynamic properties of the brain’s functional architecture. Hum. Brain Mapp. 35, 5754–5775. doi: 10.1002/hbm.22583
Kotter, R., and Wanke, E. (2005). Mapping brains without coordinates. Philos. Trans. R. Soc. Lond. B Biol. Sci. 360, 751–766. doi: 10.1098/rstb.2005.1625
Kroeger, D., and Amzica, F. (2007). Hypersensitivity of the anesthesia-induced comatose brain. J. Neurosci. 27, 10597–10607. doi: 10.1523/JNEUROSCI.3440-07.2007
Liang, Z., King, J., and Zhang, N. (2012). Intrinsic organization of the anesthetized brain. J. Neurosci. 32, 10183–10191. doi: 10.1523/JNEUROSCI.1020-12.2012
Liu, X., Zhu, X. H., Zhang, Y., and Chen, W. (2011). Neural origin of spontaneous hemodynamic fluctuations in rats under burst-suppression anesthesia condition. Cereb. Cortex 21, 374–384. doi: 10.1093/cercor/bhq105
Liu, X., Zhu, X. H., Zhang, Y., and Chen, W. (2013). The change of functional connectivity specificity in rats under various anesthesia levels and its neural origin. Brain Topogr. 26, 363–377. doi: 10.1007/s10548-012-0267-5
Logothetis, N. K., Guggenberger, H., Peled, S., and Pauls, J. (1999). Functional imaging of the monkey brain. Nat. Neurosci. 2, 555–562. doi: 10.1038/9210
Lv, Q., Yang, L., Li, G., Wang, Z., Shen, Z., Yu, W., et al. (2016). Large-scale persistent network reconfiguration induced by ketamine in anesthetized monkeys: relevance to mood disorders. Biol. Psychiatry 79, 765–775. doi: 10.1016/j.biopsych.2015.02.028
Marchant, N., Sanders, R., Sleigh, J., Vanhaudenhuyse, A., Bruno, M. A., Brichant, J. F., et al. (2014). How electroencephalography serves the anesthesiologist. Clin. EEG Neurosci. 45, 22–32. doi: 10.1177/1550059413509801
Martuzzi, R., Ramani, R., Qiu, M., Rajeevan, N., and Constable, R. T. (2010). Functional connectivity and alterations in baseline brain state in humans. Neuroimage 49, 823–834. doi: 10.1016/j.neuroimage.2009.07.028
Mhuircheartaigh, R. N., Rosenorn-Lanng, D., Wise, R., Jbabdi, S., Rogers, R., and Tracey, I. (2010). Cortical and subcortical connectivity changes during decreasing levels of consciousness in humans: a functional magnetic resonance imaging study using propofol. J. Neurosci. 30, 9095–9102. doi: 10.1523/JNEUROSCI.5516-09.2010
Nallasamy, N., and Tsao, D. Y. (2011). Functional connectivity in the brain: effects of anesthesia. Neuroscientist 17, 94–106. doi: 10.1177/1073858410374126
Nasrallah, F. A., Yeow, L. Y., Biswal, B., and Chuang, K. H. (2015). Dependence of BOLD signal fluctuation on arterial blood CO2 and O2: implication for resting-state functional connectivity. Neuroimage 117, 29–39. doi: 10.1016/j.neuroimage.2015.05.035
Paasonen, J., Salo, R. A., Huttunen, J. K., and Grohn, O. (2017). Resting-state functional MRI as a tool for evaluating brain hemodynamic responsiveness to external stimuli in rats. Magn. Reson. Med. 78, 1136–1146. doi: 10.1002/mrm.26496
Paasonen, J., Stenroos, P., Salo, R. A., Kiviniemi, V., and Grohn, O. (2018). Functional connectivity under six anesthesia protocols and the awake condition in rat brain. Neuroimage 172, 9–20. doi: 10.1016/j.neuroimage.2018.01.014
Purdon, P. L., Pierce, E. T., Mukamel, E. A., Prerau, M. J., Walsh, J. L., Wong, K. F., et al. (2013). Electroencephalogram signatures of loss and recovery of consciousness from propofol. Proc. Natl. Acad. Sci. U.S.A. 110, E1142–E1151. doi: 10.1073/pnas.1221180110
Purdon, P. L., Sampson, A., Pavone, K. J., and Brown, E. N. (2015). Clinical Electroencephalography for anesthesiologists: Part I: background and basic signatures. Anesthesiology 123, 937–960. doi: 10.1097/ALN.0000000000000841
Raichle, M. E., Macleod, A. M., Snyder, A. Z., Powers, W. J., Gusnard, D. A., and Shulman, G. L. (2001). A default mode of brain function. Proc. Natl. Acad. Sci. U.S.A. 98, 676–682. doi: 10.1073/pnas.98.2.676
Rohlfing, T., Kroenke, C. D., Sullivan, E. V., Dubach, M. F., Bowden, D. M., Grant, K. A., et al. (2012). The INIA19 template and neuromaps atlas for primate brain image parcellation and spatial normalization. Front. Neuroinform. 6:27. doi: 10.3389/fninf.2012.00027
Ryu, K. H., Song, K., Lim, T. Y., Choi, W. J., Kim, Y. H., and Kim, H. S. (2018). Does equi-minimum alveolar concentration value ensure equivalent analgesic or hypnotic potency: a comparison between desflurane and sevoflurane. Anesthesiology 128, 1092–1098. doi: 10.1097/ALN.0000000000002158
Samann, P. G., Wehrle, R., Hoehn, D., Spoormaker, V. I., Peters, H., Tully, C., et al. (2011). Development of the brain’s default mode network from wakefulness to slow wave sleep. Cereb. Cortex 21, 2082–2093. doi: 10.1093/cercor/bhq295
Scholvinck, M. L., Maier, A., Ye, F. Q., Duyn, J. H., and Leopold, D. A. (2010). Neural basis of global resting-state fMRI activity. Proc. Natl. Acad. Sci. U.S.A. 107, 10238–10243. doi: 10.1073/pnas.0913110107
Sheskin, D. J. (2003). Handbook of Parametric and Nonparametric Statistical Procedures, 3rd Edn. Boca Raton, FL: CRC Press. doi: 10.1201/9781420036268
Steriade, M., Amzica, F., and Contreras, D. (1994). Cortical and thalamic cellular correlates of electroencephalographic burst-suppression. Electroencephalogr. Clin. Neurophysiol. 90, 1–16. doi: 10.1016/0013-4694(94)90108-2
Teichert, T., Grinband, J., Hirsch, J., and Ferrera, V. P. (2010). Effects of heartbeat and respiration on macaque fMRI: implications for functional connectivity. Neuropsychologia 48, 1886–1894. doi: 10.1016/j.neuropsychologia.2009.11.026
Van Den Heuvel, M. P., and Hulshoff Pol, H. E. (2010). Exploring the brain network: a review on resting-state fMRI functional connectivity. Eur. Neuropsychopharmacol. 20, 519–534. doi: 10.1016/j.euroneuro.2010.03.008
Vanhaudenhuyse, A., Demertzi, A., Schabus, M., Noirhomme, Q., Bredart, S., Boly, M., et al. (2011). Two distinct neuronal networks mediate the awareness of environment and of self. J. Cogn. Neurosci. 23, 570–578. doi: 10.1162/jocn.2010.21488
Vincent, J. L., Patel, G. H., Fox, M. D., Snyder, A. Z., Baker, J. T., Van Essen, D. C., et al. (2007). Intrinsic functional architecture in the anaesthetized monkey brain. Nature 447, 83–86. doi: 10.1038/nature05758
Vogt, B. A., and Laureys, S. (2005). Posterior cingulate, precuneal and retrosplenial cortices: cytology and components of the neural network correlates of consciousness. Prog. Brain Res. 150, 205–217. doi: 10.1016/S0079-6123(05)50015-3
Wang, Z., Chen, L. M., Negyessy, L., Friedman, R. M., Mishra, A., Gore, J. C., et al. (2013). The relationship of anatomical and functional connectivity to resting-state connectivity in primate somatosensory cortex. Neuron 78, 1116–1126. doi: 10.1016/j.neuron.2013.04.023
White, N. S., and Alkire, M. T. (2003). Impaired thalamocortical connectivity in humans during general-anesthetic-induced unconsciousness. Neuroimage 19, 402–411. doi: 10.1016/S1053-8119(03)00103-4
Wu, T., Grandjean, J., Bosshard, S. C., Rudin, M., Reutens, D., and Jiang, T. (2017). Altered regional connectivity reflecting effects of different anaesthesia protocols in the mouse brain. Neuroimage 149, 190–199. doi: 10.1016/j.neuroimage.2017.01.074
Keywords: burst suppression, simultaneous EEG-MRI, functional connectivity, default mode network, thalamocortical connectivity, nonhuman primate
Citation: Zhang Z, Cai D-C, Wang Z, Zeljic K, Wang Z and Wang Y (2019) Isoflurane-Induced Burst Suppression Increases Intrinsic Functional Connectivity of the Monkey Brain. Front. Neurosci. 13:296. doi: 10.3389/fnins.2019.00296
Received: 02 November 2018; Accepted: 13 March 2019;
Published: 11 April 2019.
Edited by:
Peter Herman, Yale University, United StatesReviewed by:
Kai-Hsiang Chuang, The University of Queensland, AustraliaElzbieta Olejarczyk, Institute of Biocybernetics and Biomedical Engineering (PAN), Poland
Copyright © 2019 Zhang, Cai, Wang, Zeljic, Wang and Wang. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Zheng Wang, emhlbmcud2FuZ0Bpb24uYWMuY24= Yingwei Wang, d2FuZ3l3QHBlcmlvcGVyYXRpdmUtc2NpZW5jZS5jb20=
†These authors have contributed equally to this work