Associations between genetic risk, functional brain network organization and neuroticism

Neuroticism and genetic variation in the serotonin-transporter (SLC6A4) and catechol-O-methyltransferase (COMT) gene are risk factors for psychopathology. Alterations in the functional integration and segregation of neural circuits have recently been found in individuals scoring higher on neuroticism. The aim of the current study was to investigate how genetic risk factors impact functional network organization and whether genetic risk factors moderate the association between neuroticism and functional network organization. We applied graph theory analysis on resting-state fMRI data in a sample of 120 women selected based on their neuroticism score, and genotyped two polymorphisms: 5-HTTLPR (S-carriers and L-homozygotes) and COMT (rs4680-rs165599; COMT risk group and COMT non-risk group). For the 5-HTTLPR polymorphism, we found that subnetworks related to cognitive control show less connections with other subnetworks in S-carriers compared to L-homozygotes. The COMT polymorphism moderated the association between neuroticism and functional network organization. We found that neuroticism was associated with lower efficiency coefficients in visual and somatosensory-motor subnetworks in the COMT risk group compared to the COMT non-risk group. The findings of altered topology of specific subnetworks point to different cognitive-emotional processes that may be affected in relation to the genetic risk factors, concerning emotion regulation in S-carriers (5-HTTLPR) and emotional salience processing in COMT risk carriers. Electronic supplementary material The online version of this article (doi:10.1007/s11682-016-9626-2) contains supplementary material, which is available to authorized users.

For the SLC6A4, the 5-HTTLPR S/La/Lg variants were determined using PCR with Forward primer FAM-5'TGAATGCCAGCACCTAACCC-3' and Reverse primer 5-TTCTGGTGCCACCTAGACGC-3' (35 cycli of 30 seconds at 95°C, 30 seconds at 61°C and 1 minute at 72°C), and subsequent ingestion of the PCR product with the restriction enzyme Msp-I for at least 3 hours at 37 °C. The resulting restriction fragments were separated using capillary electrophoresis (ABI 3130 analyzer; Applied Biosystems, Nieuwerkerk a/d IJssel, the Netherlands), and fragment lengths were estimated using the ABI Prism® GeneMapper™ software, version 3.0 (Applied Biosystems). The La, Lg and S variants were determined by the detection of fragments of 325 base pairs (bp), 152 bp and 284 bp, respectively (validated inhouse method, (Doornbos et al. 2009) ).

Supplement 2: Overview of the full fMRI session
The full fMRI session consisted of four tasks, resting state and an anatomical scan. The following tasks/scans were presented in consecutive order: emotional face matching task (Hariri et al. 2002), mood (worry) induction paradigm (Paulesu et al. 2010), anatomical scan, resting state, interoceptive sensitivity task (Pollatos et al. 2007) and Ultimatum Game (Sanfey et al. 2003). The total duration of the fMRI session was approximately 60 minutes. The order was fixed and identical for all participants.

Supplement 3: Preprocessing steps
First, structural as well as functional images were reoriented parallel to the AC-PC plane.
Second, functional images were realigned to the first image using rigid body transformations and the mean EPI image, created during this step, was coregistered to the anatomical T1 image.
Third, structural images were corrected for bias field inhomogeneities, registered using linear transformations and segmented into grey matter (GM), white matter (WM) and cerebrospinal fluid (CSF) (MNI template space). Fourth, we used DARTEL (diffeomorphic anatomical registration through exponentiated lie algebra toolbox) (Ashburner 2007) to create a customized group template to increase the accuracy of inter-subject alignment. Individual GM and WM tissue segments were iteratively aligned to the group template in order to acquire individual deformation flow fields. Fifth, the coregistered functional images were normalized to MNI space using the customized group template and individual deformation flow fields. Furthermore, images were resampled to 2 mm 3 isotropic voxels and smoothed with an 8 mm full-width at halfmaximum (FWHM) Gaussian kernel.

Supplement 4: Scrubbing procedure
The indices framewise displacement (FD) and DVARS were calculated to indicate volumes (i.e. frames) that may be affected by motion artifacts (Power et al. 2012). FD is calculated as the root of the sum of the squared differentials per volume. Rotations were converted to translations assuming a distance of 65 mm from the origin of rotation (ArtRepair toolbox, http://cibsr.stanford.edu/tools/human-brain-project/artrepair-software.html). DVARS is calculated as the root mean square (RMS) of the derivatives of the time series across voxels included in the whole-brain mask per volume (Power et al. 2011(Power et al. , 2012. Volumes were removed when FD>0.5 mm and DVARS>mean + 3*SD. Additionally, one backward and two forward neighboring volumes were removed as well. The median of the number of scans that were removed per subject was 11.0 (IQR: 14.2). Subjects were excluded when more than one third of the volumes had to be removed.

Supplement 5: Module decomposition
A two-step procedure was applied to achieve the optimal modular structure using a threshold of 1.8% (see the next paragraph for details on the selection of this threshold). Input for this procedure was the binarized correlation matrix averaged across subjects. First, nodes were partitioned into modules using the algorithm of Blondel et al. (2008) (Blondel et al. 2008), wherein nodes are divided into groups with a maximum number of within-group edges and a minimum number of between-group edges. This calculation was repeated 500 times to increase the chance of escaping local maxima. The statistic was further optimized by applying the modularity fine-tuning algorithm of Sun et al. (2009) (Sun et al. 2009), wherein nodes are randomly assigned to other modules until modularity no further improves.

Supplement 6: Selection of the optimal threshold for module decomposition
First, correlation matrices were binarized using a range of threshold values (T=0.01-0.30, in increments of 0.01). Second, these matrices were averaged across subjects per threshold value and the entropy was calculated for each of them to indicate for which threshold value the edges showed the largest stability information-wise (lowest entropy). These results were compared to results obtained via randomized matrices (for details, see  ). The optimal threshold is the threshold where (i) the original matrix shows the largest stability across subjects (low entropy) and (ii) the difference in entropy is the largest between the original matrix and random matrix. The optimal threshold in the current study was 1.8%.   Table 2 Permutation results for the main effect of genetic group and the interaction between genetic group and neuroticism. For the main effect, the mean difference was calculated between the genetic risk and non-risk group per network measure for both polymorphisms. For the interaction effect, the difference in slope was calculated between the genetic risk and non-risk group for the association between neuroticism and a specific network measure for both polymorphisms. For the latter analyses, we only examined network measures that were related to neuroticism in our previous paper (Servaas et al. 2015). AS, affective subnetwork; COMT, catechol-O-methyltransferase; COS, cingulo-operculum subnetwork; DMS, default mode subnetwork; FPS, fronto-parietal subnetwork, SMS, somatosensory-motor subnetwork; VS, visual subnetwork. ** p-value 0.05, * p-value 00.10. Figure 1 For the main effect of the 5-HTTLPR polymorphism (participation coefficient DMS), density plots and boxplots are presented for several proportional thresholds (5%, 10%, 15%, 20%, 25% and 30%). We observed that differences were only pronounced for lower proportional thresholds (0.01-0.06). Note the different axes. DMS, default mode subnetwork; prop. thres., proportional threshold.