Effective brain network analysis with resting-state EEG data: a comparison between heroin abstinent and non-addicted subjects

Objective. Neuro-electrophysiological tools have been widely used in heroin addiction studies. Previous studies indicated that chronic heroin abuse would result in abnormal functional organization of the brain, while few heroin addiction studies have applied the effective connectivity tool to analyze the brain functional system (BFS) alterations induced by heroin abuse. The present study aims to identify the abnormality of resting-state heroin abstinent BFS using source decomposition and effective connectivity tools. Approach. The resting-state electroencephalograph (EEG) signals were acquired from 15 male heroin abstinent (HA) subjects and 14 male non-addicted (NA) controls. Multivariate autoregressive models combined independent component analysis (MVARICA) was applied for blind source decomposition. Generalized partial directed coherence (GPDC) was applied for effective brain connectivity analysis. Effective brain networks of both HA and NA groups were constructed. The two groups of effective cortical networks were compared by the bootstrap method. Abnormal causal interactions between decomposed source regions were estimated in the 1–45 Hz frequency domain. Main results. This work suggested: (a) there were clear effective network alterations in heroin abstinent subject groups; (b) the parietal region was a dominant hub of the abnormally weaker causal pathways, and the left occipital region was a dominant hub of the abnormally stronger causal pathways. Significance. These findings provide direct evidence that chronic heroin abuse induces brain functional abnormalities. The potential value of combining effective connectivity analysis and brain source decomposition methods in exploring brain alterations of heroin addicts is also implied.

heroin abuse led to abnormal functional brain organization [2][3][4]. As functional neuroimaging technologies develop, it is possible to reveal the brain functional alternations caused by heroin addiction. Functional magnetic resonance imaging (fMRI) and EEG recordings have been widely applied to study Alzheimer's disease, depression, schizophrenia [5], as well as the abnormalities of brain functional system (BFS) induced by chronic heroin abuse [2,6]. fMRI data can provide accurate imaging of abnormal regions induced by heroin addiction. Studies implied that the gray matter density was reduced in the prefrontal and temporal cortices of the heroin addicts using resting-state fMRI data [7][8][9]; the regional homogeneities of the bilateral medial orbitofrontal cortex (OFC), bilateral dorsal medial thalamus, bilateral cuneus, and lingual gyrus were diminished in heroin addicts [1]. Nevertheless, the low temporal resolution of fMRI cannot provide sufficient brain information in the time-frequency domain.
EEG uses non-invasive electrodes to measure the electrophysiological signals accompanying real-time brain activities; it can be easily administered in either research or clinical settings [10]. Functional effective network based EEG is becoming one of the main tools to study the BFS of neurological disorders and cognitive functions. Babiloni et al used directional network methods to detect the abnormal brain functions in amnesic mild cognitive impairment (MCI) and Alzheimers disease (AD); they demonstrated that directionality of parieto-to-frontal EEG synchronization was abnormal in AD and amnesic MCI [11]. Winterer et al employed EEGcoherence to evaluate the frontotemporal structural connectivity and found that schizophrenic patients showed weaker frontal-temporal coherence [12]. The study by [13] applied directional methods to evaluate the BFS of Stroop cognitive tasks; they found the anterior cingulate cortex was modeled by motor areas, and dorsal lateral prefrontal cortex was modeled by Brodmann 9/46 regions. Studies of EEG based brain network also suggested BFS features changed during a cognitive task, as well as under sleeping state [14,15]. Nevertheless, there are still few EEG-based effective network studies on heroin addiction. EEG based heroin addiction studies mainly rely on power and coherence analysis tools [2,3].
EEG power and coherence analysis are not enough to reveal the brain's functional architecture and operational principles. One limitation is that the spatial resolution is low. Due to head volume conduction, electrical signals originating from one source in the brain are detected by multiple EEG electrodes. Conversely, each electrode measures a superposition of activity from multiple sources [18]. Simply classifying sensors into groups leads to misinterpretations of the synchronization results obtained at the sensor level [19,20].
Another limitation is that coherence analysis and functional connectivity only provide correlations between brain areas without directional information [21]. However, the synchronous neuronal activities in different brain regions are widely recognized as a mechanism of directional communication within the brain [19,22].
These two shortcomings can be overcome by brain source localization and effective connectivity methods. Source localization is the process of transferring signals from the sensor level into the source space [23]. Source inversing techniques and independent component analysis (ICA) are common approaches for source localization [18]. Source inversing depends on accurate models of head anatomy and electrical properties, as well as accurate electrode locations, while ICA performs a blind decomposition of EEG channels without a head model [18,24]. Among methods extended from the ICA approach, MVARICA is a methodology based on multivariate autoregressive (MVAR) modeling and ICA, which is good at determining the temporal activation of the intracerebral EEG sources as well as their approximate locations [18,20].
Effective connectivity developed from causality analysis theory explicitly reveals the information flow from one brain region to another; it involves both correlation and directional information between brain areas [21]. A number of research groups have applied various kinds of causal measure to estimate the effective brain connectivity on cognitive control, epileptics, lexical mediation studies [25][26][27]. Among the causal measures for estimating effective connectivity, GPDC has good performance in distinguishing between direct and indirect connections. It was also good at counteracting the impact of noisy data [28,29]. To date, few heroin addiction studies have focused on the effective brain connectivity induced by heroin abuse.

