Faded Critical Dynamics in Adult Moyamoya Disease Revealed by EEG and fMRI

Abstract

Criticality is considered a dynamic signature of healthy brain activity that can be measured on the short-term timescale with neural avalanches and long-term timescale with long-range temporal correlation (LRTC). It is unclear how the brain dynamics change in adult moyamoya disease (MMD). We used BOLD-fMRI for LRTC analysis from 16 hemorrhagic ($H_{\text{MMD}}$) and 34 ischemic ($I_{\text{MMD}}$) patients and 25 healthy controls. Afterwards, they were examined by EEG recordings in the eyes-closed (EC), eyes-open (EO), and working memory (WM) states. The EEG data of 11 $H_{\text{MMD}}$ and 13 $I_{\text{MMD}}$ patients and 21 healthy controls were in good quality for analysis. Regarding the 4 metrics of neural avalanches (e.g., size ($\alpha$), duration ($\beta$), $\kappa$ value, and branching parameter ($\sigma$)), both MMD subtypes exhibited subcritical states in the EC state. When switching to the WM state, $H_{\text{MMD}}$ remained inactive, while $I_{\text{MMD}}$ surpassed controls and became supercritical ($p<0.05$). Regarding LRTC, the amplitude envelope in the EC state was more analogous to random noise in the MMD patients than in controls. During state transitions, LRTC decreased sharply in the controls but remained chaotic in the MMD individuals ($p<0.05$). The spatial LRTC reduction distribution based on both EEG and fMRI in the EC state implied that, compared with controls, the two MMD subtypes might exhibit mutually independent but partially overlapping patterns. The regions showing decreased LRTC in both EEG and fMRI were the left supplemental motor area of $H_{\text{MMD}}$ and right pre-/postcentral gyrus and right inferior temporal gyrus of $I_{\text{MMD}}$. This study not only sheds light on the decayed critical dynamics of MMD in both the resting and task states for the first time but also proposes several EEG and fMRI features to identify its two subtypes.

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Moyamoya disease (MMD) is a chronic cerebrovascular disease characterized by progressive stenosis and occlusion of the terminal portion of the bilateral internal carotid arteries and their main branches. The vascular pathology leads to widespread and continuous cerebral hypoperfusion and gradual formation of compensation from collaterals as a response [1]. The blood supply of neuron activities is often disrupted, followed by less regional neural activities and cognitive impairment [2–4]. The MMD is the only known chronic cerebrovascular disease with both hemorrhagic ($H_{\text{MMD}}$) and ischemic ($I_{\text{MMD}}$) clinical outcomes [5, 6]. Besides the general pathophysiology mentioned above in common, the two clinical subtypes are suspected to exhibit different but unknown processes, which are proved to be irrelevant with cerebral infarct or hematoma [3, 7–9]. However, their morphological and hemodynamic differences are difficult to be discerned and rarely reported [6, 9].

Healthy brain networks in the resting state are generally characterized by well-balanced excitatory and inhibitory synaptic activities [10]. These balanced brain states may sit in a dynamic state close to the criticality. Neural networks in a critical state have been found to be efficient for information transmission and other functions [11, 12]. Such critical dynamics are associated with a sustained neural activity that exhibits scale-free avalanche distribution [13]. Cognitive tasks may shift the brain activity from the critical point to a supercritical state by activating some brain regions and suppressing others [14]. The critical dynamics can be characterized quantitatively on two timescales: neural avalanches in the short term (milliseconds) and long-range temporal correlation (LRTC) in the long term (seconds to hours) [15]. During a neuronal avalanche, spontaneous activation of one neuronal group can trigger consecutive activations of other neuronal groups within just a few milliseconds, propagating cascading waves of activity. This phenomenon has been revealed by multiple neuroimaging modalities, such as electroencephalography (EEG) and fMRI in healthy subjects [16–18]. Additionally, brains in disease states such as unconsciousness and schizophrenia exhibit faded critical dynamics [19, 20].

LRTC is another notable measurement to evaluate the critical dynamics of the neuronal system [21, 22]. Although neural oscillations present themselves with large variability in both frequency and amplitude, their fluctuations reflect a tendency toward self-organized criticality [23]. The oscillatory activity at any time is influenced by previous activities, and the LRTC is built up though local interactions. Such oscillations, reflecting short-term and long-term interactions, are observed to exist throughout the entire system. Recent studies have shown that healthy brains in the resting state are associated with a relatively large value of LRTC, while decreased LRTC has been reported in various diseased brain states such as a major depressive disorder and Alzheimer’s disease [24, 25].

Clinically, rapid and accurate differential diagnosis between acute ischemia and hemorrhage is crucial for early medical and interventional treatment but the optimal time window of treatment is often missed because of its high reliance on a CT or MR scan. The MMD is expected as the promising disease template to develop a more rapid and accurate differential diagnostic tool. According to one published fMRI study of ours, a dynamic measurement of entropy was proposed as an index of the critical dynamics to describe quantitatively the spatiotemporal changes of neural communication in adult MMD [26]. It found that critical dynamics faded not only in the diseased brain but also with disease progression. Therefore, this study was performed to examine two issues as the first of its kind. One is to explore the faded critical dynamics of MMD in both the resting and task states directly through a combination of EEG (high temporal resolution) and fMRI (high spatial resolution). The other is to investigate whether the two subtypes with a similar extent of cognitive impairment exhibit different neuronal dynamics and generate several features for future rapid and bedside identification of acute cerebral ischemia and hemorrhage.

2. Materials and Methods

2.1. Participants

This study was approved by the Institutional Review Board in our hospital and conducted in accordance with the Helsinki declaration. Informed consent was signed by all the subjects of this study. From March 2017 to August 2018, 50 adult patients with MMD (16 $H_{\text{MMD}}$ and 34 $I_{\text{MMD}}$) were recruited. The inclusion criteria were as follows: (1) Chinese nationality, right handed, and age of over 18 years; (2) diagnosis through digital subtraction angiography, with a Suzuki grade of III or IV [1]; (3) no evidence of infarct and hematoma larger than 8 mm in the maximum dimension on structural brain images [2, 27, 28]; (4) physical ability to undergo cognitive testing and EEG tasks; (5) no severe systemic or other cerebrovascular diseases; and (6) no medical history of neurosurgery. Twenty-five healthy subjects with no cerebrovascular or mental diseases were enrolled as controls. The Mini-Mental State Examination (MMSE) was adopted for the neuropsychological evaluation of global cognitive states.

2.2. Data Acquisition and Preprocessing

EEG data were acquired at a sampling rate of 1000 Hz in a sound-attenuated room by using a 64-channel actiCHamp Brain Products recording system (Brain Products GmbH Inc., Munich, Germany). The impedance of all channels was below 10 KΩ. The experimental paradigm was presented in E-Prime 2.0, and preprocessing was performed using the MATLAB R2017b software plug-in EEGLAB 14.0.0. Data were filtered to the frequency range of 0.5–100 Hz. Then, the interference at a frequency of 50 Hz was removed using a notch filter. Independent component analysis and the ADJUST toolbox were used to remove the components of eyeblink and cardio ballistic artifacts. EEG recording was started with a 5-minute eyes-closed (EC) resting state, a 5-minute eyes-open (EO) resting state, and then 30 trials of a delayed-response working memory (WM) task for around 20 minutes and ended with a procedure composed of a 5-minute EC and a 5-minute EO resting state to examine the consistency of the data recording quality (Supplementary Figure S2). Ultimately, the EEG data from 24 patients (11 $H_{\text{MMD}}$ and 13 $I_{\text{MMD}}$) and 21 healthy controls were selected for good quality.

