Brain lesions disrupting addiction map to a common human brain circuit

Drug addiction is a public health crisis for which new treatments are urgently needed. In rare cases, regional brain damage can lead to addiction remission. These cases may be used to identify therapeutic targets for neuromodulation. We analyzed two cohorts of patients addicted to smoking at the time of focal brain damage (cohort 1 n = 67; cohort 2 n = 62). Lesion locations were mapped to a brain atlas and the brain network functionally connected to each lesion location was computed using human connectome data (n = 1,000). Associations with addiction remission were identified. Generalizability was assessed using an independent cohort of patients with focal brain damage and alcohol addiction risk scores (n = 186). Specificity was assessed through comparison to 37 other neuropsychological variables. Lesions disrupting smoking addiction occurred in many different brain locations but were characterized by a specific pattern of brain connectivity. This pattern involved positive connectivity to the dorsal cingulate, lateral prefrontal cortex, and insula and negative connectivity to the medial prefrontal and temporal cortex. This circuit was reproducible across independent lesion cohorts, associated with reduced alcohol addiction risk, and specific to addiction metrics. Hubs that best matched the connectivity profile for addiction remission were the paracingulate gyrus, left frontal operculum, and medial fronto-polar cortex. We conclude that brain lesions disrupting addiction map to a specific human brain circuit and that hubs in this circuit provide testable targets for therapeutic neuromodulation.


Table of contents
Supplementary methods Figure S1. Smoker cohort lesions Figure S2. Lesion locations and voxel-lesion symptom mapping (VLSM) results. Figure S3. Addiction remission network maps are reproducible when computed in each dataset independently. Figure S4. Lesion network overlap results. Figure S5. Results are robust to different cut-offs for active daily smoking Figure S6. Addiction remission network results are robust to methodological choices. Figure S7. Addiction remission network results are robust to the connectome Figure S8. Mediation analysis Figure S9. Subsamples with neuropsychological assessment are representative for the full sample Figure S10. Structural connections related to addiction remission Figure S11. Generalizability of findings to other substance use disorders Figure S12. Brainsway coil electric field (E-field) models Table S1. Demographical data for smoker lesion datasets. Table S2. Significant clusters in the addiction remission network map. Table S3. Demographical data for the smoker connectome subjects Table S4. Neuropsychological test scores Table S5. White matter tracts associated with addiction remission Table S6. Demographical data for the Vietnam Head Injury Study Table S7. Regions whose connectivity profile best matches the addiction remission network map.

Supplementary appendix
Joutsa et al.

Overview:
We leveraged data from two existing cohorts in which information on smoking addiction was collected following a focal brain lesion (Table S1). 13,15 Prior publications using these cohorts tested whether lesions intersecting the insula (based on qualitative review of structural brain scans) were associated with disruption of addiction. 13,15 For the current study, we returned to the original structural brain scans to generate precise outlines of each lesion location in atlas space. We then used a new technique termed lesion network mapping to test a new hypothesis, namely whether lesions associated with disruption of addiction map to a specific brain circuit 18,19 . This analysis was not possible when the datasets were originally collected as it requires a wiring diagram of the human brain called the human connectome, which has only recently become available. 18,19 The human connectome is generated using resting state fMRI data from a large cohort of healthy volunteers (n=1,000). 43,44 In lesion network mapping, the connectome is used to approximate connectivity with the lesion location at the time of the lesion. As such, lesion network mapping does not require connectivity imaging of the patients themselves. This is an important advantage, as the technique can be used with previously collected lesion cohorts in which only structural brain scans are available. 18,19 It is worth noting that even if resting state fMRI data had been collected in these patients after the lesion, it would not be helpful for this analysis, as connectivity with lesion locations can't be computed after the lesioned tissue is lost.