Subjects
A total of 15 heroin addicted (HA) subjects (age range 22-47 years) participated in this study. They were from the Second Compulsory Isolation Rehabilitation Center in Gansu Province, People's Republic of China. All subjects met the opioid dependence criteria in the diagnostic and statistical manual of mental disorders (DSM-V) (American Psychiatric Association, 2013), and had the history of taking heroin continuously for more than two years. Subjects with neurological and psychiatric disorders and associated genetic or family history were excluded; although having habits of smoking and drinking, the subjects selected did not have alcohol, nicotine or other drug addiction history. They did not receive methadone or alternative drug treatment, and had emotional stability and normal capacity of learning, motion, and perception, according to an interview conducted by the Compulsory Isolation Rehabilitation Center. The study was approved by the Ethics Committee of the institute of psychology of Chinese Academy of Sciences (Approval Number:H15020).
A total of 14 age-matched and education-matched nonaddicted (NA) participants were recruited from Linxia state in Gansu Province with normal intelligence. They were excluded from having mental illness or brain damage before enrollment. Participants were right-handed males, and signed the informed consent form before starting the study. Table 1 summarized the demographic information for both groups.

Preparation and EEG recording
All resting-state EEG data were acquired in the Second Compulsory Isolation Rehabilitation Center (heroin abstinent) and Linxia state (healthy controls) in Gansu Province. Before starting the experiment, the procedure was explained to the participants, who sat in a comfortable chair of a soundproof and dimly lit room. All participants were told to close their eyes, remain relaxed, and avoid movements. Data acquisition lasted 5 min.
A noninvasive scalp EEG was acquired using a 64-channel electroencephalographic recording analysis system (Brain Products, Germany) following the 10-20 electrode system of the international federation. Impedance of every participant was always below 20 kΩ. The sampling rate was 1000 Hz. The FCz electrode was used as the reference and EOG was recorded for later removal of eye artifacts. EEG data was recorded by vision recorder software version 2.0 from German Brain Products.