All MRI data were collected using a 3.0 Tesla scanner (GE Healthcare, GE Asian Hub, Shanghai, China) with a 32-channel intraoperative head coil. The fMRI data were acquired using gradient echo-planar imaging with the following parameters: 3.2 mm slice thickness, 2000 ms repetition time, 30 ms echo time, 90° flip angle, $220\times 220\text{m}{\text{m}}^2$ field of view, $64\times 64$ matrix size, and $3.438\times 3.438\times 3.2\text{m}{\text{m}}^3$ voxel size. Each scan lasted for 8 min and collected 240 volumes. Preprocessing was performed using the SPM12 (https://www.fil.ion.ucl.ac.uk/spm) and DPARSF [29]. Volumes were corrected for slice timing and head motion, and subjects were excluded if their head movement exceeded 3 mm or 3 degrees. Here, we used the Friston 24-parameter model [30] to regress out the head motion effects from the realigned data (i.e., 6 head motion parameters, 6 head motion parameters on the time point before, and the 12 corresponding squared items) that were based on recent reports that higher-order models demonstrated benefits in removing head motion effects [31]. We also obtained the mean framewise displacement (FD) derived from Jenkinson’s relative root mean square algorithm [32] for each participant to assess the effects of head motion on group-level statistics. After spatial normalization (Montreal Neurological Institute space), resampling (3 mm isotropic voxels), and spatial smoothing (4 mm, full-width, half-maximum Gaussian kernel), volumes were preprocessed using the linear trend subtraction and temporal filtering (0.01–0.08 Hz). In addition, the signals from the cerebrospinal fluid and white matter were regressed out to reduce the effect of nonneuronal BOLD fluctuation [33].

2.3. Short-Term Timescales of Critical Dynamics: Neuronal Avalanches

2.3.1. Transfer of EEG Data from Amplitude to Mean Frequency

The time sequence of each channel was translated from amplitude to mean frequency by the short-term Fourier transform (Figures 1(a) and 1(b)). The mean frequency of each 1-sec epoch was calculated as $$\begin{array}{}\text{(1)}& aveff=\frac{\sum Pf\times f}{\sum Pf},\end{array}$$with a 2 ms time shift in each channel, where $f$ and $P$ represent the frequency and power of the power spectrum of the signal, respectively. The average baseline of each channel was replaced with the group ensemble average of the channel to avoid noise interference caused by individual difference perturbations. Given that the EEG amplitude is mainly dominated by the low-frequency component of synchronized slow oscillations and some high-frequency activities are submerged in the slow-wave activities with high amplitudes, the EEG frequency may be a better sequence for reflecting the intensity of the neural activity rate.

[figures omitted; refer to PDF]

2.3.2. Overthreshold Event Detection and Neuronal Avalanche Extraction

For each channel, the threshold for events was defined as the mean $\text{value}+1$ SD (standard deviation) of the EC state averaged overall controls. It is noted that for frequency events, 1 SD fluctuation above the mean frequency rate is strong enough to be a burst event (two and three SD were also tested, but the events were very sparse). The function “findpeaks” in MATLAB 2015b was used to find the peaks higher than the threshold (with a duration over 20 ms), which were defined as the events. Neuronal avalanches were defined as continuous time bins with at least one event in each channel, beginning and ending with event-free time bins. The time bin used was 16 ms (although the conclusion held true with other time bins from 4 to 20 ms). Binarized data were used for event detection.

2.3.3. Calculation of the Avalanche Size ($\alpha$) and Duration ($\beta$)

As shown in Figures 1(c) and 1(d), we first obtained the probability density function (PDF) of the avalanche size (number of events) and duration (number of time bins) and then fitted it to the power law distribution $P_{\alpha }\chi =C_{\alpha }{\chi }^{\alpha }$ using maximum likelihood estimation to obtain the slope, denoted by the size ($\alpha$) and duration ($\beta$) [34, 35]. The exponent was estimated by calculating $L\alpha \mid x={\prod }_{i=1}^n{P}_{\alpha }{x}_i$, which means that the best-fitting exponent was calculated by maximizing the log-likelihood function. The higher absolute value of the power law exponent means the steeper fitted power law curve and indicates the higher probability of small size and duration avalanche [18]. Since the large size and duration avalanche means a rather large range and longer time of neuronal dynamics, the comparison of these parameters between the MMD and controls can be utilized to reflect the aberrant neuronal dynamics of MMD.

2.3.4. Calculation of the $\kappa$ Value

The $\kappa$ value was used to estimate how far the system was from criticality. This nonparametric measure quantifies the difference in the subjects’ cumulative density function (CDF) from the standard critical reference CDF with a power law exponent of 1.5 (Figure 1(e)) [11, 36]: $$\begin{array}{}\text{(2)}& \kappa =1+\frac{1}{m}\sum _{k=1}^m{F}^{NA}{\beta }_k-F{\beta }_k.\end{array}$$

The value of this parameter at approximately 1 indicates a critical state of the system, while values above or below 1 indicate the super- or subcritical states, respectively.

2.3.5. Calculation of the Branching Parameter ($\sigma$)

The branching parameter $\sigma$ is defined as the average number of subsequent events triggered by a single preceding event in an avalanche [34]. Thus, it was first calculated by averaging the ratio of the numbers of events in the second bin of an avalanche to the number of events in the first bin and then averaging all avalanches of a subject’s time series: $$\begin{array}{}\text{(3)}& \sigma =\frac{1}{N_{\text{av}}}\sum _{k=1}^{N_{av}}\frac{2\text{nd}\text{bin}\text{of}{k}^{\prime }\text{th}\text{avalanche}}{n_{\text{events}}1\text{st}\text{bin}\text{of}{k}^{\prime }\text{th}\text{avalanche}},\end{array}$$where $N$ is the total number of cascades in the time series and $n$ is the number of the events in a certain time bin. For the single-bin cascades, the branching parameter is equal to 0. Therefore, the branching parameter of each cascade can vary over a wide range and the number of the cascades must be at least 700 in this study.

2.4. Long-Term Timescales of Critical Dynamics: LRTC

2.4.1. LRTC Based on EEG Data

Detrended fluctuation analysis was used to calculate the Hurst exponent of LRTC [37, 38]. The EEG data were bandpass filtered (finite impulse response filter) to the delta (1–4 Hz), theta (4–8 Hz), alpha (8–12 Hz), beta (12–25 Hz), and gamma (25–100 Hz) bands using a filter order of 2/minimum frequency of the band. Given a consecutive time series $yt$, the first step was to subtract the mean signal and obtain the cumulative sum of the signal: $$\begin{array}{}\text{(4)}& xt=\sum _{k=1}^{t}yk-y,\end{array}$$where $<,>$ denotes the time average. Then, the signal was divided into a set of time windows $W$ of the same length $L$ with 50% overlap. For each window in $W$, the linear trend was removed using the least-squares method. Then, we computed the average root-mean-square fluctuation, $FL$, of all the time windows in $W$ with an identical length $L$. The time window length was logarithmically spaced with a lower bound of 4 samples and an upper bound of 10% of the signal length [38]. According to our timescale in this work, we used 29 different window lengths $L$ ranging from 100 samples to 10000 samples (0.1 ~ 10 s). The Hurst exponent $H$ was defined as the coefficient of the linear regression of the sequence $Rn/Sn$ plotted on a log-log scale. For a scale-free signal $FL\propto L^H$, in the case of $0<H<1$, the time series was correlated. When $H$ ranged from 0.5 to 1, the signal was considered to exhibit positive autocorrelation, while $H=0.5$ indicated an uncorrelated signal.

2.4.2. LRTC Based on fMRI Data

The voxel-wise Hurst exponent was adopted to evaluate the LRTC of fMRI data through the classical rescaled range (R/S) analysis [39]. The BOLD signal time series was divided into multiple shorter time series (length, $n$). The number of subtime series was $M$, and $N=M\times n$ was the length of the full-time series. For the subtime series of length $n$, $X_m=X_{m1},X_{m2},\cdots ,X_{mn}$, where $m=\mathrm{1,2,3},\cdots ,M$, we calculated the rescaled range $Rn/Sn$: $$\begin{array}{}\text{(5)}& \frac{Rn}{Sn}=\frac{1}{M}\sum _{m=1}^M\frac{R_{mn}}{Smn}.\end{array}$$

The $R_{mn}$ is the difference between the extremes of subtime series $X_m$. $CdXm$ is the cumulative deviation from $i$ to $k$ of the subtime series: $$\begin{array}{}\text{(6)}& Cd{X}_m=\sum _{i=1}^k{X}_{mi}-{\overline{X}}_m,\end{array}$$where $1\le k\le n$. $R_{mn}$ is the difference between the maximum and minimum of $Cd{X}_m$ and minimum: $$\begin{array}{}\text{(7)}& R_{mn}=\underset{1\le k\le n}{\max }\sum _{i=1}^k{X}_{mi}-{\overline{X}}_m-\underset{1\le k\le n}{\min }\sum _{i=1}^k{X}_{mi}-{\overline{X}}_m.\end{array}$$

The $S_{mn}$ is the standard deviation of subtime series $X_m$: $$\begin{array}{}\text{(8)}& S_{mn}=\sqrt{\frac{1}{n}\sum _{i=1}^n{X_{mi}-{\overline{X}}_m}^2}.\end{array}$$

Calculating all the subtime series, we obtained the sequence $Rn/Sn$ corresponding to the sequence $n$, where $n$ was all the possible values of the subtime series length.

2.5. Statistical Analysis

For EEG data analysis, one-way ANOVA was used and all the pairwise comparisons were assessed using Tukey-Kramer’s multiple comparison method. For fMRI data analysis, one-sample $t$-tests were first performed on individual Hurst exponent maps for each group in a voxel-wise manner and to explore the within-group LRTC patterns using Rest version 1.8 ($Z>3.48$ for voxel-wise $p<0.001$ and cluster-wise $p<0.05$, corrected through Gaussian random field theory) [40]. Afterwards, the averaged whole-brain Hurst exponent values were generated and compared among the three groups using the one-way ANOVA. Then, two-sample $t$-tests were performed between each pair of subgroups on individual Hurst exponent maps ($Z>3.29$ for voxel-wise $p<0.001$ and cluster-wise $p<0.05$).

3. Results

3.1. Clinical Information

Table 1 shows the clinical information of the involved subjects. The demographic differences among the three groups were not significant ($p>0.05$). However, the mean MMSE score of the healthy controls was significantly higher than that of the patient groups ($p<0.001$). Among the subjects who were selected for the EEG analysis, differences of the demographics and cognitive and task performances among the three groups were not significant ($p>0.05$).

Table 1

Demographic features of the 3 groups.

Index	$H_{\text{MMD}}$	$I_{\text{MMD}}$	Controls	$F/{\chi }^2$ value ($p$ value)
fMRI cohort
Subjects, $n$	16	34	25	/
Age (years)	$$42.31\pm 7.87$$	$$40.41\pm 9.71$$	$$41.32\pm 11.38$$	0.815 (0.309)
Male (%)	6 (37.5)	16 (47.1)	12 (48.0)	0.509 (0.775)
Education (years)	$$9.06\pm 5.13$$	$$9.38\pm 3.99$$	$$8.96\pm 4.37$$	0.729 (0.694)
MMSE	$$22.75\pm 4.86$$	$$25.44\pm 3.82$$	$$28.2\pm 1.71$$	22.076 (<0.001)
EEG cohort
Subjects, $n$	11	13	21	/
Age (years)	$$40.91\pm 9.58$$	$$41.23\pm 14.60$$	$$43.38\pm 9.06$$	0.326 (0.858)
Male (%)	4 (36.4)	5 (38.5)	12 (57.1)	2.057 (0.357)
Education (years)	$$8.84\pm 2.60$$	$$8.72\pm 2.45$$	$$8.33\pm 2.39$$	0.408 (0.815)
MMSE	$$24.92\pm 1.44$$	$$24.82\pm 1.17$$	$$25.16\pm 2.39$$	2.305 (0.316)
WM accuracy	$$0.773\pm 0.217$$	$$0.759\pm 0.150$$	$$0.827\pm 0.146$$	0.670 (0.7165)

$H_{\text{MMD}}$: hemorrhagic moyamoya disease; $I_{\text{MMD}}$: ischemic moyamoya disease; MMSE: Mini-Mental State Examination.

3.2. Faded Critical Dynamics of Neuronal Avalanches

3.2.1. Cascade Size ($\alpha$) and Duration ($\beta$)

Figure 2 indicates that the three groups exhibited significant differences in both $\alpha$ (EC, $F=14.11$, $p=2.6\times {10}^{-5}$; EO, $F=68.06$, $p=2.75\times {10}^{-13}$; and WM, $F=5.07$, $p=0.0115$) and $\beta$ (EC, $F=9.88$, $p=0.0004$; EO, $F=10.37$, $p=0.0003$; and WM, $F=5.08$, $p=0.0114$). In the subgroup analysis, $H_{\text{MMD}}$ and $I_{\text{MMD}}$ exhibited significant differences in $\alpha$ in the EO (${\overline{\alpha }}_H=1.51\pm 0.01,$${\overline{\alpha }}_I=1.47\pm 0.005,$$p=1.03\times {10}^{-9}$) and WM (${\overline{\alpha }}_H=1.48\pm 0.01,$${\overline{\alpha }}_I=1.42\pm 0.01$, $p=0.0082$) states and in $\beta$ in the WM state (${\overline{\beta }}_H=1.90\pm 0.04,$${\overline{\beta }}_I=1.77\pm 0.03$, $p=0.0088$). The control and $I_{\text{MMD}}$ groups showed significant differences in both $\alpha$ (${\overline{\alpha }}_C=2.08\pm 0.02,$${\overline{\alpha }}_I=1.68\pm 0.02$, $p=3.35\times {10}^{-5}$) and $\beta$ (${\overline{\beta }}_C=2.08\pm 0.02,$${\overline{\beta }}_I=2.31\pm 0.04$, $p=5.23\times {10}^{-4}$) in the EC state. Additionally, the control and $H_{\text{MMD}}$ groups exhibited significant differences in $\alpha$ in both the EC (${\overline{\alpha }}_C=2.08\pm 0.02$, ${\overline{\alpha }}_H=1.63\pm 0.02$, $p=0.0059$) and EO (${\overline{\alpha }}_C=1.46\pm 0.007,{\overline{\alpha }}_H=1.51\pm 0.01$, $p=9.57\times {10}^{-10}$) states and in $\beta$ in the EO state (${\overline{\beta }}_C=1.85\pm 0.01,$${\overline{\beta }}_H=2.00\pm 0.04$, $p=1.53\times {10}^{-4}$).

[figures omitted; refer to PDF]

3.2.2. Deviations from Criticality with κ and the Branching Parameter of $\sigma$

The three groups exhibited significant differences in both $\kappa$ (EC, $F=15.51$, $p=1.18\times {10}^{-5}$; EO, $F=66.34$, $p=4.02\times {10}^{-13}$; and WM, $F=5.85$, $p=0.0063$) and $\sigma$ (EC, $F=23.01$, $p=2.84\times {10}^{-7}$; EO, $F=11.35$, $p=0.0001$; and WM, $F=5.33$, $p=0.0094$) in all three states. $H_{\text{MMD}}$ and $I_{\text{MMD}}$ showed significant differences in $\kappa$ in the EO (${\overline{\kappa }}_H=0.97\pm 0.006,$${\overline{\kappa }}_I=0.99\pm 0.003,$$p=1.26\times {10}^{-9}$) and WM (${\overline{\kappa }}_H=0.99\pm 0.007,$${\overline{\kappa }}_I=1.03\pm 0.008$, $p=0.0044$) states and in $\sigma$ in the WM state (${\overline{\sigma }}_H=1.14\pm 0.045,$${\overline{\sigma }}_I=1.32\pm 0.04$, $p=0.0027$). The control and $I_{\text{MMD}}$ groups showed significant differences in both $\kappa$ (${\overline{\kappa }}_C=0.95\pm 0.004,$${\overline{\kappa }}_I=0.91\pm 0.005$, $p=2.03\times {10}^{-5}$) and $\sigma$ (${\overline{\sigma }}_C=0.96\pm 0.21,$${\overline{\sigma }}_I=0.70\pm 0.005$, $p=1.43\times {10}^{-6}$) in the EC state. Additionally, the control and $H_{\text{MMD}}$ showed significant differences in $\kappa$ in both the EC (${\overline{\kappa }}_C=0.95\pm 0.004,$${\overline{\kappa }}_H=0.93\pm 0.006$, $p=0.0022$) and EO (${\overline{\kappa }}_C=1.01\pm 0.004,$${\overline{\kappa }}_H=0.97\pm 0.006,$$p=9.57\times {10}^{-10}$) states and in $\sigma$ in both the EC (${\overline{\sigma }}_C=0.96\pm 0.21,$${\overline{\sigma }}_H=0.74\pm 0.03$, $p=4.42\times {10}^{-5}$) and EO (${\overline{\sigma }}_C=1.19\pm 0.02,$${\overline{\sigma }}_H=1.04\pm 0.02$, $p=1.13\times {10}^{-4}$) states. The results kept the same tendency when utilizing a different parameter to calculate.

3.2.3. State Transition

The three groups exhibited a similar transition tendency between different states in all parameters (Figure 2(d)). For example, considering $\alpha$, both subtypes of MMD presented with larger values than controls in the EC state, implying a subcritical state. Similarly, MMD exhibited less activity than controls when switching to the EO state. Interestingly, when changing into the WM state, the $H_{\text{MMD}}$ group remained inactive but the $I_{\text{MMD}}$ group exhibited supercritical dynamics that were even more active than those of the controls.

3.3. Faded Critical Dynamics of LRTC

3.3.1. LRTC Based on EEG Data

The alpha band signal was taken for Hurst exponent analysis. Figure 3(a) indicates that the amplitude envelopes of the MMD groups were more analogous to the random noise than the control amplitude envelope was (the EC state, with a central-parietal channel as an example). The averaged Hurst exponents of the $H_{\text{MMD}}$ and $I_{\text{MMD}}$ groups were approximately 0.63 and 0.70, respectively, while that of the controls was 0.87 (Figure 3(b)).

[figures omitted; refer to PDF]

The Hurst exponents of all the recording channels were calculated and mapped in Figure 3(c). Then, channels with significant differences in Hurst exponents among the three groups were mapped (Figure 3(d), Supplementary Table S1, ANOVA, $p<0.05$). Afterwards, the mean Hurst exponent value of the generated channels was calculated and compared among the three groups (Figure 3(e)). In the EC state, the control group exhibited significantly higher values ${\overline{H}}_C=0.85\pm 0.02$ in these channels than either $H_{\text{MMD}}$ (${\overline{H}}_H=0.78\pm 0.02$, $p=0.0089$) or $I_{\text{MMD}}$ (${\overline{H}}_I=0.76\pm 0.02$, $p=0.005$). In the EO state, more channels presented a significant group difference and the $I_{\text{MMD}}$ showed a significantly higher value ${\overline{H}}_I=0.76\pm 0.01$ than the controls ${\overline{H}}_C=0.66\pm 0.02$ in these channels ($p=0.0062$). Furthermore, even more channels showed significant group differences in Hurst exponents in the WM state and the controls showed significantly lower values ${\overline{H}}_C=0.66\pm 0.02$ than the $H_{\text{MMD}}$ group (${\overline{H}}_H=0.76\pm 0.02$, $p=0.0019$) or the $I_{\text{MMD}}$ group (${\overline{H}}_I=0.74\pm 0.02$, $p=0.0080$).

When switching from the EC to EO state, 18 channels had significant Hurst exponent changes in the $I_{\text{MMD}}\text{group}$ compared to the control group (Figures 4(a), $p<0.05$). Interestingly, the Hurst exponents of these 18 channels all decreased in the controls, while the $I_{\text{MMD}}$ group exhibited a chaotic pattern of change in these channels. In comparison, only 4 channels showed significant changes in the $H_{\text{MMD}}$ group compared to controls (Figure 4(b), $p<0.05$). Similarly, the Hurst exponents of these 4 channels all decreased in the controls, while the $H_{\text{MMD}}$ also showed a chaotic changing pattern in these channels. When switching from the EC to WM state, 27 channels in the $I_{\text{MMD}}$ and 28 channels in the $H_{\text{MMD}}$ showed significant changes in their Hurst exponents compared to those in the controls (Figures 4(c) and 4(d), $p<0.05$). Similarly, all these channels exhibited a large decline in Hurst exponents in the controls, while presenting with no regular pattern of change in either $I_{\text{MMD}}$ or $H_{\text{MMD}}$.

[figures omitted; refer to PDF]

3.3.2. LRTC Based on fMRI Data

The Hurst exponent patterns of the three groups are presented in Figures 5(a)–5(c). Visual inspection indicated that in all three groups, the bilateral orbital frontal gyrus (OFG) and left precuneus (PCu) exhibited high values, while the bilateral fusiform gyrus (FFG), left inferior temporal gyrus(ITG), and left insular gyri (IG) showed low values. In addition, the bilateral supplemental motor area (SMA), precentral gyrus (PreCG), and postcentral (PoCG) gyrus of the controls; the left medial superior frontal gyrus (SFGmed) and right PCu of $H_{\text{MMD}}$ patients; and the bilateral SFGmed, dorsolateral prefrontal gyrus (DLPFC), SMA, and right PCu of $I_{\text{MMD}}$ patients all exhibited high values. In addition, the bilateral caudate nucleus (CN), hippocampus (HIP), parahippocampal gyrus (PHG), and thalamus (THA) of both controls and $I_{\text{MMD}}$ showed low values. Similar to the EEG results in the EC state (Figure 3(e)), the controls exhibited the highest mean Hurst exponent value among the three groups (Figure 5(d), ANOVA, $F=3.397$, $p=0.038$), showing as controls $0.786\pm 0.018>I_{\text{MMD}}0.779\pm 0.016>H_{\text{MMD}}0.774\pm 0.011$.

[figures omitted; refer to PDF]

Compared with controls, significant decreases in the Hurst exponent in the $H_{\text{MMD}}$ group were found in the left SMA, left DLPFC, left PCu, left superior parietal gyrus (SPG), and left middle occipital gyrus (MOG). Afterwards, significant decreases in the $I_{\text{MMD}}$ group compared to controls were found in the bilateral DLPFC, left SMA, right PreCG, right PoCG, and right ITG. No regions in MMD exhibited significantly higher Hurst exponent values than the controls. In addition, the left PreCG exhibited a significantly lower value in $H_{\text{MMD}}$ than in $I_{\text{MMD}}$ but no region showed a significant difference in the opposite direction (Figure 6, Table 2).

[figures omitted; refer to PDF]

Table 2

Regional LRTC differences between each pair of the three groups.

Brain regions			MNI coordinates (mm)
	BA	Vol (mm3)	$x$	$y$	$z$	Maximum $Z$
$H_{\text{MMD}}$ vs controls
Left MOG	19/39/37	783	−30	−75	24	4.627
Left SMA	6	567	−12	−6	78	4.533
Left PCu	7/5	324	−15	−63	63	4.389
Left SPG	5/7	405	−18	−63	63	4.402
Left DLPFC	6	918	−15	−3	78	4.311
$I_{\text{MMD}}$ vs controls
Left DLPFC	6	675	−18	−6	78	5.058
Left SMA	6	459	−12	0	78	4.808
Right PoCG	4/3	540	15	−30	75	4.606
Right DLPFC	6/8	972	33	0	63	4.458
Right ITG	19/37	702	48	−69	−6	4.232
Right PreCG	4/6	405	15	−27	75	4.137
$H_{\text{MMD}}$ vs $I_{\text{MMD}}$
Left PreCG	6	540	−51	3	27	4.156

The $x$, $y$, and $z$ coordinates represent the primary peak location in the MNI space. BA: Brodmann’s area; SMA: supplementary motor area; DLPFC: dorsolateral prefrontal gyrus; MNI: Montreal Neurological Institute; MOG: middle occipital gyrus; PCu: precuneus; SPG: superior parietal gyrus; PoCG: postcentral gyrus; ITG: inferior temporal gyrus; PreCG: precentral gyrus.

Since head micromovements could introduce artefactual interindividual differences in resting-state fMRI metrics [41, 42], we also measured the difference of head motion among the three groups. Although the $I_{\text{MMD}}$ exhibited the highest FD (Welch’s ANOVA, $p=0.045$, Supplementary Figure S5), the mean Hurst exponent was not correlated significantly with mean FD across participants ($p>0.05$, Supplementary Figure S6). The impact of head motion on the group-level differences was then tested by adding mean FD as a covariate, and results confirmed that head motion was not responsible for the differences among the three groups (Supplementary Figure S7 and Table S2).

3.3.3. LRTC Colocalization Patterns Based on EEG and fMRI Data in the EC State

The EEG channel placement was projected onto the cortical surface and converted to the Talairach Stereotactic System based on a published Brodmann’s area (BA) atlas [43]. Afterwards, we assessed the colocalization of regions with LRTC abnormities based on EEG data and fMRI data in the EC state. Compared to controls, regions with a significant Hurst exponent decrease in the $H_{\text{MMD}}$ group based on overlapping data were found in the left SMA (BA6, Figure 7(a)). In addition, regions with a significant decrease in the $I_{\text{MMD}}$ group based on overlapping data compared to controls were found in the right PreCG and PoCG (BA4 and BA6) and right ITG (BA19, Figure 7(b)).

[figures omitted; refer to PDF]

4. Discussion

The criticality theory provides a novel insight into the neuronal dynamics underlying brain disorders. This study was the first to apply multiscale critical dynamics analysis to examine multimodal dynamical features in two moyamoya subtypes as compared to healthy controls. The neuronal avalanches on both fast and slow timescales were analyzed during rest and task performance, and several critical EEG features were derived. Both hemorrhagic and ischemic MMD exhibited particularly low EEG frequency activity and distinct subcritical dynamics, which could be distinguished easily from those of healthy controls. In addition, the decreased long-term correlations revealed in both high temporal (EEG) and spatial (fMRI) resolution were observed to reflect distinct neurophysiological processes associated with abnormal vascular network patterns in hemorrhagic and ischemic brains. Besides, this study provided clues for further rapid differential diagnosis between acute stroke and hemorrhage at the very early phase by use of EEG instead of CT and MR, which had greater advantages of rapidness, convenience, low cost, and radiation safe. Undoubtedly, time is the brain in treatment of acute ischemic stroke [44].

Previous investigations have suggested that the healthy brain in the resting state is usually characterized by well-balanced excitatory and inhibitory synaptic activities. These balanced levels of excitation and inhibition drive irregular spontaneous firing activities that exhibit scale-free avalanche distributions in the brain. Such a scale-free state can be effectively described by criticality [13]. Any input stimulus could effectively drive the brain into a supercritical state with additional excitatory activity, while relaxed low-signal states such as sleep can slow down the activity of the brain and shift it into a subcritical state. We noted that both subtypes of MMD exhibited subcritical states in the EC state and these suppressed dynamics prevented adaptive switching of brain function from introspective to extrospective states. This phenomenon might result from serious metabolic decrease and low neural activity rates caused by the chronic steno-occlusive angiopathy of MMD.

When switched to the EO state, both healthy controls and MMD presented with more neural activity and the MMD group remained less active than the controls. However, the ischemic moyamoya brains demonstrate stronger neural activities than hemorrhagic ones in the EO state. When switched to the WM state, all three groups exhibited more neural activities than that in the EO state, as was expected. Interestingly, the ischemic subtype surpassed the controls, while the hemorrhagic subtype still remained the least active. All parameters of neuronal avalanches exhibited similar results and were mutually verified. For healthy subjects, this phenomenon is reasonable because working memory is a behavioral state and requires effective neuronal activity to accomplish tasks [45]. However, we note that patients with ischemic MMD need to maintain a supercritical status to achieve similar scores to patients with hemorrhagic MMD and healthy subjects in certain WM tasks. In addition, these EEG features of neuronal avalanches during tasks may provide valuable clues for understanding the different neurophysiological processes of $I_{\text{MMD}}$ and $H_{\text{MMD}}$.

For healthy brains in the EC state, spatially distributed neuronal activity may oscillate in phase with each other and result in high LRTC values [46–48]. Compared with controls, the MMD group exhibited a reduced LRTC value close to 0.5 in the EC state, implying less correlated and more random brain activity. When subjects switched to EO and WM states, however, the LRTC value decreased sharply in the controls but both the $I_{\text{MMD}}$ and $H_{\text{MMD}}$ groups remained chaotic. The controls exhibited the lowest LRTC value in both the EO and WM states. In addition, spatial patterns of LRTC differences among the three groups were mapped based on EEG channels in both resting and WM states and under state transition. The healthy brain breaks the overall non-task-related dynamic balance into a set of subnetworks, such as the default-mode network, to deal with input signals effectively during state transition [49, 50]. However, the abnormal LRTC values and distribution of MMD in this study imply a completely different breakdown process for long-range neuronal fluctuations.

To further locate the regions with significant LRTC abnormities in MMD, we also examined BOLD fluctuations on fMRI due to the high spatial resolution of this modality. The results indicate that in the EC state, the patterns of the LRTC decreases in the two moyamoya subtypes are mutually independent but overlap in the left DLPFC of the executive control network and the left SMA of the salience network. Nevertheless, all regions of these patterns are key nodes involved in planning or direct control of movement, language, and visual information [51]. Referring to the pathophysiological nature of MMD, chronic stenosis/occlusion of the anterior circulation (bilateral internal carotid arteries and their main branches) is often followed by the collaterals from bilateral external carotid arteries and posterior circulation (bilateral vertebrobasilar arteries). Thus, the mismatch of the anterior circulation degradation and collateral development often results in a seemingly random and individualized cerebral hypoperfusion [7, 27, 52]. However, previous fMRI studies of MMD all revealed that patterns of functional deterioration are not random and key nodes of brain network such as the DLPFC, left SMA, are always involved [3, 4, 26, 28, 53]. Thus, this paper provides a crucial evidence that to output a similar extent of cognitive impairment, the neurophysiological processes of the two moyamoya subtypes may be mutually independent but overlap in some key nodes of the brain networks. Furthermore, we wondered whether there were potential links between spatial delay and temporal decay of neuronal oscillations in MMD and we attempted to trace their anatomical basis through both EEG and fMRI in the EC state. The identified regions are believed to play key roles in the neurophysiological processes of cognitive impairment in both ischemic and hemorrhagic MMD.

This study has several limitations that must be addressed. First, the EEG and fMRI data were not acquired at the same time. In order to generate a more stable and reliable result, the simultaneous EEG-fMRI technology should be used in future studies. Second, the study is based on a small sample size because completing EEG tasks is difficult for moyamoya patients with executive dysfunction. More patients are in need not only to increase the statistical power but to involve patients with Suzuki grading I–II and V–VI so as to obtain more knowledge of disease progression. Nevertheless, this study is the first of its kind to characterize the variability of brain dynamics in MMD on both short-term and long-term timescales and to show different neurophysiological features of its hemorrhagic and ischemic subtypes.

Authors’ Contributions

Y Mao, YX Gu, and YG Yu provided the idea of the article. Y Lei, YZ Li, and LC Yu wrote the article. LZ Xu, GX Zheng, and XY Qi analyzed the data. X Zhang, L Chen, and W Zhang collected the data. Yu Lei, Yuzhu Li, and Lianchun Yu contributed equally to this work.

References

[1] J. Suzuki, A. Takaku, "Cerebrovascular “moyamoya” disease. Disease showing abnormal net-like vessels in base of brain," Archives of Neurology, vol. 20 no. 3, pp. 288-299, DOI: 10.1001/archneur.1969.00480090076012, 1969.

[2] S. H. Su, J. Hai, L. Zhang, F. Yu, Y. F. Wu, "Assessment of cognitive function in adult patients with hemorrhagic moyamoya disease who received no surgical revascularization," European Journal of Neurology, vol. 20 no. 7, pp. 1081-1087, DOI: 10.1111/ene.12138, 2013.

[3] Y. Lei, Y. Li, W. Ni, H. Jiang, Z. Yang, Q. Guo, Y. Gu, Y. Mao, "Spontaneous brain activity in adult patients with moyamoya disease: a resting-state fMRI study," Brain Research, vol. 1546, pp. 27-33, DOI: 10.1016/j.brainres.2013.12.022, 2014.

[4] K. Kazumata, K. K. Tha, H. Narita, I. Kusumi, H. Shichinohe, M. Ito, N. Nakayama, K. Houkin, "Chronic ischemia alters brain microstructural integrity and cognitive performance in adult moyamoya disease," Stroke, vol. 46 no. 2, pp. 354-360, DOI: 10.1161/STROKEAHA.114.007407, 2015.

[5] G. M. Burke, A. M. Burke, A. K. Sherma, M. C. Hurley, H. H. Batjer, B. R. Bendok, "Moyamoya disease: a summary," Neurosurgical Focus, vol. 26 no. 4, article E11,DOI: 10.3171/2009.1.FOCUS08310, 2009.

[6] R. M. Scott, E. R. Smith, "Moyamoya disease and moyamoya syndrome," The New England Journal of Medicine, vol. 360 no. 12, pp. 1226-1237, DOI: 10.1056/NEJMra0804622, 2009.

[7] Y. Lei, Y. J. Li, Q. H. Guo, X. D. Liu, Z. Liu, W. Ni, J. B. Su, H. Yang, H. Q. Jiang, B. Xu, Y. X. Gu, Y. Mao, "Postoperative executive function in adult moyamoya disease: a preliminary study of its functional anatomy and behavioral correlates," Journal of Neurosurgery, vol. 126 no. 2, pp. 527-536, DOI: 10.3171/2015.12.JNS151499, 2017.

[8] S. Kang, X. Liu, D. Zhang, R. Wang, Y. Zhang, Q. Zhang, W. Yang, J. Z. Zhao, "Natural course of moyamoya disease in patients with prior hemorrhagic stroke," Stroke, vol. 50 no. 5, pp. 1060-1066, DOI: 10.1161/STROKEAHA.118.022771, 2019.

[9] M. Zhang, J. Tang, N. Liu, Y. Xue, X. Ren, J. Fu, "Postoperative functional outcomes and prognostic factors in two types of adult moyamoya diseases," Journal of Stroke and Cerebrovascular Diseases, vol. 29 no. 8, article 104846,DOI: 10.1016/j.jstrokecerebrovasdis.2020.104846, 2020.

[10] S. S. Poil, R. Hardstone, H. D. Mansvelder, K. Linkenkaer-Hansen, "Critical-state dynamics of avalanches and oscillations jointly emerge from balanced excitation/inhibition in neuronal networks," The Journal of Neuroscience, vol. 32 no. 29, pp. 9817-9823, DOI: 10.1523/JNEUROSCI.5990-11.2012, 2012.

[11] W. L. Shew, H. Yang, T. Petermann, R. Roy, D. Plenz, "Neuronal avalanches imply maximum dynamic range in cortical networks at criticality," The Journal of Neuroscience, vol. 29 no. 49, pp. 15595-15600, DOI: 10.1523/JNEUROSCI.3864-09.2009, 2009.

[12] W. L. Shew, H. Yang, S. Yu, R. Roy, D. Plenz, "Information capacity and transmission are maximized in balanced cortical networks with neuronal avalanches," The Journal of Neuroscience, vol. 31 no. 1, pp. 55-63, DOI: 10.1523/JNEUROSCI.4637-10.2011, 2011.

[13] O. Shriki, J. Alstott, F. Carver, T. Holroyd, R. N. A. Henson, M. L. Smith, R. Coppola, E. Bullmore, D. Plenz, "Neuronal avalanches in the resting MEG of the human brain," The Journal of Neuroscience, vol. 33 no. 16, pp. 7079-7090, DOI: 10.1523/JNEUROSCI.4286-12.2013, 2013.

[14] E. D. Fagerholm, R. Lorenz, G. Scott, M. Dinov, P. J. Hellyer, N. Mirzaei, C. Leeson, D. W. Carmichael, D. J. Sharp, W. L. Shew, R. Leech, "Cascades and cognitive state: focused attention incurs subcritical dynamics," The Journal of Neuroscience, vol. 35 no. 11, pp. 4626-4634, DOI: 10.1523/JNEUROSCI.3694-14.2015, 2015.

[15] A. Zhigalov, G. Arnulfo, L. Nobili, S. Palva, J. M. Palva, "Relationship of fast- and slow-timescale neuronal dynamics in human MEG and SEEG," The Journal of Neuroscience, vol. 35 no. 13, pp. 5385-5396, DOI: 10.1523/JNEUROSCI.4880-14.2015, 2015.

[16] G. Solovey, K. J. Miller, J. G. Ojemann, M. O. Magnasco, G. A. Cecchi, "Self-regulated dynamical criticality in human ECoG," Frontiers in Integrative Neuroscience, vol. 6, 2012.

[17] E. Tagliazucchi, P. Balenzuela, D. Fraiman, D. R. Chialvo, "Criticality in large-scale brain FMRI dynamics unveiled by a novel point process analysis," Frontiers in Physiology, vol. 3, 2012.

[18] L. Cocchi, L. L. Gollo, A. Zalesky, M. Breakspear, "Criticality in the brain: a synthesis of neurobiology, models and cognition," Progress in Neurobiology, vol. 158, pp. 132-152, DOI: 10.1016/j.pneurobio.2017.07.002, 2017.

[19] C. Meisel, E. Olbrich, O. Shriki, P. Achermann, "Fading signatures of critical brain dynamics during sustained wakefulness in humans," The Journal of Neuroscience, vol. 33 no. 44, pp. 17363-17372, DOI: 10.1523/JNEUROSCI.1516-13.2013, 2013.

[20] G. J. Yang, J. D. Murray, G. Repovs, M. W. Cole, A. Savic, M. F. Glasser, C. Pittenger, J. H. Krystal, X. J. Wang, G. D. Pearlson, D. C. Glahn, A. Anticevic, "Altered global brain signal in schizophrenia," Proceedings of the National Academy of Sciences of the United States of America, vol. 111 no. 20, pp. 7438-7443, DOI: 10.1073/pnas.1405289111, 2014.

[21] K. Linkenkaer-Hansen, V. V. Nikouline, J. M. Palva, R. J. IImoniemi, "Long-range temporal correlations and scaling behavior in human brain oscillations," The Journal of Neuroscience, vol. 21 no. 4, pp. 1370-1377, DOI: 10.1523/JNEUROSCI.21-04-01370.2001, 2001.

[22] K. Linkenkaer-Hansen, D. J. Smit, A. Barkil, T. E. van Beijsterveldt, A. B. Brussaard, D. I. Boomsma, A. van Ooyen, E. J. de Geus, "Genetic contributions to long-range temporal correlations in ongoing oscillations," The Journal of Neuroscience, vol. 27 no. 50, pp. 13882-13889, DOI: 10.1523/JNEUROSCI.3083-07.2007, 2007.

[23] V. V. Nikulin, T. Brismar, "Long-range temporal correlations in electroencephalographic oscillations: relation to topography, frequency band, age and gender," Neuroscience, vol. 130 no. 2, pp. 549-558, DOI: 10.1016/j.neuroscience.2004.10.007, 2005.

[24] K. Linkenkaer-Hansen, S. Monto, H. Rytsälä, K. Suominen, E. Isometsä, S. Kähkönen, "Breakdown of long-range temporal correlations in theta oscillations in patients with major depressive disorder," The Journal of Neuroscience, vol. 25 no. 44, pp. 10131-10137, DOI: 10.1523/JNEUROSCI.3244-05.2005, 2005.

[25] T. Montez, S.-S. Poil, B. F. Jones, I. Manshanden, J. P. A. Verbunt, B. W. van Dijk, A. B. Brussaard, A. van Ooyen, C. J. Stam, P. Scheltens, K. Linkenkaer-Hansen, "Altered temporal correlations in parietal alpha and prefrontal theta oscillations in early-stage Alzheimer disease," Proceedings of the National Academy of Sciences of the United States of America, vol. 106 no. 5, pp. 1614-1619, DOI: 10.1073/pnas.0811699106, 2009.

[26] Y. Lei, B. Song, L. Chen, J. Su, X. Zhang, W. Ni, Y. Yu, B. Xu, L. Yu, Y. Gu, Y. Mao, "Reconfigured functional network dynamics in adult moyamoya disease: a resting-state fMRI study," Brain Imaging and Behavior, vol. 14 no. 3, pp. 715-727, DOI: 10.1007/s11682-018-0009-8, 2020.

[27] P. Karzmark, P. D. Zeifert, T. E. Bell-Stephens, G. K. Steinberg, L. J. Dorfman, "Neurocognitive impairment in adults with moyamoya disease without stroke," Neurosurgery, vol. 70 no. 3, pp. 634-638, DOI: 10.1227/NEU.0b013e3182320d1a, 2012.

[28] Y. Lei, J. Su, H. Jiang, Q. Guo, W. Ni, H. Yang, Y. Gu, Y. Mao, "Aberrant regional homogeneity of resting-state executive control, default mode, and salience networks in adult patients with moyamoya disease," Brain Imaging and Behavior, vol. 11 no. 1, pp. 176-184, DOI: 10.1007/s11682-016-9518-5, 2017.

[29] C. G. Yan, X. D. Wang, X. N. Zuo, Y. F. Zang, "DPABI: data processing & analysis for (resting-state) brain imaging," Neuroinformatics, vol. 14 no. 3, pp. 339-351, DOI: 10.1007/s12021-016-9299-4, 2016.

[30] K. J. Friston, S. Williams, R. Howard, R. S. Frackowiak, R. Turner, "Movement-related effects in fMRI time-series," Magnetic Resonance in Medicine, vol. 35 no. 3, pp. 346-355, DOI: 10.1002/mrm.1910350312, 1996.

[31] C. G. Yan, B. Cheung, C. Kelly, S. Colcombe, R. C. Craddock, A. Di Martino, Q. Li, X. N. Zuo, F. X. Castellanos, M. P. Milham, "A comprehensive assessment of regional variation in the impact of head micromovements on functional connectomics," NeuroImage, vol. 76, pp. 183-201, DOI: 10.1016/j.neuroimage.2013.03.004, 2013.

[32] M. Jenkinson, P. Bannister, M. Brady, S. Smith, "Improved optimization for the robust and accurate linear registration and motion correction of brain images," NeuroImage, vol. 17 no. 2, pp. 825-841, DOI: 10.1006/nimg.2002.1132, 2002.

[33] C. Yan, Y. Zang, "DPARSF: a MATLAB toolbox for “pipeline” data analysis of resting-state fMRI," Frontiers in Systems Neuroscience, vol. 4, 2010.

[34] J. M. Beggs, D. Plenz, "Neuronal avalanches in neocortical circuits," The Journal of Neuroscience, vol. 23 no. 35, pp. 11167-11177, DOI: 10.1523/JNEUROSCI.23-35-11167.2003, 2003.

[35] A. Clauset, C. R. Shalizi, M. E. J. Newman, "Power-law distributions in empirical data," SIAM Review, vol. 51 no. 4, pp. 661-703, DOI: 10.1137/070710111, 2009.

[36] H. Yang, W. L. Shew, R. Roy, D. Plenz, "Maximal variability of phase synchrony in cortical networks with neuronal avalanches," The Journal of Neuroscience, vol. 32 no. 3, pp. 1061-1072, DOI: 10.1523/JNEUROSCI.2771-11.2012, 2012.

[37] C. K. Peng, S. V. Buldyrev, S. Havlin, M. Simons, H. E. Stanley, A. L. Goldberger, "Mosaic organization of DNA nucleotides," Physical Review. E, Statistical Physics, Plasmas, Fluids, and Related Interdisciplinary Topics, vol. 49 no. 2, pp. 1685-1689, DOI: 10.1103/PhysRevE.49.1685, 1994.

[38] R. Hardstone, S. S. Poil, G. Schiavone, R. Jansen, V. V. Nikulin, H. D. Mansvelder, K. Linkenkaer-Hansen, "Detrended fluctuation analysis: a scale-free view on neuronal oscillations," Frontiers in Physiology, vol. 3, 2012.

[39] B. B. Mandelbrot, J. R. Wallis, "Robustness of the rescaled range R/S in the measurement of noncyclic long run statistical dependence," Water Resources Research, vol. 5 no. 5, pp. 967-988, DOI: 10.1029/WR005i005p00967, 1969.

[40] X. Chen, B. Lu, C. G. Yan, "Reproducibility of R-fMRI metrics on the impact of different strategies for multiple comparison correction and sample sizes," Human Brain Mapping, vol. 39 no. 1, pp. 300-318, DOI: 10.1002/hbm.23843, 2018.

[41] J. D. Power, K. A. Barnes, A. Z. Snyder, B. L. Schlaggar, S. E. Petersen, "Spurious but systematic correlations in functional connectivity MRI networks arise from subject motion," NeuroImage, vol. 59 no. 3, pp. 2142-2154, DOI: 10.1016/j.neuroimage.2011.10.018, 2012.

[42] J. D. Power, K. A. Barnes, A. Z. Snyder, B. L. Schlaggar, S. E. Petersen, "Steps toward optimizing motion artifact removal in functional connectivity MRI; a reply to Carp," NeuroImage, vol. 76, pp. 439-441, DOI: 10.1016/j.neuroimage.2012.03.017, 2013.

[43] L. Koessler, L. Maillard, A. Benhadid, J. P. Vignal, J. Felblinger, H. Vespignani, M. Braun, "Automated cortical projection of EEG sensors: anatomical correlation via the international 10-10 system," NeuroImage, vol. 46 no. 1, pp. 64-72, DOI: 10.1016/j.neuroimage.2009.02.006, 2009.

[44] J. T. Kim, B. H. Cho, K. H. Choi, M. S. Park, B. J. Kim, J. M. Park, K. Kang, S. J. Lee, J. G. Kim, J. K. Cha, D. H. Kim, H. W. Nah, T. H. Park, S. S. Park, K. B. Lee, J. Lee, K. S. Hong, Y. J. Cho, H. K. Park, B. C. Lee, K. H. Yu, M. S. Oh, D. E. Kim, W. S. Ryu, J. C. Choi, J. H. Kwon, W. J. Kim, D. I. Shin, M. J. Yeo, S. I. Sohn, J. H. Hong, J. S. Lee, J. Lee, H. J. Bae, K. H. Cho, "Magnetic resonance imaging versus computed tomography angiography based selection for endovascular therapy in patients with acute ischemic stroke," Stroke, vol. 50 no. 2, pp. 365-372, DOI: 10.1161/STROKEAHA.118.023173, 2019.

[45] A. Roli, M. Villani, A. Filisetti, R. Serra, "Dynamical criticality: overview and open questions," Journal of Systems Science and Complexity, vol. 31 no. 3, pp. 647-663, DOI: 10.1007/s11424-017-6117-5, 2018.

[46] A. Eke, P. Herman, L. Kocsis, L. R. Kozak, "Fractal characterization of complexity in temporal physiological signals," Physiological Measurement, vol. 23 no. 1, pp. R1-38, DOI: 10.1088/0967-3334/23/1/201, 2002.

[47] A. Renart, J. de la Rocha, P. Bartho, L. Hollender, N. Parga, A. Reyes, K. D. Harris, "The asynchronous state in cortical circuits," Science, vol. 327 no. 5965, pp. 587-590, DOI: 10.1126/science.1179850, 2010.

[48] B. J. He, "Scale-free properties of the functional magnetic resonance imaging signal during rest and task," The Journal of Neuroscience, vol. 31 no. 39, pp. 13786-13795, DOI: 10.1523/JNEUROSCI.2111-11.2011, 2011.

[49] J. M. Palva, A. Zhigalov, J. Hirvonen, O. Korhonen, K. Linkenkaer-Hansen, S. Palva, "Neuronal long-range temporal correlations and avalanche dynamics are correlated with behavioral scaling laws," Proceedings of the National Academy of Sciences of the United States of America, vol. 110 no. 9, pp. 3585-3590, DOI: 10.1073/pnas.1216855110, 2013.

[50] D. Krzemiński, M. Kamiński, A. Marchewka, M. Bola, "Breakdown of long-range temporal correlations in brain oscillations during general anesthesia," NeuroImage, vol. 159, pp. 146-158, DOI: 10.1016/j.neuroimage.2017.07.047, 2017.

[51] S. Tu, J. Qiu, U. Martens, Q. Zhang, "Category-selective attention modulates unconscious processes in the middle occipital gyrus," Consciousness and Cognition, vol. 22 no. 2, pp. 479-485, DOI: 10.1016/j.concog.2013.02.007, 2013.

[52] K. Kazumata, K. Tokairin, M. Ito, H. Uchino, T. Sugiyama, M. Kawabori, T. Osanai, K. K. Tha, K. Houkin, "Combined structural and diffusion tensor imaging detection of ischemic injury in moyamoya disease: relation to disease advancement and cerebral hypoperfusion," Journal of Neurosurgery, vol. 3, 2020.

[53] K. Kazumata, K. K. Tha, H. Uchino, M. Ito, N. Nakayama, T. Abumiya, "Mapping altered brain connectivity and its clinical associations in adult moyamoya disease: a resting-state functional MRI study," PLoS One, vol. 12 no. 8, article e0182759,DOI: 10.1371/journal.pone.0182759, 2017.