Introduction
Heroin addiction is a chronic disease characterized by compulsive drug seeking and use. In clinical practice, the maintenance of pharmacotherapy is a strong determinant of the outcome, but adherence to treatment remains problematic due to that relapse is common (Volkow, Jones, Einstein, & Wargo, Reference Volkow, Jones, Einstein and Wargo2019). A deep understanding of the mechanism of addiction is important and may provide an opportunity of seeking new pharmacotherapeutic targets.
By means of neurobiological imaging, the procedure from hedonic drug use to compulsive drug seeking has been well described and elaborated with a four-circuit addiction model composed of reward circuit, motivation/drive circuit, memory/learning circuit, and control circuit (Volkow, Fowler, & Wang, Reference Volkow, Fowler and Wang2003). The reward circuit mainly comprises the ventral tegmental area (VTA) and nucleus accumbens (NAc); the motivation/drive circuit comprises the orbitofrontal cortex (OFC), thalamus, caudate and putamen; the memory/learning circuit comprises hippocampus, amygdala; the control circuit comprises anterior cingulate cortex (ACC), inferior prefrontal cortex (IPC), dorsolateral and prefrontal cortex (dlPFC). It was proved that between the four circuits existed functional unbalance interplay, such as overactivating of reward circuit and motivation/drive circuit and lack of interaction of control circuit with reward circuit and motivation/drive circuit (Volkow et al., Reference Volkow, Fowler and Wang2003; Volkow, Wang, Tomasi, & Baler, Reference Volkow, Wang, Tomasi and Baler2013). Recently, studies supplemented that insular cortex was an essential structure for addiction behavioral maintenance representing the interoceptive effects of drug taking and making this information available to conscious awareness, memory and executive functions (Naqvi, Gaznick, Tranel, & Bechara, Reference Naqvi, Gaznick, Tranel and Bechara2014). Overall, the brain addiction-related topological organization model has been proposed and it, as a network comprising a number of modules made of nests sub-modules (circuits), is mediated by inter-modular connections of hubs.
Although the four-circuit model provides a relatively comprehensive neural framework for understanding the main neurobiological processes through which biological and sociocultural factors contribute to resilience against or vulnerability to drug use, it cannot completely explain the remarkable differences in the outcomes among addicts with pharmacotherapy. Methadone maintenance treatment (MMT), the most effective and accepted therapy for heroin dependence, for example, has situations reported like illicit drug use or retention rates altering widely among clinics (Gauthier, Eibl, & Marsh, Reference Gauthier, Eibl and Marsh2018; Zhou & Zhuang, Reference Zhou and Zhuang2014). Clinical neuropsychology studies demonstrated that multiple factors, e.g. age, gender, genetic variation, methadone dose, social experience, stress, trauma and environmental factors, may individually or combinedly contribute to the gap of MMT outcomes and further suggested that those factors could be used as predictors for MMT outcomes to a certain extent (Andersen, Reference Andersen2019; Coller, Barratt, Dahlen, Loennechen, & Somogyi, Reference Coller, Barratt, Dahlen, Loennechen and Somogyi2006; Darker, Ho, Kelly, Whiston, & Barry, Reference Darker, Ho, Kelly, Whiston and Barry2016; Jaremko, Sterling, & Van Bockstaele, Reference Jaremko, Sterling and Van Bockstaele2015; Lin, Hung, Peng, Chao, & Lee, Reference Lin, Hung, Peng, Chao and Lee2015; Lister, Brown, Greenwald, & Ledgerwood, Reference Lister, Brown, Greenwald and Ledgerwood2019; Luo et al., Reference Luo, Guo, Dai, Wang, Zhu, Miao and Gong2017; Proctor et al., Reference Proctor, Copeland, Kopak, Hoffmann, Herschman and Polukhina2015; Zhou, Li, Wei, Li, & Zhuang, Reference Zhou, Li, Wei, Li and Zhuang2017). However, few of the predictors or strategies was applied into clinical practice because of diversity, inconsistency and subjectivity of the studies' results.
Considering that the fundamental neural activity pattern underlies the relapse risk gap among MMT patients, scientists applied neuroimaging techniques to search the neuronal substrates associated with an eventual outcome to MMT. Li et al., studied the relationship between the percentage of positive urine drug screens and the functional connectivity (FC) strength between salience network, default mode network (DMN) and executive control network in MMT patients, and suggested the disrupted coupling between the salience network and DMN and between the left executive control network and DMN was associated with relapse behavior (Li et al., Reference Li, Liu, Wang, Wang, Li, Chen and Wang2018). With independent component analysis, Li and her colleagues compared the FC of DMN between the relapsers and abstainers during MMT, revealing the potential predictive value of DMN concerning heroin relapse under MMT (Li et al., Reference Li, Li, Wang, Xiao, Liu, Shi and Wang2015b). Some other neuroimaging indexes were also addressed for the association with relapse to heroin use in MMT, such as fractional anisotropy and axial diffusivity in the posterior limb of the internal capsule (Li et al., Reference Li, Zhu, Li, Ye, Chen, Liu and Wang2016c), regional homogeneity value in caudate (Chang et al., Reference Chang, Li, Li, Chen, Zhu, Ye and Wang2016), and cue-induced activation in NAc/subcallosal cortex (Li et al., Reference Li, Li, Wang, Wang, Zhang, Zhu and Liu2015a).
Although several relapse-related brain regions and neuroimaging indices were identified by studying the neuronal activity characteristics in the neuroimaging researches mentioned above, the practicality of these results was restricted from their retrospective evaluation design which means that participants' relapse outcomes were known in advance. Thus, a well-designed prospective cohort study using neuroimaging techniques at the condition of being completely blind to relapse outcomes was necessarily needed. Furthermore, an addiction-related brain network should be a suitable index of analysis for predicting patients' response to MMT, because the pathological alteration led by heroin use is wide at whole-brain level.
According to prior researches, the topological patterns in addiction-related networks were made of a set of submodules-neural circuits, consisting of some key nodes densely connected to each other and inter-modular connector hubs mentioned above. Degree centrality (DC), as topological features of graphs meaning the number of functional relationships of a given voxel with the rest of the network, can describe the addiction-related graphical models of the brain better (Buckner et al., Reference Buckner, Sepulcre, Talukdar, Krienen, Liu, Hedden and Johnson2009; Guo et al., Reference Guo, Liu, Hou, Wei, Liu and Chen2016; Luo et al., Reference Luo, Guo, Dai, Wang, Zhu, Miao and Gong2017; Wang, Jiao, Zhang, & Lin, Reference Wang, Jiao, Zhang and Lin2017; Zhou et al., Reference Zhou, Wang, Rao, Liang, Chen, Zheng and Li2014; Zuo et al., Reference Zuo, Ehmke, Mennes, Imperati, Castellanos, Sporns and Milham2012). And it was conceived as a stable and sensitive parameter in test−retest in the face of BOLD signal noise that was difficult to overcome by other methods (Zuo et al., Reference Zuo, Ehmke, Mennes, Imperati, Castellanos, Sporns and Milham2012). So it has been widely applied to examine node characteristics of intrinsic connectivity networks and the topology of brain hub changes with the disease at the global voxel level, for example, to investigate the association of alcohol dependence with changed intrinsic functional hubs and to reveal the attribution of several neural systems to the impaired motor behavior (Luo et al., Reference Luo, Guo, Dai, Wang, Zhu, Miao and Gong2017). Even DC was suggested as an index for attentional impulsivity in a recent study on codeine-containing cough syrups dependence (Hua et al., Reference Hua, Wang, Li, Li, Ma, Li and Jiang2018).
In this current study, based on the key nodes and inter-circuits brain hubs disclosed by prior studies on addiction and those using DC value, we constituted the addiction-related brain network of MMT patients and identified the relapse-related brain network organization pattern with cluster analysis. Our aim was to prospectively identify the poor MMT responders with neuroimaging data-driven method and reveal the brain regions contributing to relapse behavior. It was to be expected that the prospective identification of poor responders to MMT with cluster analysis based on neuroimaging data would be useful in modifying intervention strategies for improving the outcomes of MMT patients in future.
Methods
Participants
This study was approved by the Institution Board of the Fourth Military Medical University and conducted in accordance with the Declaration of Helsinki. Sixty-two patients diagnosed with heroin dependence receiving MMT (MMT patients) (male; mean age, 35.8 years; age range, 22–53 years) and 30 matched male healthy controls (HC) (mean age, 34.8 years; age range, 19–48 years) were recruited in this study. All participants were fully informed about the details of the experiment and signed the written consents for their involvement. Recruited from the outpatient of Xi'an Methadone Substitution Treatment Center, all the MMT patients met the DSM-IV criteria for heroin dependence, received MMT for no less than 6 months, and had at least 1 month of stable daily methadone dose before entering the study. Other inclusion criteria included: (1) right handed; (2) aged 18–55 years old; (3) heroin use history of more than 12 months. They were excluded if they had current or past psychiatric medical illness other than heroin dependence. HC were recruited from local community by advertising and they were free of DSM-IV-TR Axis I disorders and had no history of drug use. The exclusion criteria for all participants included: (1) current or past major medical illnesses or current use of prescription medications; (2) dependence other than psychoactive substances (except nicotine); (3) head trauma or neurological illness history; (4) MRI examination contraindication.
Evaluation measures
The demographics information (age, educational level), heroin use history (lifetime use and use dosage per day) and methadone treatment information (duration and methadone daily dosage) of each participant were collected during the enrollment interview. Before the MRI scan, all participants were asked to complete the questionnaire survey under psychiatrist guidance, including Beck depression inventory (BDI), Hamilton Anxiety Scale (HAMA). Furthermore, subjective heroin craving was assessed using a 0–10 visual analog scale, asking, “To what extent do you feel the urge to use heroin?” (0 indicating the least craving and 10 the strongest).
According to the reports from MMT fix-clinics in different cities of China, the MMT retention rate varied greatly (Cao et al., Reference Cao, Wu, Rou, Li, Lin and Wang2014; Jiang et al., Reference Jiang, Han, Du, Wu, Zhang, Zhang and Zhao2014; Zhang et al., Reference Zhang, Chow, Zhuang, Liang, Wang, Tang and Wilson2013). It was notable that studies on retention in MMT among opioid-dependent patients in Xi'an revealed the dropout and relapse mostly occurred between 12 and 60 months after MMT initiation (Wei et al., Reference Wei, Wang, Wang, Li, Li and Jia2013; Zhou & Zhuang, Reference Zhou and Zhuang2014). Thus, the follow-up duration in this study was set as 26 months.
After MRI scan, the MMT participants were followed for 26 months. During each appointment, the participant underwent a monthly structured interview assessing illicit drug use and a urine drug test. As opiates (mainly heroin, 38.1%) and synthetic drugs (mainly methamphetamine, 60.5%) were the main drugs of abuse in China (China Anti-Drug Network, annual report in 2017), relapse was defined as any use of heroin or methamphetamine identified by positive urine drug test (Morphine/Methamphetamine Diagnostic Kit, Guangzhou Jianlun Biological Technology Co., Ltd). In addition, the participants who missed the appointments without any contact with the study coordinator or who were loss of follow-up were considered to be relapsed. Finally, the accumulative total amount of positive urine tests and a missed interview was collected to calculate MMT relapse rate. Our coordinator contacted MMT participants 3 days before each appointment and encouraged them to continue the MMT program. There were 3 MMT patients lost in the follow-up period, two patients lost to follow up after MRI scan and one lost after the first clinical structural interview. There were 451 relapse counts totally, among them 216 were due to missing appointments and 139 due to loss to follow-up.
Brain imaging methodology and data analysis
Magnetic resonance imaging acquisition
MRI data were acquired with a 3.0T GE-Signa HDxt MRI scanner using an eight-channel head coil (GE Healthcare, Milwaukee, USA) in the Department of Radiology, Tangdu Hospital, the Fourth Military Medical University. Alcohol, tea, caffeine, any drug or medicine were prohibited in the 12 h before MRI scanning. All the MRI scanning were conducted 4 h after the MMT patients took their daily dosage (methadone peak plasma level) to avoid drug withdrawal symptoms (Wolff, Hay, Raistrick, & Calvert, Reference Wolff, Hay, Raistrick and Calvert1993).
Lying supine with the head fixed by a belt and foam pads, participants were instructed to keep their heads still, close their eyes and not to think anything specific. A routine structure MRI scanning was conducted to exclude gross cerebral pathology. After that, the resting-state fMRI (rs-fMRI) data were obtained using a gradient-echo planar imaging (GRE-EPI) pulse sequence, settings as: 30 axial slices, TR = 2000 ms, TE = 30 ms, flip angle = 90°, FOV = 256 mm × 256 mm, slice thickness = 5.0 mm, skip = 0 mm, matrix = 64 × 64. This session lasted for 5 min and 10 s, including a 10-s dummy scan period at the beginning of the scan. For each participant, 150 echo-planar volumes were collected during the rs-fMRI scan. And then a high-resolution structure session was performed with contiguous slices to cover the whole brain in a steady state using an axial fast spoiled gradient recalled echo (3D-FSGPR) for spatial normalization of the data sets to a standard atlas, scanning parameters: TR = 7.8 ms, TE = 3.0 ms, TI = 450 ms, FOV = 256 mm × 192 mm, slice thickness = 1.0 mm, skip = 0 mm, matrix = 256 × 256. After scanning, all participants reported they were awake during all scans.
Imaging data processing
Data preprocessing was conducted with SPM 8 software (http://www.fil.ion.ucl.ac.uk/spm) and DPABI (http://rfmri.org/dpabi). Prior to data preprocessing, the first five volumes were discarded to improve signal stabilization. Then, the remaining rs-fMRI data were corrected for slice timing and realigned for motion correction. The standard Montreal Neurological Institute (MNI) template was used for spatial normalization with a resampling voxel size of 3 mm × 3 mm × 3 mm. Participants with head motion more than 2 mm of translation and 2° of rotation in any of the x, y and z axes or mean frame-displacement more than 0.2 were excluded. Next, nuisance signals (including Friston 24-head motion parameters, the white matter and cerebrospinal fluid signals) were extracted and regressed out from the time series of every voxel to reduce the effects of nonneuronal signals. After a linear trend of the time course was removed, a band-pass filter (0.01–0.1 Hz) was applied. Then, scrubbing signal spikes based on the method of Power and colleagues was conducted (Power, Barnes, Snyder, Schlaggar, & Petersen, Reference Power, Barnes, Snyder, Schlaggar and Petersen2012). Specifically, a threshold of Frame-Wise Displacement = 0.5 mm was used and the bad volume and the volumes including 2 before and 1 after were removed. Two MMT participants and one HC were excluded from the next analysis because of excessive head motion or mean frame displacement.
Degree centrality (DC) calculation
Weighted DC measures were calculated using the “REST-DC” toolkit in the REST V1.8 packages (Zuo et al., Reference Zuo, Ehmke, Mennes, Imperati, Castellanos, Sporns and Milham2012) as previously described (Di Martino et al., Reference Di Martino, Zuo, Kelly, Grzadzinski, Mennes, Schvarcz and Milham2013; Liu et al., Reference Liu, Liu, Wei, Sun, Yang, Meng and Qiu2015). Briefly, to obtain each participant's graph, Pearson correlation coefficients were computed between the time series of all pairs of brain voxels within gray matter mask. Each voxel represented a node in the graph, and each significant functional connection (i.e. Pearson correlation) between any pair of voxels were the edges. As a result, a N × N undirected adjacency matrix (N was the count of all voxels within the gray matter mask) was obtained to construct the whole-brain FC matrix for each participant, each element of the matrix represents the Pearson correlation coefficient between pairs of voxels. To eliminate the possible spurious connectivity (Li et al., Reference Li, Ma, Huang, Li, Tian, Wen and Jiang2016b), the adjacency matrix was further thresholded with each Pearson correlation coefficient at r > 0.25. The weighted DC of a voxel was calculated as the sum of the connections' strength between a given brain voxel and all other voxels.
The formula is: $DC_i = \mathop {\sum {a_{ij}} }\limits_{j = 1}^N$
.
Where the DC i represents the weighted DC value for given voxel i, and a ij means the Pearson correlation coefficient between voxel i and voxel j.
Furthermore, standardized weighted DC (DC Z-score, zDC) maps were acquired by subtracting the mean value, and then divided by the standard deviation within the whole gray matter mask (Takeuchi et al., Reference Takeuchi, Taki, Nouchi, Sekiguchi, Hashizume, Sassa and Kawashima2015; Zuo et al., Reference Zuo, Ehmke, Mennes, Imperati, Castellanos, Sporns and Milham2012) according to the Z-score standardization formula below:
Where Z i represents the Z-scored DC value of voxel i, μ and σ are the mean and standard deviation of the DC measures across all N voxels. Finally, the resulting DC maps were spatially smoothed with a 4-mm FWHM Gaussian kernel.
Addiction-related network constitution
Based on the current addiction neurobiological theoretical framework, 10 pairs of brain nodes and hubs were chosen from “four circuits model” and addiction theory of insula and set as regions of interest (ROIs) to constitute brain addiction-related network for MMT patients. These regions included the bilateral NAc, amygdala, ACC, caudate, OFC, hippocampus, insular, putamen, thalamus, and DLPFC. The ROIs were determined with FSL https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/FSL (Jenkinson, Beckmann, Behrens, Woolrich, & Smith, Reference Jenkinson, Beckmann, Behrens, Woolrich and Smith2012; Smith et al., Reference Smith, Jenkinson, Woolrich, Beckmann, Behrens, Johansen-Berg and Matthews2004; Woolrich et al., Reference Woolrich, Jbabdi, Patenaude, Chappell, Makni, Behrens and Smith2009) according to Harvard–Oxford cortical and subcortical structural atlases based on MNI 152 1 mm template. The threshold probability was set as >50%. The DLPFC ROI was obtained by merging BA9 and BA46 of WFU PickAtlas http://fmri.wfubmc.edu/software/PickAtlas (Lancaster et al., Reference Lancaster, Woldorff, Parsons, Liotti, Freitas, Rainey and Fox2000; Maldjian, Laurienti, Kraft, & Burdette, Reference Maldjian, Laurienti, Kraft and Burdette2003). All ROI masks were resampled into 3 mm functional space. For characterizing the topological properties of brain networks, the DC value was extracted from the 10 pairs of ROIs of the 60 MMT patients after regressed out total methadone consumption, and used to make a 60 × 20 matrix M, representing the addiction-related network for subjects.
Cluster analysis and cluster validation and stability test
The clustering analysis was used to classify MMT patients' matrix M into subgroups with R (https://www.R-project.org/) and several R packages (factoextra, NbClust). The procedure included six primary steps: (1) assessing clustering tendency of the matrix M with the Hopkins' statistic and a visual approach; (2) identifying the optimal number K of clustering with Euclidean distance of ward D2 aggregation algorithm; (3) grouping the MMT participants according to K-means clustering algorithm; (4) clustering efficiency test with Silhouette analysis (Sun et al., Reference Sun, Lui, Yao, Deng, Xiao, Zhang and Gong2015); (5) clustering stability test with Fowlkes–Mallows (FM) index; (6) grouping stability test with the method provided by Hening C (online Supplementary-Materials).
Statistical analysis
To investigate the DC features coupling with subgroups, the differences in DC between MMT subgroups and HC were conducted with two-sample t test on the whole brain level (TFCE corrected p < 0.05 and cluster size K > 10). The difference in relapse rate between subgroups was calculated via Pearson's χ2 test with Yates' continuity correction. The differences in age, education level, smoking habit, heroin use information, methadone use history and behavioral score between MMT subgroups were calculated using two-sample t test or analysis of variance (>2 subgroups). The results were considered as significant if p value were less than 0.05.
Regression analysis
To explore which brain region contributed to relapse significantly, a sphere ROIs with 3 mm radius was made according to the coordination of peak voxel of brain region with significantly different inter-subgroup DC value. Robust Poisson regression (Tsou, Reference Tsou2006) was run to predict the relapse rate based on the z-scored DC signal strength for each ROI in 60 MMT participants, in which relapse count was set as dependent variable and only one brain DC value as independent variable in the model at each time. The statistical threshold was adjusted to 0.05/n after Bonferroni correction (n = the number of brain regions with significantly different inter-subgroup DC values).
Results
Clustering analysis
The results of clustering are shown as an ordered dissimilarity image (ODI) of 60 × 20 matrix M (Fig. 1). According to the index of clustering tendency, Hopkins statistic value for M was 0.38, far less than 0.5, indicating the clusterability of the matrix M. Besides, it was also implied by the ODI that the matrix M was of significant clusterability. The optimal clustering number K was determined as 2 according to Wards' algorithm and the calculations of 7 clustering indices out of 26. Clustering efficiency and stability tests also suggested the 2-cluster solution was the most stable and efficient one (online Supplementary-Materials). Finally, the 60 MMT participants were classified into two subgroups, 29 MMT patients (48.3%) in subgroup-1 and 31 (51.7%) in subgroup-2.
Ordered dissimilarity image of matrix M. The color level is proportional to the value of the dissimilarity between observations. Objects belonging to the same cluster are displayed in consecutive order. The dissimilarity matrix image confirmed that there is a cluster structure in the HD participants' data set. Two main subgroups (subgroup1 and subgroup2) were identified. Red: High similarity, Blue: Low similarity.
Differences in demographic data, clinical information between MMT subgroups
There was a significant difference in the relapse rate between subgroup-1 and subgroup-2 [χ2 = 22.578, df = 1, p = 2.018 × 10−06; ϕ = 0.13, 95% confidence interval (CI) = (0.07–0.017)], subgroup-1 higher than subgroup-2 (Fig. 2). Besides, subgroup-1 was significantly higher than subgroup-2 in the total heroin consumption [p = 0.047; Cohen's d = 0.58, 95% CI = (0.06–1.11)]. No significant differences existed in age, education level, methadone use information, BDI, HAMA, craving score, smoking duration, and daily cigarette number between the two subgroups (p > 0.05) (Table 1).
Violin plot of relapse rate of HD subgroups. HD subgroup1 had a significant higher relapse rate than subgroup2 (p < 0.05).
Demographic and clinical information of study subjects and subgroups
Note: HC, healthy control; HD, patients diagnosed with heroin dependence; S1, subgroup-1; S2, subgroup-2; NA, not applicable; s.d., standard deviation; *, statistically significant; HD-HC, HD v. HC, and so on.
Differences in DC value between MMT subgroups and HC group
As a result of circular analysis, the areas with significant differences in DC value between MMT subgroups and HC group reported were almost the ROIs used for cluster analysis. Compared with HC, MMT subgroup-1 had more and wider distribution of brain areas with changed DC, including the bilateral caudate, putamen, thalamus, left amygdala, hippocampus, vACC, and OFC with decreased DC, the bilateral paracentral gyrus and right postcentral gyrus with increased DC. On the contrary, subgroup-2 showed a less brain regions with altered DC, the bilateral calcarine and lingual gyrus with decreased DC value. When compared to subgroup-2, subgroup-1 had 17 brain regions with decreased DC, including the bilateral amygdala, caudate, OFC, thalamus, hippocampus, NAc, putamen, vACC, and left insular (TFCE corrected p < 0.05 and cluster size K > 10) (Fig. 3).
DC differences between HD and control and between subgroup 1 and subgroup 2 (TFCE corrected p < 0.05 and cluster size K > 10). Subgroup1 had a widespread altered density of connectivity involving in cortical–striatal–thalamic circuit at the whole-brain level, whereas a limited alteration in the subgroup2. As a result of circular analysis, the areas with significant differences in DC value between MMT subgroups and HC group reported were almost the ROIs used for cluster analysis (the upper, middle, and lower rows, respectively).
Regression analysis
The DC value of the bilateral vACC and NAc, hippocampus and right amygdala of MMT patients had a significant negative relationship with relapse rate. For every unit increased in the DC value of left vACC, right vACC, left NAc, right NAc, left hippocampus, right hippocampus and right amygdala, relapse rate decreased by 53.57, 58.11, 55.47, 43.41, 42.69, 50.89 and 41.4%, respectively (p = 6 × 10−06, Exp(B) = 0.4643; p = 0, Exp(B) = 0.4189; p = 1 × 10−06, Exp(B) = 0.4453; p = 3 × 10−06, Exp(B) = 0.4659; p = 0.001 Exp(B) = 0.5731; p = 5.3 × 10−05 Exp(B) = 0.4911; p = 2 × 10−04, Exp(B) = 0.526) (Fig. 4).
Estimates for signal of six brain ROIs in the models showed by Robust Poisson regression coefficient plot with 95% confidential intervals.
Discussion
This present study constructed the addiction-related brain network for MMT patients and at the condition of being completely blind to relapse rate in advance, we, by clustering the data set of the network, identified two MMT subgroups with different relapse-level defined by distinct patterns of addiction-related brain network abnormalities. The subsequent statistical analysis of clinical information confirmed that no differences were found in the demographic characteristics, methadone use or BDI and HAMA score between the two subgroups, and that subgroup-1 was the poor responder to MMT because it had a higher relapse rate. Here subgroup-1 was called as high relapse group (HRG) and subgroup-2 low relapse group (LRG). We think the prospective identification of HRG of MMT patients is the highlight of this study, because it may provide an opportunity for improved treatment tailoring and/or based on it, healthcare providers can adjust more aggressive adjunct therapies for the high-risk populations. In the following discussion part, we do not discuss the areas with significant DC value between the MMT subgroups and HC group to avoid emphasizing circular analysis.
Our study demonstrated that HRG had a widespread altered density of connectivity involving a cortical-striatal-thalamic circuit at the whole-brain level, whereas a limited one in the LRG. This finding implied that the brain network pattern of HRG identified is characterized by alteration of the importance of hubs extensively. Two factors may account for such global deficits.
The first one might be related to chronic heroin use. Heroin exposure cause brain structural injuries, including cell apoptosis (Tramullas, Martinez-Cue, & Hurle, Reference Tramullas, Martinez-Cue and Hurle2008), mitochondrial dysfunction (Feng et al., Reference Feng, Jia, Su, Wang, Lv, Xu and Yao2013), synaptic defects (Garcia-Fuster et al., Reference Garcia-Fuster, Ramos-Miguel, Rivero, La Harpe, Meana and Garcia-Sevilla2008), spongiform leukoencephalopathy (Pirompanich & Chankrachang, Reference Pirompanich and Chankrachang2015), disturbance of neurogenesis development (Bayer et al., Reference Bayer, Franke, Ficker, Richter, Lessig, Buttner and Weber2015), neuronal loss (Feng et al., Reference Feng, Jia, Su, Wang, Lv, Xu and Yao2013). In line with these findings, gray matter nucleus volume loss (Muller et al., Reference Muller, Truebner, Schiltz, Kuhn, Mawrin, Dobrowolny and Steiner2015), abnormal cortical thickness (Li et al., Reference Li, Tian, Zhang, Qiu, Wen, Ma and Huang2014) and white matter structural connectivity damage (Liu et al., Reference Liu, Li, Hao, Cao, Xu, Rohrbaugh and Liu2008; Wollman et al., Reference Wollman, Alhassoon, Stern, Hall, Rompogren, Kimmel and Perez-Figueroa2015) were discovered in morphological studies. Correspondingly, disrupted activities of neuronal clusters, brain regions, neuronal circuits and brain networks in heroin hijacked brain were revealed in functional imaging studies (Chang et al., Reference Chang, Li, Li, Chen, Zhu, Ye and Wang2016; Jiang et al., Reference Jiang, Qiu, Zhang, Han, Lv, Li and Tian2011; Li et al., Reference Li, Li, Li, Zhang, Wang, Zhu and Liu2016a; Wang et al., Reference Wang, Zhu, Li, Li, Wu, Zheng and Wang2013; Zhang et al., Reference Zhang, Li, Wen, Cai, Li, Tian and Wang2017). Actually, the pathological alterations resulted from heroin use is not equally distributed over the brain, but preferentially affect the hub regions (Volkow & Morales, Reference Volkow and Morales2015; Volkow, Wang, Fowler, & Tomasi, Reference Volkow, Wang, Fowler and Tomasi2012; Volkow, Wang, Fowler, Tomasi, & Telang, Reference Volkow, Wang, Fowler, Tomasi and Telang2011). Thus, the brain hubs with prominent reduced DC in HRG may reflect the neuropathological features associated with chronic heroin use.
From the view of neurodevelopment, the large-scale changed DC in hubs of high relapse risk network pattern could be interpreted as the consequence of complex interactions between genetic, epigenetic and environmental factors. Genetic factor accounts for approximately half of the risk for addiction (Volkow & Li, Reference Volkow and Li2005), and a series of gene variants are associated with addiction vulnerability by impacting on synaptic plasticity (Oliver et al., Reference Oliver, Brigman, Bolognani, Allan, Neisewander and Perrone-Bizzozero2018; Randesi et al., Reference Randesi, van den Brink, Levran, Yuferov, Blanken, van Ree and Kreek2018), receptor binding affinity (Shi et al., Reference Shi, Hui, Xu, Wang, Huang and Hu2002), dopamine receptor density (Noble, Reference Noble2000), neuron cell sensitivity to the substance of abuse (Alia-Klein et al., Reference Alia-Klein, Parvaz, Woicik, Konova, Maloney, Shumay and Goldstein2011), cortical neurotransmission and neurogenesis (Crews & Vetreno, Reference Crews and Vetreno2011), and gray matter volume shrinkage in chronic substance use (Gitik et al., Reference Gitik, Srivastava, Hodgkinson, Shen, Goldman and Meyerhoff2016). In addition, the epigenetic factors also contribute to the risk for addiction. People who experienced early life stress are at significantly high risk for the development of addiction caused by stress-induced neuro-plastic changes (Nestler, Reference Nestler2014). Other factors associated with addiction vulnerability include neurodevelopment imbalance between circuits supporting reward-seeking behaviors and self-control during adolescence (Giedd, Reference Giedd2008), early exposure to drugs of abuse (Nestler, Reference Nestler2014) and psychological traits of impulsivity, excitement seeking and stress reactivity (Bickel, Jarmolowicz, Mueller, Gatchalian, & McClure, Reference Bickel, Jarmolowicz, Mueller, Gatchalian and McClure2012; Holmes, Hollinshead, Roffman, Smoller, & Buckner, Reference Holmes, Hollinshead, Roffman, Smoller and Buckner2016; Jasinska, Stein, Kaiser, Naumer, & Yalachkov, Reference Jasinska, Stein, Kaiser, Naumer and Yalachkov2014). Consistent with the results, the longer heroin use duration in HRG indicated a younger heroin use age relative to LRG, which might imply the influence of gene and environment (Su et al., Reference Su, Li, Du, Jiang, Chen, Sun and Zhao2015).
The areas with significantly reduced DC in the primary and secondary visual cortex in the HRG might reflect the target of action of methadone. The previous study has demonstrated that methadone can depress visual function by acting on the visual parts, including midbrain and thalamic visual nuclei (Rothenberg, Peck, Schottenfeld, Betley, & Altman, Reference Rothenberg, Peck, Schottenfeld, Betley and Altman1979). Then, the relationship between the disrupted neural activity of occipital and long-term methadone use is confirmed by electroencephalograph techniques (McGlone et al., Reference McGlone, Mactier, Hamilton, Bradnam, Boulton, Borland and McCulloch2008; Wang, Kydd, & Russell, Reference Wang, Kydd and Russell2016). Although we did not observe the similar DC change of visual cortex in the LRG, we still think it is worth exploring the possible methadone-associated alteration of visual cortex function with imaging technology in the future. The increased DC value of paracentral and postcentral lobe in the HRG is inconsistent with those of prior opiate study, such as decreased glucose metabolism, reduced response to heroin-related cue and GO/NOGO task in paracentral and/or postcentral lobe (London et al., Reference London, Broussolle, Links, Wong, Cascella, Dannals and Rippetoe1990; Mei, Zhang, & Xiao, Reference Mei, Zhang and Xiao2010; Ye et al., Reference Ye, Li, Zhang, Li, Zhu, Chen and Wang2018). It might be a reflection of methadone playing its therapeutic role in addicted brain by increasing regional blood cerebral flow (Danos et al., Reference Danos, Kasper, Grunwald, Klemm, Krappel, Broich and Moller1998; Jiang et al., Reference Jiang, Qiu, Zhang, Han, Lv, Li and Tian2011). In an early study performed in heroin addicts, postcentral gyrus showed a reduced response to drug-related cue along with a decreased craving after buprenorphine dose (Mei et al., Reference Mei, Zhang and Xiao2010). Another study on positive therapeutic effect of acupuncture on alcohol use disorder revealed that several brain regions including postcentral gyrus were activated (Yang et al., Reference Yang, Choi, Kim, Ryu, Lim, Kim and Lee2017). Thereupon, we speculate the increased DC could be a reflection that treatment works. However, this phenomenon of increased DC of paracentral and postcentral region was not observed in LRG, which still calls for further investigations to understand the neuropathology underlying it.
Our regression analysis revealed that several hubs of addiction-related circuits, such as NAc, vACC, hippocampus and amygdala, were closely related to relapse. It raises the possibility of the potential therapeutic targets for addiction treatment.
NAc as a crucial hub in the network of reward, motivation and craving has rich and extensive interconnections with the VTA, amygdala, thalamus, prefrontal cortex, etc. Because of the key role in the process from voluntary drug use to compulsive drug-seeking behavior, NAc has been identified as a predictor for abstinence period following cession and a common target of pharmacological and cognitive-based interventions (Konova, Moeller, & Goldstein, Reference Konova, Moeller and Goldstein2013; Luigjes et al., Reference Luigjes, van den Brink, Feenstra, van den Munckhof, Schuurman, Schippers and Denys2012; Owens et al., Reference Owens, MacKillop, Gray, Beach, Stein, Niaura and Sweet2017). In an earlier clinical study, a prominent improvement of 5-year retention rate was observed in former heroin addicts with bilateral NAc ablation stereotactic surgery (Li et al., Reference Li, Wang, Wang, Chang, Ge, Gao and Gao2013). Furthermore, animal studies revealed that the possible way to attenuate drug reinstatement by manipulating NAc is local activation/or activation of GABAergic interneurons in the medial prefrontal cortex via antidromic stimulation of cortico-accumbal afferents (Vassoler et al., Reference Vassoler, White, Hopkins, Guercio, Espallergues, Berton and Pierce2013).
Hippocampus is involved in the formation of associations between drug use and special events and the regulation of reinstatement of drug-reinforced response (Francis, Chaudhury, & Lobo, Reference Francis, Chaudhury and Lobo2014; Koob & Volkow, Reference Koob and Volkow2010; Pascoli et al., Reference Pascoli, Terrier, Espallergues, Valjent, O'Connor and Luscher2014). Exposure to drug-associated contextual cue contributes significantly to relapse, which is mediated by drug-induced synaptic plasticity change in the ventral hippocampus (Alvandi, Bourmpoula, Homberg, & Fathollahi, Reference Alvandi, Bourmpoula, Homberg and Fathollahi2017; Borjkhani, Bahrami, & Janahmadi, Reference Borjkhani, Bahrami and Janahmadi2018). Therefore, there has been diverse and efficacy interventions targeting the hippocampus to prevent recurrent relapse in drug-dependent animal models. Wright, V. L. successfully used alpha7 nicotinic acetylcholine receptor antagonist to block the rats' reinstatement of morphine-conditioned place preference which was mediated by the retrieval of associative drug memories in ventral hippocampus (Wright et al., Reference Wright, Georgiou, Bailey, Heal, Bailey and Wonnacott2018). Another experiment proved the reinstatement induced by stress or combination of stress and primer dose of morphine was attenuated by antagonist administration targeting D1- and D2-like receptor in CA1 region of the hippocampus, addressing the role of hippocampal structures in treatment (Nazari-Serenjeh, Rezaee, Zarrabian, & Haghparast, Reference Nazari-Serenjeh, Rezaee, Zarrabian and Haghparast2018). In addition, several studies implicated that the drug-induced increase in synaptic plasticity in hippocampus projections to NAc shell played role in the contextual-mediated associations with drug-taking, and manipulation on hippocampus influenced drug reinstatement and drug-seeking behavior (Britt et al., Reference Britt, Benaliouad, McDevitt, Stuber, Wise and Bonci2012; Pascoli et al., Reference Pascoli, Terrier, Espallergues, Valjent, O'Connor and Luscher2014).
vACC is associated with relapse vulnerability and therapeutic outcomes (Forster, Dickey, & Forman, Reference Forster, Dickey and Forman2018). A study demonstrated that the abstinent time-related FC pattern of vACC and frontopolar cortex in a group of recently abstinent stimulant use disorder patients had a shift from strengthened FC first to a greatly reduced FC afterwards, implicating the patients' relapse risk (Camchong et al., Reference Camchong, Macdonald, Mueller, Nelson, Specker, Slaymaker and Lim2014). And it was found in cocaine dependents that their increasing impulsivity and higher relapse rate was associated with the aberrant FC between ventral caudate and subgenual ACC (Contreras-Rodriguez et al., Reference Contreras-Rodriguez, Albein-Urios, Perales, Martinez-Gonzalez, Vilar-Lopez, Fernandez-Serrano and Verdejo-Garcia2015). It was confirmed that the drug-induced neuroadaptation causing an imbalance between emotional processing and cognitive processing did also exist in the vACC (Shapira-Lichter et al., Reference Shapira-Lichter, Strauss, Oren, Gazit, Sammartino, Giacobbe and Lozano2018; Wang et al., Reference Wang, Wang, Qin, Yuan, Tian, Li and Guo2010). In view of this, it is reasonable to presume that vACC could be a potential medical interference target for addiction treatment.
Compelling evidence demonstrated that the amygdala, specifically the basolateral nucleus appears to play a critical role in relapse related to exposure to reminders of drug consumption, drug withdraw, fear or stress (Goode & Maren, Reference Goode and Maren2019; Koob, Reference Koob2009). The neuronal substrates underlying it is the afferents from the VTA to the basolateral nucleus that increase its activity and the reciprocal connections with the NAc involved with the information of seeking behavior (Ambroggi, Ishikawa, Fields, & Nicola, Reference Ambroggi, Ishikawa, Fields and Nicola2008; Ford, Mark, & Williams, Reference Ford, Mark and Williams2006; Grace & Rosenkranz, Reference Grace and Rosenkranz2002). As the potential mediator of the conditioned cue reinstatement, the amygdala has been chosen as the target of medication in the context of relapse. Disruption of neural activity in the basolateral nucleus or central nucleus of the amygdala in a rat model of relapse after chronic cocaine or ethanol self-administration blocks the initial acquisition of cue conditioning and the subsequent relapse, as well as withdrawal-induced drinking (Bale & Vale, Reference Bale and Vale2004; See, Fuchs, Ledford, & McLaughlin, Reference See, Fuchs, Ledford and McLaughlin2003). These findings are also strong evidence supporting the opinion about the negative correlation between amygdala and relapse rate proposed by our study.
Although we successfully divided the MMT patients into two subgroups with statistically different relapse rate according to the distinct addiction-related network patterns, the classification could not be as good enough as we had expected. Several MMT patients might be arranged into the wrong subgroup. We think the limitations of this study may contribute to this.
The first limitation should be originated from the selection strategy of ROIs for constructing Matrix M representing addiction-related network. We were mainly based on classical addiction neurobiological theories to determine the ROIs. Despite it is relatively simpler, the ROIs we selected, perhaps, could not be the best for constructing an addiction-related brain network. It is believed that the more specific the brain network reflecting neuronal activity associated with addiction, the better the clustering analysis results. There could be other methods to select the ROIs more sensitive to relapse, which deserves our more attemptation in future. The second one that might influence the clustering results is the sensitivity of MRI index describing ROIs neuronal activity. The reason we chose DC as the MRI index was that it had stability and sensitivity when as a topological feature of graphics and few MRI indicators were suitable for clustering analysis reportedly. As for whether DC is the most suitable for clustering analysis, it is still unclear and needs more investigations. To explore the impact of other MRI parameters on cluster results is valuable. The third one is the relatively longer interval of a structured interview and a urine drug test for assessing illicit drug use, which may lead to a bias of relapse rate. Moreover, we confess that the final results could be impacted by the limitation of only having male subjects. It is necessary to use a more appropriate interview plan and gender-matched group to verify the reliability and repeatability of our findings in future.
In summary, we prospectively classified the former heroin addicts with stable MMT into a high or low relapse risk group by grouping their brain neural activity patterns with a neuroimaging data-driven method. The two groups shared the same demographic information, methadone use history and smoking use except total heroin consumption and relapse rate. In addition, it was found that HRG was more extensive in the neuropathological changes in the brain, and that NAc, vACC, hippocampus and amygdala had a role in the contribution to illegal drug use. These brain hubs may be selected as potential psychological and pharmacotherapeutic targets in future. Importantly this new strategy of identifying HRG with a data-driven approach measuring addiction-related network pattern rather than clinical phenomenology would suggest a great possibility of medication pathway and benefit the MMT patients in the optimizing treatment plan and improving therapeutic effect.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/S0033291721003937
Acknowledgement
This work was supported by grants from the National Natural Science Foundation of China (82071497, 81471648 and 81771883) and the General Project (Youth project) of the Basic Natural Science Research Program of Shaanxi Province (2020JQ-883).
Conflict of interest
No conflict of interest was declared.