EEG analysis.
After checking the data quality by visual inspection, we found during the process of EEG signal recording, several subjects produced some actions, such as swallowing, coughing, and muscle shake, which resulted in large artifacts. To avoid polluting further data analysis, the corresponding bad segments were removed; for sample matching, EEG data lengths for all subjects were equalized. Finally, for each subject, 4 min lengths of data were cropped from 5 min. The ECG channel was unused in the experiment; 63 channels of EEG signals were used. The whole brain average was computed as a re-reference for all participants [31]. Since reliable source decomposition relied on high signal to noise ratio (SNR), we generated filtered datasets (1-45 Hz) for both subject groups. Since filtering might affect the original phases [32], we generated additional unfiltered datasets of both the HA and NA groups for effective connectivity analysis.
The EEGLAB toolbox in MATLAB was used to process and analyze continuous resting-state data [33]. For the preprocessing of filtered datasets, the data was filtered through a bandpass filter with low and high cutoffs of 1-45 Hz and a 50 Hz notch filter. For the unfiltered datasets preprocessing, only a 0.5 Hz high-pass filter and a 50 Hz notch filter were used to remove the voltage drift and linear noises. The filtered datasets and unfiltered datasets had the same artifact rejections processing for eye movements and environmental noise. Eye movements of sensor data were detected and removed by performing independent component analysis (ICA), after that the clean data was recomposed from the remaining components [34,35].
Individual de-noised EEG data had 4 min * 60 s * 1000 Hz = 240 000 samples. To save time and space in calculating resources, the sampling rate was down-sampled to 200 Hz, the long-term individual recordings were cropped into 24 segments with a time window of 10 sec, corresponding to a matrix (segments = 24, channels = 63, samples = 2000). All the individual recordings were arranged into HA and NA groups. This work assumed resting-state EEG data were similar within each group. Then, individual EEG data was normalized using a z-score method with corresponding group ensemble mean and standard deviation along the sample dimension. Subsequently, group matrices for further source decomposition and connectivity analysis were constructed. There were two kinds of matrices: filtered group and unfiltered group matrices. Each group matrix had the shape as (subjects * segments, channels, samples).

Brain source decomposition.
Filtered EEG data with high SNR was more suitable for brain sources estimation. SCoT toolbox implemented the blind source decomposition and VAR model fitting with the MVARICA approach [18,36]. Following the work of Billinger et al ( [18]), EEG signals were assumed as the linear mixture of source activations; the relationship was as defined in equation (1), the mixing matrix M transformed every sample n of the sources activations s n into the sample n of the filtered group EEG signals x n . The detail of source activations was as presented in equation (2), s n was modeled using previous p samples (also named as the model order) multiplied by corre sponding MVAR coefficient matrices A (k) and added with one innovation processes e, e was assumed to be an independent non-Gaussian noise processes of the fitted MVAR model.
The source decomposition worked in four steps (as shown in figure 1). Firstly, principal component analysis (PCA) was applied on EEG recording x n to reduce components with least contribution to total EEG variance and limit the number of sources in equation (3), in which C is the PCA transformation matrix. Secondly the MVAR model of y n was fitted using the previous samples k (from 1 to p) of principal components (PCs) multiplied by B (k) and added with white noise residuals r in equation (4) From equations (1)-(4), we got the relationship between two white-noise residuals in equation (5): Thirdly, a whiteness test was adopted to check whether residuals were non-white: when the p-value estimated from the whiteness test was higher than 0.05, the residuals of the fitted MVAR models were deemed white [18,37], otherwise the amount of PCs or previous samples had to be adjusted to fit the MVAR model well.
Fourthly, the residuals of the PC model were decomposed by ICA in order to obtain the estimated Ĉ M [38]. Finally, the unmixing matrix U was obtained from equation (6), EEG signals x were transferred into estimated source activations s in equation (7): The topography of the decomposed sources were estimated by combining s and EEG electrodes positions [18,20].
Previous heroin addiction studies indicated abnormality in temporal, frontal, parietal and occipital cortices [1,16]. In this work, we took four, six and eight as the candidates for the number of principal components (PCs); we expected the decomposed PCs were from the abnormal brain areas induced by heroin addiction. Based on the above procedure, the unmixing matrix of HA and NA groups was undertaken as follows:(1) PCA was applied to HA and NA group EEG signals separately (2) two VAR models were fitted using HA and NA group principal components; (3) the whiteness of the residuals of the two VAR models was tested; (4) the two white residuals of the two VAR models were concatenated and ICA run on the concatenated residuals; (5) the un-mixing matrix was generated, and variance sources between the two groups estimated.