Iowa Lesion Cohort
Active daily smoking was defined as smoking more than five cigarettes per day for more than two years at the time of lesion onset. 13 From the original 69 patients in the University of Iowa cohort, 13 two patients were excluded (one because the original scan was not found and one without a clearly identifiable focal brain lesion). Of these 67 patients, lesions from 47 patients were manually traced onto a template brain using the MRI or CT scan as a guide. 45,46 Because lesion tracing methods at the University of Iowa were changed in 2006, the other 20 patients had their lesions manually traced on native T1weighted scans with FSL and then transformed to the MNI152 atlas using nonlinear registration and lesion masking techniques available in ANTs software. 47 Lesion masking techniques included enantiomorphic normalization, which replaces the lesion volume with the voxel intensities from its non-damaged homologue, and applying a cost function to mask to the lesion volume for bilateral lesions. 48 The anatomical accuracy of all lesion tracings was reviewed in both native and MNI space and edited as needed by a neurologist (A.D.B.) blinded to clinical outcome data.

Rochester Lesion Cohort
Active daily smoking was defined as smoking at least one cigarette per day and at least 100 in their life-time at the time of lesion onset. 15 From the 156 patients in the University of Rochester cohort, 15 brain imaging data was available for 66 patients. Four patients were excluded because there was not a well-defined focal brain lesion (three had extensive hemorrhages and one lacked a clearly identifiable focal brain lesion). The remaining 62 brain lesions were manually traced on the MNI template based on the original brain MRI or CT by a neurologist (K.M.) blinded to clinical outcome data.
Both Cohorts: Disruption of smoking addiction was defined identically in both cohorts using previously established criteria: 13,15 the patient reported quitting smoking less than a day after their brain lesion, rated difficulty of quitting 1 or 2 on a scale of 1-7, reported not starting smoking again since their brain injury, and reported that they felt no urge to smoke since quitting. We will continue to use the term "smoking addiction" as it has been specifically defined for these cohorts, however updated terms including "nicotine dependence" or "substance/tobacco use disorder" can be considered analogous. For simplicity, we will also use the term "addiction remission" to refer to patients meeting the above criteria "disruption of smoking addiction".
Voxel-based lesion symptom mapping VLSM analyses were conducted using standard methods and NiiStat software (https://github.com/neurolabusc/NiiStat). 17 The analysis was conducted voxel-by-voxel across the whole brain. The contrast of interest was addiction remission vs. not quitting. Positive values in this VLSM map identify voxels more likely to be lesioned in addiction remission cases while negative values identify voxels more likely to be lesioned in non-quitters. Statistical significance was set at p<0.05 using

Supplementary appendix
Joutsa et al. Freedman-Lane permutation with 2000 permutations. Voxels not affected by any of the lesions were excluded from the analyses. Any significant findings were confirmed by repeating the analysis but excluding voxels affected in < 5% of the cases (n < 7). To assess for reproducibility across our two smoking datasets, the main analysis was repeated on each dataset individually, generating two independent VLSM maps for addiction remission (Fig S2C). We tested for similarity between these two maps using spatial correlation across all non-zero voxels, as in prior work from our group. Spatial correlation was selected to avoid any potential bias associated with arbitrary thresholds or weighting of certain regions, because it provides a simple metric of how much variance in one map matches another, and because if was used in our prior work. 49-51 This spatial correlation value was compared to a null distribution derived from randomly permuting clinical outcome and repeating the analysis (1000 iterations). We also tested whether clinical outcome in one dataset could be determined using the VLSM results from the other independent dataset. For this analysis, we computed the intersection of each lesion location with the VLSM map derived from the other dataset, summing the VLSM voxel values within the lesion mask. This generates a single number for each lesion that reflects the probability of addiction remission based on the VLSM results from the other dataset. Values for the addiction remission lesions were compared to the values from non-quitters using an independent sample twotailed t-test.

Lesion network mapping
To ensure that the results are independent of the subjects included in the connectome, lesion network mapping was conducted separately using two different connectomes: a connectome derived from current smokers (n=126) and connectome derived from healthy volunteers (n=1,000) using practically identical methodology. Detailed connectome robustness analyses comparing the two connectomes are described in "Robustness of smoker versus normative connectome to split-half replication".
Lesion network mapping was used to investigate functional connectivity between each lesion location and all other brain voxels using the same methods described extensively in prior work. 18,19,44 Briefly, lesion locations were used as seed regions for resting state functional connectivity MRI (rs-fcMRI) analysis using the two connectomes in separately. Pearson correlation maps of connectivity were transformed to a normal distribution using Fisher r-to-z transformation and averaged across all subjects in the connectome. This method identifies the network of brain regions functionally connected to each lesion location without the need for functional imaging data from the patients themselves.

Smoker connectome used for lesion network mapping
For the smoker connectome, participants were required to be right-handed, 18 to 60 years of age, free of active drug or alcohol abuse/dependence (other than nicotine), report no current psychiatric or neurological disorders, and have no contraindications for magnetic resonance imaging (MRI). Written informed consent was obtained in accordance with the National Institute on Drug Abuse Intramural Research Program Institutional Review Board. Data from a total of 126 participants recruited over four separate study protocols were included in the analyses. The demographics for the participants are summarized in Table S3. MRI data acquisition parameters: Whole-brain echo planar images were acquired on a 3T Siemens Allegra (n=81), Siemens Trio (n=35) and a Siemens Prisma scanner (n = 10). High-resolution oblique-axial T1weighted structural images were acquired using a 3D magnetization-prepared rapid gradient-echo (MPRAGE) sequence (repetition time (TR)=2000 ms (n=45), 2500 ms (n=81); echo time (TE)=3.5 ms (n=45), 4.38 ms (n=81); inversion time (TI) = 900 ms (n=45), 1000 ms (n=81), flip angle (FA)= 9° (n=45), 8° (n=81), voxel size=1 mm3). The resting state fMRI data were collected as oblique axial slices, using a T2*-weighted, single-shot gradient echo, echo planar imaging sequence sensitive to blood oxygenation level-dependent (BOLD) effects (TR=2000 ms; TE=27 ms; FA= 80° (n=69), 78° (n=45), 77° (n=12); slice thickness= 4 mm (n=75), 3.5 mm (n=51); field of view 220 × 220 mm2; image matrix 64 × 64; run length=180 frames/6 minutes). Structural data were processed using FreeSurfer version 6.0.1 to generate cortical surface models and subcortical segmentations that were then used for the registration and regression steps of the functional data preprocessing below. 52,53 Resting state functional data were processed with Dr. Thomas Yeo's Computational Brain Imaging Group's fMRI preprocessing pipeline (https://github.com/ThomasYeoLab/Standalone_CBIG_fMRI_Preproc2016) with the aid of wrapper scripts that allow the application of this pipeline to BIDS-formatted data (https://github.com/bchcohenlab/BIDS_to_CBIG_fMRI_Preproc2016) in similar fashion to the primary connectome data described in the main text. Preprocessing included the following steps: 1) removal of the first 4 frames/TRs of data, 2) slice time correction with FSL, 54,55 3) motion correction using FSL's mcflirt, and motion outliers were identified using fsl_motion_outliers at thresholds of framewise displacement > 0.2mm or DVARS > 50, with segments less than 5 frames also Supplementary appendix Joutsa et al.

4
being excluded 56 , 4) alignment of structural and functional images using boundary-based registration 57 from the FsFast software package (http://surfer.nmr.mgh.harvard.edu/fswiki/FsFast), and 5) nuisance removal via regression of whole brain/global signal, twelve motion correction parameters, averaged ventricular signal, averaged white matter signal, as well as their temporal derivatives (18 regressors in total). The flagged outlier volumes were ignored during the regression procedure. 6) The data were interpolated across censored frames using least squares spectral estimation of the values at censored frames, 58 and 7) finally, a band-pass filter (0.009 ≤ f ≤ 0.08 Hz) was applied. The preprocessed data was then projected to the 2mm version of the MNI152 NLIN 6 th generation atlas packaged with FSL using mri_vol2vol and smoothed with a FWHM of 6mm.
Functional connectivity maps were then created using AFNI (https://afni.nimh.nih.gov/ version 20.2.06): AFNI's 3dmaskave was used to extract the averaged time series from the residual files for each of the regions of interest (ROIs) for every participant. These time series were then used as seed regions in a whole-brain seed-based functional connectivity (FC) analysis using the 3dDeconvolve command. The correlation coefficients were then converted to z-scores using the Fisher's rto-Z transform. Finally, for each of the seed ROIs, a group t-statistic map versus zero was calculated over the 126 participants using AFNI's 3dttest++.

Normative connectome used for lesion network mapping
The normative resting state functional connectivity data used in this study is derived from the Brain Genomics Superstruct Project (GSP) dataset of 1,570 healthy young adult participants. 43 Versions of this dataset have been used for functional connectivity analysis for over a decade 59,60 and prior lesion network analyses from our group. 19,61 We have made a fully preprocessed version of this connectome publicly available [https://doi.org/10.7910/DVN/ILXIKS], along with all of the code and parameters used to process the data [https://doi.org/10.5281/zenodo.4905738].
The GSP data was acquired on matched Siemens 3T MAGNETOM Tim Trio MRI systems (Erlangen, Germany) at Harvard University and Massachusetts General Hospital using the vendor-supplied 12-channel phased-array head coil. Structural data included a high-resolution (1.2 mm isotropic) multi-echo T1-weighted magnetization-prepared gradient-echo image (multiecho MPRAGE). Functional imaging data were acquired using a gradient-echo echo-planar imaging (EPI) sequence sensitive to blood oxygenation level-dependent (BOLD) contrast. BOLD runs consisted of 47 interleaved slices (foot-head; 1, 3, 5 … 45, 47, 2, 4, 6 …, 44, 46). One hundred and twenty-four measurements were collected for each BOLD run (TR=3000 msec; 4 initial TRs collected to allow for T1-stablization and 120 valid measurements). During BOLD data collection participants were instructed to remain still, stay awake, and keep their eyes open while blinking normally. One or two BOLD runs were acquired per subject (72.6% of sessions included two runs). Sequences, parameters, and instructions were unchanged throughout the collection process. However, not all subjects were acquired on the same scanner, as five different scanners were used. In addition, during the scanning period, the scanner console changed from B13 to B15 to B17. The scanner (Scanner_Bin) and console version (Console) for each imaging session are available within the CSV files included in the original data release (GSP_list_140630.csv and GSP_retest_140630.csv). GSP subjects included in this connectome were chosen based on a combined 3-way score of normalized DVARS, normalized Entropy Focused Criterion (EFC), and normalized Framewise Displacement (FD) values generated from mriqc (https://mriqc.readthedocs.io/en/stable/) consistent with the goal of minimizing the effects of motion without "scrubbing". 62 The 1000 GSP subjects with the best combined scores were selected for connectome inclusion. Additionally, if a subject contains 2 resting-state BOLD runs, then the 2 motion quality values are averaged to produce a single value. Anatomical surfaces were created with FreeSurfer v4.5.
Functional preprocessing was then performed with a modified version of Thomas Yeo's Computational Brain Imaging Group (CBIG) fMRI preprocessing pipeline to generate the BOLD runs used in the connectome that allows for separate subcortical, cerebellar, and cortical spatial smoothing, as was done in, 59 but is not currently enabled in the latest versions of CBIG Pipeline. Our complete code and wrapper script to allow application of the CBIG pipeline to BIDS formatted data is publicly available at [https://doi.org/10.5281/zenodo.4905738]. The CBIG pipeline allows for a configurable sequence of steps for functional connectivity preprocessing and the exact parameters used here are available in our 'config file' online [https://github.com/bchcohenlab/BIDS_to_CBIG_fMRI_Preproc2016/blob/1.0.0/configs/Legacy_GSP.config]. Briefly, preprocessing included the following steps: 1) removal of the first 4 frames/TRs of data, 2) slice time correction with FSL, 54,55 3) motion correction using FSL's mcflirt, and motion outliers were identified using fsl_motion_outliers at thresholds of framewise displacement > 0.2mm or DVARS > 50, with segments less than 5 frames also being excluded 56 , 4) alignment of freesufer surfaces and structural and functional images using boundary-based registration 57 from the FsFast software package (http://surfer.nmr.mgh.harvard.edu/fswiki/FsFast), and 5) nuisance removal via regression of whole brain/global signal, twelve motion correction parameters, averaged ventricular signal, averaged white matter signal, as well as their temporal derivatives (18 regressors in total). The flagged outlier volumes were ignored during the regression procedure. 6) Supplementary appendix Joutsa et al.

5
Next, a band-pass filter (0.0001 ≤ f ≤ 0.08 Hz) was applied. 7) Finally, the preprocessed data was then projected to the 2mm version of the MNI152 NLIN 6 th generation atlas packaged with FSL using mri_vol2vol and smoothed with a FWHM of 7mm. Data from multiple runs, where available, was then concatenated.

Connectivity profile of addiction remission / addiction remission network
At each voxel, a general linear model was constructed to relate lesion connectivity to clinical outcome (addiction remission, quitting without remission, not quitting). Study site was included as a covariate and connections associated with addiction remission vs. not quitting were identified. The analysis was conducted across the entire brain and statistical significance was set at FWE-corrected p<0.05 using threshold-free cluster enhancement (tfce) implemented in FSL randomise. 34

Robustness of Smoker versus Normative Connectome to Split-Half Replication
Our primary results were nearly identical using the smoker versus the normative connectome. To guide our decision of which connectome to use in our many subsequent analyses (e.g. comparison with lesions associated with reduced alcoholism risk) we tested the robustness of each connectome to split-half replication. The robustness of the connectomes was assessed by splitting the connectome into two randomly selected halves (n=500 for the normative connectome and n=63 for the smoker connectome) and obtaining two randomly selected smaller samples from the normative connectome matching the number of subjects in split halves of the smoker connectome (n=63). We then computed the spatial correlation coefficient (r) between the resulting lesion connectivity maps within each connectome pair. The test-retest reliability within each connectome was then investigated by comparing the explained variance (r 2 ) between the connectomes (Fig. S7).

Reproducibility of the addiction remission network
Reproducibility was assessed by repeating the above analysis independently in the Iowa and Rochester lesion datasets, identical to the reproducibility testing performed for the VLSM analysis (Fig S3). The spatial correlation across non-zero voxels in each map was computed. This correlation value was compared to a null distribution derived from randomly permuting group labels (1000 iterations). To test for clinical predictive value, within each dataset, spatial correlation was computed between each subject's lesion connectivity and the lesion network map generated from the other dataset. These spatial correlations were compared between addiction remission vs. not quitting across both datasets using independent samples t-test (two-tailed). To demonstrate that our findings are specific for connectivity, the same reproducibility analyses were conducted using the lesion network maps and the VLSM maps.
Determining which group drives the addiction remission network To investigate which group drove the connectivity difference, we conducted two analyses. First, we determined the similarity of the connectivity maps within each group by computing the spatial correlation between each individual map and the average map for that group. The resulting correlation coefficients were compared between groups using Mann-Whitney Utest. Second, to investigate commonly connected regions within each group, we conducted a lesion network overlap analysis (Fig S4). Individual 'lesion networks' were created from the voxelwise connectivity maps, which were thresholded with |t| ≥ 5 corresponding to uncorrected P<10 -6 and whole brain voxel-level FWE-corrected P<0.05. The thresholded maps were binarized then overlaid to show brain regions connected to the greatest number of lesion locations. This analysis was performed separately for positive and negative connectivity and for lesions associated with addiction remission and not quitting. Higher lesion network overlap indicates greater consistency in the lesion networks.

Robustness of the addiction remission network to methodological choices
To ensure our results were robust to methodological choices, we repeated our analysis multiple times varying our inclusion criteria, covariates, cutoffs, and even our connectome. For all analyses, similarity to the original addiction remission network was assessed through spatial correlation.
1. First, we repeated our analysis varying our cutoff for classifying a subject as an "active daily smoker" (Fig S5). Note that this had been defined differently in the original studies of the Iowa and Rochester cohorts. To exclude the possibility that these different definitions impacted our results, we recomputed the addiction remission map after excluding eight patients

Supplementary appendix
Joutsa et al. 6 from the Rochester cohort who smoked less than five cigarettes per day (n = 121). We also repeated our analysis including only patients who smoked at least 20 cigarettes per day (n = 86).
2. Second, we repeated our analysis including lesion size (Fig S6B), proportion of gray/white matter in the lesion mask ( Fig  S6C) and age (Fig S6D) as a covariate. Although lesion size is usually included as a covariate in traditional lesion analyses, 17 it is not usually included in lesion network mapping analyses. 18,44,63-68 However, we sought to verify that the results were not influenced by lesion size and to make our lesion network mapping results directly comparable to the traditional VLSM results. Proportion of gray/white matter in the lesion masks was defined by extracting the number of voxels falling in the FSL gray and white matter tissue prior maps (threshold >0.5 for both) and calculated the difference of their proportions from the total number of voxels in the lesion, resulting in values from -1 (fully white matter lesion) to +1 (fully gray matter lesion).
3. Third, to verify that our results are not dependent on the lesion etiology or imaging modality, we replicated the main analysis in strokes only, which constituted 118/129 of our lesions (Fig S6E), and lesion masks defined with MRI only (Fig  S6F).
5. Fifth, to assess the effect of global signal regression on our results, we repeated our main analysis using a connectome processed without global signal regression. 69 Specifically, we processed the first 100 subjects of the 1,000 subject healthy volunteer connectome using aCompCor and the Conn Toolbox (www.nitrc.org/projects/conn) (Fig S7). 70

Mediation Analysis
We conducted a mediation analysis to test whether our connectivity results were mediated by the previously reported lesion overlap with the insula 13,15 (Fig S8). To maintain consistency with the original reports, we used categorical intersection with the insula (yes / no) as previously described and recorded for each lesion 13,15 . We next computed the similarity of each patient's lesion connectivity to our addiction remission map using spatial correlation, converted to a normal distribution using Fisher's r to z transform. This results in a single number for each lesion that reflects connectivity of the lesion location to our addiction remission network. Using the MATLAB Mediation Toolbox by Tor Wager (https://github.com/canlab/MediationToolbox), 71 we conducted a mediation analysis with lesion connectivity as the primary predictor, addiction remission as the primary outcome, and insula overlap as the mediator. We hypothesized that the effect of connectivity on addiction remission would not be mediated by insula overlap. For completeness, we also tested whether connectivity to our addiction remission network might mediate the previously reported association between insula damage and addiction remission.

Behavioral effects in the Iowa smoker cohort
After occurrence of the lesion, a subset of the patients in the Iowa cohort underwent neuropsychiatric testing, including Wechsler Adult Intelligence Scale III, Minnesota Multiphasic Personality Inventory and Trail Making Test Part B. To investigate if addiction remission was related to other neuropsychiatric changes, such as frontal lobe dysfunction or depression, we compared patients who remitted to patients who did not quit smoking using independent samples t-test (Table  S4). To ensure that these subsets were representative of the whole sample, we re-computed our addiction remission map using only the subset of patients that underwent any (n = 57) or all (n = 34) neuropsychological assessments. Similarity with the addiction remission map computed using these subsamples and the original addiction remission map computed using the full sample (n = 129) was assessed using spatial correlation (Fig S9).

Generalizability: Alcoholism risk scores
Brain lesions (n = 186) were collected from the Vietnam Head Injury Study 41 for patients that had completed the MacAndrew Alcoholism Scale (MAC) ( Table S6), 72 a subscale of the Minnesota Multiphasic Personality Inventory which measures personality characteristics that predispose to alcoholism. 73 Patients only completed this scale once, years after their focal brain injury. 41 All brain lesions had been previously traced based on CT scans and transformed to the MNI atlas, as described earlier. 69 Connectivity with each lesion location was computed as described above in "lesion network mapping." Connections correlated with MAC score were identified on a voxel-wise basis using Pearson's correlation. This process yielded a map of connections with the lesion location correlated with decreased risk of alcoholism. This correlation map was 7 converted to a normal distribution using Fisher's r to z transform. This alcoholism risk map was then compared with our addiction remission map using spatial correlation. To assess if the spatial similarity was higher than expected by chance, the spatial correlation coefficient was compared to a null distribution of spatial correlation coefficients calculated using permutation. 50 Specifically, our addiction remission map was re-generated with randomly permuted group labels (addiction remission versus non-quitters), our alcoholism risk map was regenerated with randomly permuted MAC scores, and the spatial correlation between the two maps was computed (1,000 iterations). The permutation method has been described in more detail earlier. 50 To determine if similar results could be obtained using lesion location alone (i.e. without using connectivity), we repeated this permutation analysis using VLSM maps rather than the lesion networks.
To ensure that the results were specific to reduced addiction risk, we repeated the analysis using the 10 main scales of the MMPI. 73 We again used spatial correlation and permutation to test whether the maps generated using each of these MMPI scales were more similar to our addiction remission network than would be expected by chance. As a direct test of specificity, we also compared the spatial correlation obtained using our alcoholism risk map to the spatial correlation obtained using these other 10 MMPI scales. The MMPI is designed such that these 10 scales are independent of one another, however we also orthogonalized the data to ensure that this was the case in our sample (results were similar with or without orthogonalization). Spatial correlations were converted to a normal distribution using Fisher's r to z transform. We then compared the spatial correlation with the addiction risk map to the distribution of spatial correlations with these other 10 maps using a one-sample two-tailed t-test.
To evaluate specificity beyond the domains measured by the MMPI, we analyzed 27 other behavioral variables available in this dataset including the Beck Depression Inventory (BDI), mini mental status exam (MMSE), and 25 subscales of the neurobehavioral rating scale (2 subscales were excluded because more than 90% of the values were zero). For each variable, we repeated the same analysis we performed for addiction risk, again assessing spatial similarity using permutation. We the qualitatively compared the spatial correlation coefficients for these 27 variables to the spatial correlation with alcoholism risk. Note that due to collinearity of these 27 measurements (e.g. depression scores from neurobehavioral rating score and BDI), we could not perform a direct statistical specificity analysis similar to the above analysis using the different MMPI scales (which are independent of one another).
Generalizability: Case reports of lesions resulting in addiction remission to other substances A literature search for case reports of lesions resulting in remission of substance use disorders was conducted in PubMed at the end of 2017 using the following searches: "substance use disorder" AND "lesion" AND "brain" which resulted in 321 hits, and "substance use disorder" AND "injury" AND "brain" which resulted in 1577 hits. Resulting entries were then screened for cases where brain lesions resulted in addiction remission (abstinence from drug seeking and craving) in patients with non-nicotine substance use disorders. This search identified three cases where an acute brain lesion resulted in remission of addiction to substances other than nicotine (Fig. S11). [74][75][76] The lesion location and connectivity in each case was defined as described above. The spatial correlation between these connectivity maps and the maps derived from smoking patients were computed and statistically compared between lesions resulting in addiction remission and not quitting using Mann-Whitney U-test.

Electric field models of the Brainsway H4 and H7 TMS Coils
The Brainsway TMS device using H4 coil was cleared by the US Food and Drug Administration (FDA) for short-term smoking cessation in late 2020. This coil was designed to target multiple brain regions. TMS acts by creating a rapidly changing electric field (E-field) in the brain via electromagnetic induction. 11,12 The E-field distribution in the brain is determined by several factors, including the coil shape and head anatomy, and can be computationally modelled. 12 The 3D Efield model of the FDA approved coils in the MNI space was obtained from the manufacturer (Brainsway, NJ, USA). The field maps were overlaid on the MNI brain and surface projection of the volumetric map is created with Mango (version 4.1, https://rii.uthscsa.edu/mango/) (Fig. S12). The E-field intensity is displayed as a percentage relative to the peak intensity in the brain.

Supplementary appendix
Joutsa et al. Lesions that resulted in addiction remission for smoking (A, n = 34). Lesions in patients who did not quit smoking (B, n = 69). Lesions in patients who quit but did not remit (C, n = 26). Each slice represents a different patient. There was no difference in lesion size between the groups (one-way ANOVA, F = 1.78, p = 0.17).

Fig. S6. Addiction remission network results are robust to methodological choices.
Our addiction remission network map from our primary analysis (A) remained unchanged when including lesion size (B), proportion of gray/white matter in the lesion (C) or patient age (D) as a covariate, excluding lesions other than strokes (E), or only including lesion patients with MRI (F). Maps are shown unthresholded to facilitate comparison. Spatial correlation coefficients between our original addiction remission network map (A) and the maps in B, C, D, E and F were r=0.99, r=0.99, r=0.999, r=0.98 and r=0.94, respectively.

Fig. S7. Addiction remission network results are robust to the connectome.
Our addiction remission network map from our primary analysis using an rs-fcMRI dataset from 1,000 healthy volunteers (A) remained unchanged compared to using rs-fcMRI dataset from 126 smokers with spatial correlation coefficient r=0.94 (B), or to healthy volunteer connectome processed without global signal regression (r=0.94, C). The significant clusters surviving correction for multiple comparisons with threshold free cluster enhancement (PFWE<0.05) were also almost identical (D-F). Test-retest reliability was significantly better for the normative connectome [mean (SD) r 2 0.98 (0.01)] compared to the smoker connectome [0.81 (0.10), two-sided paired t-test p<0.001] and normative connectome size-matched to the smoker connectome [0.90 (0.05), two-sided paired t-test p<0.001].

Fig. S12. Brainsway coil electric field (E-field) models.
E-field of the H4 coil for smoking cessation (A) and H7 coil for treatment of alcoholism. For clarity and to highlight the entire spectrum of the E-field intensities, a different color scale was chosen for this supplementary figure compared to Fig. 4 (bottom row). Note that the E-field models do not cover the back of the brain. However, the E-field intensities in the brain regions outside the E-field models are likely to be very low. A threshold of 20% of the maximum E-field is used for the illustration.

Supplementary appendix
Joutsa et al.  White matter tracts showing higher disconnection probability in patients with addiction remission compared to patients who did not quit smoking (Aspin-Welch's one-sided test, uncorrected p < 0.05). Family-wise error (FWE) correction for multiple comparisons was calculated using permutation analysis of linear models (PALM) analogous to the voxelwise analyses. Mean (SD) values of tract disconnection probability are shown. n.s. = not significant (p > 0.05). Table S7. Regions whose connectivity profile best matches the addiction remission network map.

Supplementary tables
The clusters are defined at r > 0.75. Coordinates are reported in MNI space (mm) and cluster sizes in voxels (2 x 2 x 2mm) for clusters with >100 contiguous voxels. The peaks with highest positive and negative r values are shown in Fig. 4. R = right, L = left, B = bilateral.