Effective connectivity analysis.
Filtering methods following causality analysis might influence the reliability of the causality results [32], so unfiltered group EEG matrices were used for effective network modeling. As long as sources were spatially stationary, the un-mixing matrix U estimated from filtered EEG data could serve as a spatial filter for extracting the source signals from unfiltered group EEG matrices [18,39]. The stationariness of unfiltered sources could be checked by testing the constant of the mixing matrix M; the same model order p was used to construct MVAR model x = Ms, the eigenvalues of M were tested when the model order was varied from 1 to p. If all the eigenvalues were <1, the unfiltered sources were deemed suitable for the same un-mixing matrix U [18,40].
Applying the same un-mixing matrix U, data from two different groups could use the same source populations, which facilitated comparing brain networks of the two subject groups [18]. Once unfiltered sources were extracted from these decomposed regions, GPDC function was applied for Ultimately, this matrix along with the transformation C was used to obtain estimated s [18]. causality analysis between the source signals within each group. GPDC was a modified version of partial directed coherence (PDC). PDC provided a description of the mutual interaction between pairs of signal sources while deducting the effect of other simultaneous sources [41]. Compared with PDC, GPDC introduced weight of the noise variance of each involved time series, which made it suitable for distinguishing between direct and indirect connections. Meanwhile it was more robust for unfiltered time series [28,29].
As shown in figure 2, there were five steps for effective connectivity analysis of HA and NA groups: (1) unfiltered and z-scored EEG data (n_subjects*24, 63, 10 000) were decomposed into source matrices (n_subjects*24, PCs, 10 000) for each group, n_subjects was the number of subjects of each group, and PCs was the number of source regions or principal components; (2) then decomposed sources of each group were fitted to the MVAR model with the same model order p, the whiteness and stationariness of the fitted model were tested-if the tests were failed, the un-mixing matrix U would be re-estimated again using higher model order; (3) the fitted MVAR model was transformed into the frequency domain by Z-transform approach [28]; (4) GPDC was performed for each group, and effective networks of both the HA group and NA group were constructed. (5) Statistical connectivity differences of the two groups were obtained by performing bootstrap estimates [18].

Statistics.
In this work, we applied the shuffling phases strategy to estimate the significance of the causal interactions between decomposed source signals within each group. The bootstrap method was applied to reveal network differences between two group subjects [18]. The causal interactions were mainly inferred from the phase relationship between actual source signals, so surrogate data were derived by randomly shuffling phases without changing the magnitude of sources. The phase shuffling process destroyed the original causal interactions while keeping the spectral structure of the actual sources. Subsequently, GPDC was applied to the surrogate data. Corresponding to each group of effective networks, the shuffling and causality estimation procedures were repeated 1000 times, then a distribution of GPDC values was generated under the null hypothesis that no causal interactions existed between source regions. This work took the 99.99 percentile of the surrogate GPDC values as the significant threshold of the actual causal interactions [27,28,42,43].
Group connectivity differences were statistically evaluated by the bootstrap method. The present work applied the Bootstrap function from the SCoT toolbox for the abnormal effective connectivity estimates. Bootstrap samples were drawn at the trial level, for 29 subjects (HA: 15 subjects, NA: 14 subjects) * 24 trials = 696 trials in total. The distribution of differences in effective connectivity was obtained by bootstrapping both HA and NA groups. The repeats of Bootstrap estimates were 10 000. For each bootstrap estimate, 696 samples were generated from all of the 696 trials using the random sampling with replacement method. The ratio of falsely detected significant differences was guaranteed to be less than 0.01 [18].

Decomposed sources
Source decomposition and effective network modeling were based on MVAR models. Optimized model order was selected from 50 to 120. These values were used to construct MVAR models using the procedures outlined in figures 1 and 2 respectively. Each potential MVAR model was evaluated by whiteness and stable statistical tests, and model order selection mainly based on the whiteness and stable index. After testing model orders from 50 to 120, it was determined that model orders larger than 100 passed the whiteness test. Based on the rule that the minimum qualified model order should be taken as the optimized model order, p = 100 was selected for source decomposition and effective network modeling [44].
To reveal the effective connectivity differences between the HA group and NA group, brain sources reflecting the EEG variances of both groups were identified. We tested three PC number candidates: four, six and eight; three corresponding source region sets were obtained using the MVARICA approach. With the model order 100, all the fitted VAR models passed the whiteness and stability tests. As illustrated in figure 3, when EEG signals were decomposed into four or eight sources, frontal and temporal regions were both too redundant. However when EEG signals were decomposed into six source regions, the six decomposed regions were clearly different from each other. They were mainly consistent with previous heroin addiction studies [1,16], in that abnormal regions induced by heroin Under this procedure, causal interactions of each pair of sources were derived from unfiltered EEG data. Unfiltered stable EEG signals x were decomposed by the unified un-mixing matrix U. Subsequently estimated source signals s fitted a MVAR model using the model order p, w n was the additive white noise vector and A k was the MVAR model coefficients, stable and whiteness tests were performed on the constructed MVAR model. Then Z-transform was applied to transfer the stable VAR model into the frequency representation. Furthermore, effective networks of two groups were estimated by GPDC method, the causal influence from a source signal j to a target signal i was estimated by GPDC measure, (C ii is the ith diagonal element of the noise covariance matrix C, M is the number of decomposed sources, A ij ( f ) is the spectral transfer matrix from j to i). Finally, statistical connectivity differences between the two groups were estimated. addiction covered the frontal, temporal, parietal and occipital regions. So effective connectivity analysis would conduct between these six unified source regions.

Effective network construction
Under the same spatial MVARICA filter, unfiltered sources were extracted from the selected six regions. With the model order p = 100, the fitted MVAR models of HA and NA groups passed both the whiteness and stability tests. Corresponding to each group, GPDC was conducted on the fitted MVAR model, the significant causal thresholds were estimated using shuffling phase statistics. Since previous findings of the addictive electrophysiological activities were mainly in the range of 1-45 Hz, so this work only focused on GPDC spectral values within this frequency domain. As shown in figure 4, all the decomposed regions involved the resting-state brain network, and most of the spectral causal interactions were significant relative to the thresholds of surrogates. While from figure 4, it was hard to clarify the abnormal network features of the HA group, compared with the NA group.

Group network comparisons
To statistically check the significant differences of effective connectivity between HA and NA groups, the bootstrap estimate was adopted to estimate the abnormal effective connectivity in the HA group; this estimate was repeated 10 000 times. Figure 5 shows the estimates of significant network differences between two groups (p < 0.01), and in this work the false discovery rate (FDR) was applied to correct the multiple testing results. This figure indicated chronic heroin addiction induced alterations of the BFS. We further analyzed these differences to form table 2, outflows and inflows of each decomposed region are presented in table 3. As table 2 shows, effective network differences between the HA and NA groups were widely present between six unified regions, in the six well studied bands: δ, θ, α, β1, β2 and γ. HA < NA means the strength of causal influences in HA group was abnormally weaker than the NA group, and HA > NA means the strength of causal influences in the HA group was abnormally stronger than that in the NA group.
As shown in table 2, the number of weaker causal influences (green arrows) was larger than the number of stronger effective connectivity (red arrows) in the HA group, compared with the NA group; the parietal (P) region was a dominant hub in the abnormally weaker effective connectivity, while the left occipital (lO) region was a stable region in the abnormally stronger effective connectivity across all six frequency bands (summarized in table 3); the alterations of BFS induced by heroin addiction were presented as parts of causal interactions becoming stronger or weaker; in the low frequency domain (1-13 Hz), common weaker causal pathways of the HA group were rO → P ↔ rT, while there were no common stronger causal pathways in the HA group; in the high frequency domain , common weaker causal pathways of the HA group were P → rT, lO → lT and lO → F, while common stronger causal pathways of the HA group only had lO → rT; overall, the strength of the frontal-to-occipital circuits was broadly decreased in the HA group.

Discussion
This study presents the first effective connectivity analysis of the BFS alterations induced by heroin addiction based on EEG signals. It proposes one strategy to analyze brain networks of HA and NA subjects using the advanced tools MVARICA and GPDC on the resting-state EEG data. Differently from previous studies of heroin addiction [7-9, 16, 17], [45][46][47], this work studies the BFS alterations in the source space, and effective connectivity is estimated in the frequency domain.

Studies of EEG with MRI or functional MRI (fMRI)
Previous studies showed that combining EEG with MRI or fMRI helped in increasing the spatial resolution of electromagnetic source imaging, and meanwhile kept tracing rapid information interactions between source regions [6]. Alan et al revealed the relationship between cognitive tasks and related cortical activities in theta and alpha bands [48]. The review of [49] showed EEG-fMRI to be a powerful tool in detecting epilepsy dysfunctions. There were still few heroin addiction studies based on this combination [50]. EEG-fMRI/ MRI provided the opportunity to study the spatiotemporal  mechanisms of BFS. However during EEG recording, simultaneous MRI acquisition would pollute EEG signals due to the switching of the magnetic field gradients [51], and EEG-fMRI/MRI required complex processing methods to avoid these kinds of electromagnetical artifacts, and get accurate models of head anatomy and electrical properties [18]. The MVARICA approach was more convenient: as it relied only on EEG signals, it could inverse sensor data into the source space without relying on a head model or electrical properties. This work expected to adopt an easy-to-use approach to the study of BFS alterations induced by chronic heroin addiction.

Source regions related with heroin addiction
It has been widely accepted that chronic heroin abuse induces abnormalities in active brain regions [2,46]. Previous studies using fMRI data showed that bilateral medial orbitofrontal cortex (OFC), bilateral dorsal medial thalamus, bilateral cuneus, and lingual gyrus were diminished in heroin addicts compared with normal controls [1]. Since these accurate and deep brain regions are difficult to detect by EEG, there have been some EEG studies to verify their corresponding or nearby cortical surface regions: frontal, temporal, parietal and occipital regions [16,46,52].
EEG power spectral measurement is one tool to evaluate the abnormal brain regions of heroin addicts. Polunina et al studied the resting-state EEG power spectral properties of 33 heroin addicts with different lengths of addictive periods versus 13 healthy controls. They found that α2 (10-13 Hz) frequency shifts from frontal and central locations, especially in the right hemisphere, were associated with the length of addictive history. Heroin craving was also associated with α1 (8-10 Hz) frequency shifts at central, temporal and occipital regions [46,53]. The study of [53] also explored the relationship between Wechsler adult intelligence scale (WAIS) and resting-state EEG power spectral shifts; it indicated that mean frequencies of δ (1-4 Hz), θ2 (6)(7)(8), α (8)(9)(10)(11)(12)(13), and β (13-30 Hz) of right hemispheric frontal/temporal sensors were associated with cognitive dysfunctions of heroin abusers. The study by [16] found that occipital, right parietal, temporal and frontal regions had abnormal activities in chronic opioid addicts, compared with normal controls.
In this work, individual de-noised and filtered resting EEG data were normalized using the z-score method, and arranged into the HA/NA group data. Available MVARICA function from SCoT toolbox [18] was applied for source decomposition. The fitted models were checked by whiteness and stability tests in order to confirm the decomposed sources were reliable. By testing three candidates of source region numbers-four, six and eight-we found six decomposed source regions could present frontal, parietal, left/ right temporal, and left/right occipital regions clearly. The source regions identified in this work were mainly consistent with previous studies. Significant differences between the effective connectivity of HA and NA groups. Each subplot represented the causal spectral value within 1-45 Hz from its corresponding source in the column to its corresponding sink in the row. Gray histograms were the significant abnormal spectral causal fields in HA group (all p < 0.01, false discovery rate (FDR) corrected). Blue plots were the bootstrap estimates of the HA group connectivity, and green plots were the bootstrap estimates of the NA group connectivity. X axis represents the frequency domain 1-45 Hz, y axis presented the scale of GPDC values 0-1.0, and the diagonal represents the sources involving effective connectivity.
In this work, the abnormal relationship between frontal and temporal interhemispheric frequency bands were verified from the perspective of causality analysis with GPDC. The causal measure GPDC had been proved reliable, accurate and less affected by artifacts in [41,54]; for this reason, it was selected in this work to reveal directional interactions between source regions. Effective connectivity was conducted between six decomposed cortical regions, which covered the critical areas: frontal, temporal, parietal and occipital regions.
Studies in [2,46] pointed to abnormal α frequency shifts in the central, temporal, frontal and occipital regions of heroin addicts. In this study, we found complex abnormally weaker causal interactions in the α, these abnormal effective connectivity values might result in the abnormal α frequency shifts identified by other studies. Coullautvalera et al found drug addicts had higher coherence between frontal and posterior cortical regions in θ (4-8 Hz) band by investigating the resting-state EEG data of 21 drug abusers versus 20 normal controls [17]. From the perspective of effective connectivity, this work found the causal relationship between frontal to parietal regions in θ, it was P → F. There was also one stronger causal pathway from F to P in α; there was a bidirectional causal interactions between F and P in γ.
On the basis of verifying previous studies, new unique causal interactions were also found in the HA group. To observe the stable BFS alterations of HA, we summarized the causal interactions of the six frequency bands into two types: Comparative effective connectivity analysis between the HA and NA groups. In the 'HA < NA' column, green arrows represent abnormally weaker causal interactions of the HA group, and in the 'HA > NA' column, red arrows represent abnormally stronger causal interactions of the HA group (all p < 0.01, FDR corrected). The 'LF' row was the common weaker or stronger causal pathways across the low frequency domain (1-13 Hz), the 'HF' row shows the common weaker or stronger causal pathways across the high frequency domain   common causal interactions of the low frequency domain (1-13 Hz) and common causal interactions of the high frequency domain . Previous studies showed that the BFS were diminished in heroin addicts compared with normal controls [1]. We found common abnormally weaker causal interactions were more frequent than stronger causal interactions. Most of the causal interactions between abnormal regions were remote pathways, so it could also be explained by previous studies that the remote connectivity decreased in heroin addicts [16,45]. Additionally, by counting the outflows and inflows of each source region, we found parietal and left occipital regions were the two dominant hubs related with the BFS alterations induced by heroin addiction; these findings were supported by [17] and [55], which revealed the parietal and left occipital regions had abnormal activities related with HA subjects.
Since this work specially focused on the exploring of BFS alterations of heroin addicts, there were several limitations to be improved in the future work. First, the spatial resolution of brain regions was low-more advanced source decomposition approaches will need to be developed in the future. Furthermore, there are many other advanced tools for effective network study, such as nonparametric causality analysis and dynamical causality modeling, the introduction of which to them into heroin addiction research would be meaningful. Finally, only male addicts were included, so we could not determine whether these findings generalize to the female. We plan to include female subjects in future studies.

Conclusion
This paper has described a means to characterize the restingstate effective brain network of heroin abstinent individuals relative to non-addicted individuals using EEG data. To the best of our knowledge, this study reports the directional circuits abnormality of heroin addicts in the frequency domain, implying the frontal/temporal/parietal/occipital effective network, for the first time. It would be of great interest to further investigate how heroin addiction alters the effective brain network. Prospective studies were recommended to confirm the associations between these newly explored alterations of effective connectivity, and heroin addiction.