Introduction
Major depressive disorder (MDD) is a common mental disorder and a leading cause of disability worldwide (GBD 2015 Disease and Injury Incidence and Prevalence Collaborators, 2016). Both genetic and environmental factors play important roles in the etiology of MDD, but underlying mechanisms remain poorly understood. In recent decades, considerable effort has been made to understand these mechanisms through genetic studies of MDD (Cohen-Woods et al., Reference Cohen-Woods, Craig and Mcguffin2013). However, early genome-wide association studies (GWAS), including a mega-analysis by Psychiatric Genomics Consortium (PGC), failed to show significant findings (Muglia et al., Reference Muglia, Tozzi, Galwey, Francks, Upmanyu, Kong, Antoniades, Domenici, Perry, Rothen, Vandeleur, Mooser, Waeber, Vollenweider, Preisig, Lucae, Muller-Myhsok, Holsboer, Middleton and Roses2010; Rietschel et al., Reference Rietschel, Mattheisen, Frank, Treutlein, Degenhardt, Breuer, Steffens, Mier, Esslinger, Walter, Kirsch, Erk, Schnell, Herms, Wichmann, Schreiber, Jockel, Strohmaier, Roeske, Haenisch, Gross, Hoefels, Lucae, Binder, Wienker, Schulze, Schmal, Zimmer, Juraeva, Brors, Bettecken, Meyer-Lindenberg, Muller-Myhsok, Maier, Nothen and Cichon2010; Terracciano et al., Reference Terracciano, Tanaka, Sutin, Sanna, Deiana, Lai, Uda, Schlessinger, Abecasis, Ferrucci and Costa2010; Shi et al., Reference Shi, Potash, Knowles, Weissman, Coryell, Scheftner, Lawson, Depaulo, Gejman, Sanders, Johnson, Adams, Chaudhury, Jancic, Evgrafov, Zvinyatskovskiy, Ertman, Gladis, Neimanas, Goodell, Hale, Ney, Verma, Mirel, Holmans and Levinson2011; Wray et al., Reference Wray, Pergadia, Blackwood, Penninx, Gordon, Nyholt, Ripke, Macintyre, Mcghee, Maclean, Smit, Hottenga, Willemsen, Middeldorp, de Geus, Lewis, Mcguffin, Hickie, van den Oord, Liu, Macgregor, Mcevoy, Byrne, Medland, Statham, Henders, Heath, Montgomery, Martin, Boomsma, Madden and Sullivan2012; Major Depressive Disorder Working Group of the Psychiatric GWAS Consortium, 2013). Some significant loci have been reported more recently with increased sample sizes (CONVERGE Consortium, 2015; Hyde et al., Reference Hyde, Nagle, Tian, Chen, Paciga, Wendland, Tung, Hinds, Perlis and Winslow2016; Power et al., Reference Power, Tansey, Buttenschon, Cohen-Woods, Bigdeli, Hall, Kutalik, Lee, Ripke, Steinberg, Teumer, Viktorin, Wray, Arolt, Baune, Boomsma, Borglum, Byrne, Castelao, Craddock, Craig, Dannlowski, Deary, Degenhardt, Forstner, Gordon, Grabe, Grove, Hamilton, Hayward, Heath, Hocking, Homuth, Hottenga, Kloiber, Krogh, Landen, Lang, Levinson, Lichtenstein, Lucae, Macintyre, Madden, Magnusson, Martin, Mcintosh, Middeldorp, Milaneschi, Montgomery, Mors, Muller-Myhsok, Nyholt, Oskarsson, Owen, Padmanabhan, Penninx, Pergadia, Porteous, Potash, Preisig, Rivera, Shi, Shyn, Sigurdsson, Smit, Smith, Stefansson, Stefansson, Strohmaier, Sullivan, Thomson, Thorgeirsson, Van der Auwera, Weissman, Breen and Lewis2017; Major Depressive Disorder Working Group of the Psychiatric Genomics Consortium, 2018) and in a recent study by the MDD working group of PGC, 44 loci were identified (Major Depressive Disorder Working Group of the Psychiatric Genomics Consortium, 2018). Despite a rapid increase in GWAS focusing on common variants (minor allele frequency, i.e. MAF >1 or 5%), studies on rare variants (usually MAF <1%) in MDD were lagged behind. Because it is estimated that only a quarter of the heritability of MDD is due to the common genetic variants (Cross-Disorder Group of the Psychiatric Genomics Consortium, 2013), other types of variants such as rare variants are likely to contribute to the genetic risk of the disorder. It is therefore recommended that genetics of MDD should be studied using ‘a balanced portfolio of methods’ (Sullivan et al., Reference Sullivan, Daly and O'DONOVAN2012).
When studying the rare variants, which tend to be more deleterious and have larger effect sizes than their common counterparts (Gibson, Reference Gibson2012), whole exome sequencing (WES) is the most commonly used technique. This has been applied to a number of psychiatric disorders, especially schizophrenia (Purcell et al., Reference Purcell, Moran, Fromer, Ruderfer, Solovieff, Roussos, O'DUSHLAINE, Chambert, Bergen, Kahler, Duncan, Stahl, Genovese, Fernandez, Collins, Komiyama, Choudhary, Magnusson, Banks, Shakir, Garimella, Fennell, Depristo, Grant, Haggarty, Gabriel, Scolnick, Lander, Hultman, Sullivan, Mccarroll and Sklar2014; Singh et al., Reference Singh, Walters, Johnstone, Curtis, Suvisaari, Torniainen, Rees, Iyegbe, Blackwood, Mcintosh, Kirov, Geschwind, Murray, Di Forti, Bramon, Gandal, Hultman, Sklar, Palotie, Sullivan, O'DONOVAN, Owen and Barrett2017) and autism (Tammimies et al., Reference Tammimies, Marshall, Walker, Kaur, Thiruvahindrapuram, Lionel, Yuen, Uddin, Roberts, Weksberg, Woodbury-Smith, Zwaigenbaum, Anagnostou, Wang, Wei, Howe, Gazzellone, Lau, Sung, Whitten, Vardy, Crosbie, Tsang, D'ABATE, Tong, Luscombe, Doyle, Carter, Szatmari, Stuckless, Merico, Stavropoulos, Scherer and Fernandez2015). By contrast, few studies with limited sample sizes have focused on effects of rare variants in MDD. An earlier study applied WES to 10 patients with an extreme pharmacological phenotype of MDD (Tammiste et al., Reference Tammiste, Jiang, Fischer, Magi, Krjutskov, Pettai, Esko, Li, Tansey, Carroll, Uher, Mcguffin, Vosa, Tsernikova, Saria, Ng, Eller, Vasar, Nutt, Maron, Wang and Metspalu2013). A more recent study sequenced the exome of 23 brain samples of suicide victims with MDD and 21 controls (Tombacz et al., Reference Tombacz, Maroti, Kalmar, Csabai, Balazs, Takahashi, Palkovits, Snyder and Boldogkoi2017). WES has been applied to two population cohorts to study depressive symptoms measured using the Epidemiological Studies Depression Scale in the inbreeding Erasmus Rucphen Family (Amin et al., Reference Amin, Belonogova, Jovanova, Brouwer, van Rooij, van den Hout, Svishcheva, Kraaij, Zorkoltseva, Kirichenko, Hofman, Uitterlinden, van Ijcken, Tiemeier, Axenovich and van Duijn2017a) and the outbreeding Rotterdam Study (Amin et al., Reference Amin, Jovanova, Adams, Dehghan, Kavousi, Vernooij, Peeters, de Vrij, van der Lee, van Rooij, van Leeuwen, Chaker, Demirkan, Hofman, Brouwer, Kraaij, Willems, Hankemeier, van Ijcken, Uitterlinden, Niessen, Franco, Kushner, Ikram, Tiemeier and van Duijn2017b). Despite reports that the genetic architecture of common variants in lifetime MDD overlap with current depressive symptoms in the general population (Major Depressive Disorder Working Group of the Psychiatric Genomics Consortium, 2018), it remains unclear whether this is the case for rare variants.
Given the considerable heterogeneity of the psychiatric disorders, researchers have taken advantage of endophenotypes, measurable components along the pathway between disease and distal genotype (Gottesman and Gould, Reference Gottesman and Gould2003) in genetic studies. Gray matter volume (GMV) is one of the reasonable endophenotypes for MDD due to the following reasons: (1) it is heritable to a large extent (60–80%) (Jansen et al., Reference Jansen, Mous, White, Posthuma and Polderman2015); (2) GMV abnormality, including reduction of GMV in the prefrontal cortex, temporal cortex, limbic area has been mapped in MDD (Lener and Iosifescu, Reference Lener and Iosifescu2015); (3) GMV of some regions are trait-related in MDD (Zeng et al., Reference Zeng, Shen, Liu, Fang, Liu and Hu2015); and (4) some brain regions, e.g. the insular, have reduced GMV, but to a lesser extent when compared with patients, in unaffected relatives (Matsubara et al., Reference Matsubara, Matsuo, Harada, Nakano, Nakashima, Watanuki, Egashira, Furukawa, Matsunaga and Watanabe2016). By integrating both genetic and imaging components in parallel into a multivariate model using the parallel Independent Component Analysis (pICA) approach (Liu et al., Reference Liu, Pearlson, Windemuth, Ruano, Perrone-Bizzozero and Calhoun2009), we can simultaneously investigate imaging–genetic associations and identify genetic factors underpinning GMV variation of MDD. In the multivariate approach, pICA captures clusters of variants or specific brain features co-varied across individuals. Previous studies of schizophrenia (Liu et al., Reference Liu, Pearlson, Windemuth, Ruano, Perrone-Bizzozero and Calhoun2009), bipolar disorder (Meda et al., Reference Meda, Ruano, Windemuth, O'NEIL, Berwise, Dunn, Boccaccio, Narayanan, Kocherla, Sprooten, Keshavan, Tamminga, Sweeney, Clementz, Calhoun and Pearlson2014), and attention-deficit hyperactivity disorder (ADHD) (Khadka et al., Reference Khadka, Pearlson, Calhoun, Liu, Gelernter, Bessette and Stevens2016) have shown encouraging results using pICA.
Based on this knowledge, the current study was designed to explore both the effect of rare variants on MDD and investigate pathways underlying GMV changes in MDD. We firstly targeted on rare coding variants and performed burden tests on different levels of variant sets, adopting a hypothesis-free strategy to avoid artificial bias. pICA was then used to explore correlations between the imaging and genetic modalities, where GMV of brain regions were included as potential endophenotypes of MDD. Based on the genetic–imaging correlation, we characterized genes and pathways enriched in the genetic component that might relate to the GMV changes in MDD. Finally, we investigated the relationships between clinical symptoms and cognitive assessment and the identified genetic and imaging components.
Materials and methods
Subjects
Patients with MDD (n = 77) who had not taken antidepressant in the past 3 months were recruited from the Mental Health Center in West China Hospital of Sichuan University. They were interviewed and assessed using the Structured Clinical Interview for DSM-IV Patient Edition (SCID-P) (First et al., Reference First, Spitzer, Williams and Gibbons1995a) by a trained psychiatrist. Diagnoses were assigned according to the criteria for MDD in DSM-IV. Healthy controls (n = 250) were recruited from the local community and screened using the SCID Non-patient Edition (SCID-NP) (First et al., Reference First, Spitzer, Williams and Gibbons1995b). All participants were Han Chinese and right-handed. More details of inclusion and exclusion criteria are listed in online Supplementary material. Demographic and psychopathological data are given in Table 1.
Demographic and psychopathological information of patients with MDD and HCs
Mean ± s.d. (standard deviation) or number of subjects are displayed in the table; MDD, major depressive disorder; HCs, healthy controls; HAMD, Hamilton Depression Rating Scale; IQ, intelligence quotient.
a The p value was obtained by χ2 test.
b The p value was obtained by two sample t test.
Clinical and cognitive data acquisition
The 17-item Hamilton Depression Rating Scale (HAMD-17) was used to access the severity of depressive symptoms for all patients. Only patients with a total score ⩾18 were included. We also assessed general cognitive functioning using the short version of the Wechsler Adult Intelligence Scale Revised in China (WAIS-RC) with seven subtests including information, arithmetic, digital symbol, digital span test, block design, picture completion, and similarities. Verbal and performance intelligence quotient (IQ) scores were calculated according to the WAIS-RC manual (Ryan et al., Reference Ryan, Weilage and Spaulding1999) and the detail was described in online Supplementary material.
Genetic data acquisition and processing
Sequence processing
DNA was collected from peripheral blood and extracted using Illumina Paired-End Sample Preparation kits. Agilent SureSelect Exome Enrichment Kit V5 was used to capture the exonic DNA for all cases (n = 77) and a subgroup of the controls (n = 97), while Illumina TrueSeq Exome Enrichment Kit was used for the remaining controls (n = 153). To check the consistency of sequencing in two batches with different capture kits, we sequenced two controls twice using different capture kits. Paired-end sequencing was performed on the Illumina HiSeq 2000. GATK Best Practices were adopted for data processing (online Supplementary material) (McKenna et al., Reference Mckenna, Hanna, Banks, Sivachenko, Cibulskis, Kernytsky, Garimella, Altshuler, Gabriel, Daly and Depristo2010). Rare variants in identified candidate genes and pathways were validated by Sanger sequencing.
Quality control
FastQC was used to check the overall quality of raw sequencing reads and then VeryfyBamID (Jun et al., Reference Jun, Flickinger, Hetrick, Romm, Doheny, Abecasis, Boehnke and Kang2012) was adopted to detect contamination for all subjects. We also used PLINK (Purcell et al., Reference Purcell, Neale, Todd-Brown, Thomas, Ferreira, Bender, Maller, Sklar, de Bakker, Daly and Sham2007) to check the sequencing data for population stratification, subject relatedness, and sex consistency. GATK (McKenna et al., Reference Mckenna, Hanna, Banks, Sivachenko, Cibulskis, Kernytsky, Garimella, Altshuler, Gabriel, Daly and Depristo2010) and KGGSeq (Li et al., Reference Li, Li, Li, Pan, Hsu, Liu, Zhan, Wang, Song and Sham2017) were then used for quality control (QC) checking at the variant and genotype level. The specific parameters are presented in online Supplementary material.
Annotation
The variants were annotated using KGGSeq (Li et al., Reference Li, Li, Li, Pan, Hsu, Liu, Zhan, Wang, Song and Sham2017), where they were mapped to 17 feature groups in accordance with the Human Genome Variation Society recommendations. Only non-synonymous exonic and splicing variants were considered in our analysis. MAF of the variants were retrieved from three public databases, namely Exome Aggregation Consortium (ExAC r0.3.1), East Asian panel of 1000 Genomes Project released in 2013 May (1kgeas), and South Asian panel of 1000 Genomes Project released in 2013 May (1kgsas). We then predicted the potentially damaging consequences for each variant by Combined Annotation Dependent Depletion (CADD) score (Kircher et al., Reference Kircher, Witten, Jain, O'ROAK, Cooper and Shendure2014).
MRI data acquisition and processing
Overall, 75 cases and 106 controls were scanned on a 3.0 T MR scanner (Achieva; Philips, Amsterdam, the Netherlands) using an eight-channel phased-array head coil in the department of radiology. T1 images were processed using Diffeomorphic Anatomical Registration Through Exponentiated Lie algebra toolbox in Statistical Parametric Mapping 8, and followed by the standard approach of Ashburner (Ashburner, Reference Ashburner2007). Detailed procedures are presented in online Supplementary material.
Statistical analyses of rare variants
We carried out set-based Multivariate and Collapsing (CMC) burden tests adjusting for age, sex, and the first two principal components (PCs) using RvTest (Zhan et al., Reference Zhan, Hu, Li, Abecasis and Liu2016). CMC collapsed rare variants into a single unit to increase power. We carried out set-based tests under two frequency categories (MAF ⩽1% and MAF ⩽0.1%) and two annotation categories (non-synonymous variants or potentially damaging variants), for a total of four tests. Any gene with at least three variants was included. All curated gene sets from Molecular Signatures Database (MSigDB) (Liberzon et al., Reference Liberzon, Birger, Thorvaldsdottir, Ghandi, Mesirov and Tamayo2015) were checked (refer to online Supplementary materials for details).
Parallel Independent Component Analysis
pICA extracts components (a weighted combination of variables) for genetic and imaging modalities in parallel based on infomax ICA, and enhances the inter-modality association through a constraint on the correlation (Liu et al., Reference Liu, Pearlson, Windemuth, Ruano, Perrone-Bizzozero and Calhoun2009). The algorithm for pICA can be found in online Supplementary material. Practical implementation of pICA is explained in detail in a recent review paper (Khadka et al., Reference Khadka, Pearlson, Calhoun, Liu, Gelernter, Bessette and Stevens2016). To achieve the best performance of pICA, the sample size to number of variables ratio should be larger than 0.02 (Liu et al., Reference Liu, Demirci and Calhoun2008). In this study, we only included variants with MAF ⩾5% in ExAC database and p value ⩽0.1 in logistic regression of MDD adjusting for age, sex, and top two PCs. The component number was estimated using minimum description length criteria for imaging modality and component consistency (Chen et al., Reference Chen, Calhoun and Liu2012a) for genetic modality with details shown in online Supplementary material. PICA was then performed through the Fusion ICA Toolbox2 within Matlab8.3 (FIT, http://mialab.mrn.org/software/fit/index.html). The weights of components were transformed into Z-score. To avoid false-positive associations, we constrained the significance to be only one pair of components. Then 10-fold cross-validation and 1000-run permutation test were performed to assess the validity of identified correlation (online Supplementary material).
Subsequent analyses after pICA
We extracted loading coefficients (LCs), which reflect the influence of components on subjects, in correlated component pairs. Two-sample t test was used to examine group differences of LCs. Partial correlation of LCs in two modalities controlling for age, sex, first two PCs, and MDD diagnosis was computed. Bonferroni correction was used for multiple comparisons. We also assessed respective linear correlations between pICA-derived components and clinical, cognitive measures. After the significant finding of association between LCs of GMV4 and IQ, we performed mediation analysis with MDD diagnosis as the outcome variable, IQ-related variables as the mediators, LCs of GMV4 as the exposure of interest, and years of education as a covariate, to investigate the direct and indirect causal relationships between exposure and outcome. A method specifically developed for common binary outcomes, which exploits the relationship between logit and probit models, was used (Gaynor et al., Reference Gaynor, Schwartz and Lin2018). Odds ratio (OR) was then calculated and bootstrapping with 500 random bootstrap samples was performed to obtain confidence intervals (CI). More details are available in online Supplementary material. The statistical analyses were performed in R3.3.2.
In the significant component pair, brain regions with thresholds |Z-score| >2.0 and volume >1.0 cm3 are selected to represent the imaging component. For the genetic component, we extracted genes with top-weighted variants and computed their overlaps with gene sets in curated gene sets of MSigDB. False discovery rate (FDR) was used to correct for multiple comparison.
Results
Demographic and clinical information
Summary for WES quality
Following extensive QC, the final data set was composed of 77 cases and 245 controls and had an average depth of 52×. Five controls (three with sample contamination; two related with another subject) but no case were removed. Overlapping regions of two capture kits consisted of 47.7 Mbp, with 93.7% achieving at least 10× depth. In total, 243 038 single nucleotide variants were identified. Consistency between the genotypes of two biological replicates captured by two kits in separate batches was 99.7%. No population stratification or batch effect has been observed using PC analysis as shown in online Supplementary Fig. S1.
CMC burden test of rare variants in MDD
In gene-based CMC burden tests, two genes, i.e. CSMD1 (CUB And Sushi Multiple Domains 1; p = 5.32×10−6) and CNTNAP5 (Contactin Associated Protein-Like 5; p = 1.82×10−6) were found associated with MDD. TMEM190 (Transmembrane Protein 190; p = 5.92×10−6) got a p value close to the significance threshold after Bonferroni correction (Fig. 1 and Table 2, refer to online Supplementary Table S1 and Fig. S2 for details of the variants). All variants were validated by Sanger sequencing.
Quantile–quantile plot (Q–Q plot) of gene-based association. In category A, gene CSMD1 achieved statistical significance after correcting for 4920 genes; (B) in category B, gene TMEM190 got a p value close to the significance threshold after correcting for 8647 genes; (C) in category C, no significance have been found after correcting for 3152 genes; (D) in category D, gene CNTNAP5 achieved significance after correcting for 6087 genes (category A: non-synonymous variants with MAF ⩽1% and predicted as potentially damaging; category B: non-synonymous variants with MAF ⩽1%; category C: non-synonymous variants with MAF ⩽0.1% and predicted as potentially damaging; category D: non-synonymous variants with MAF ⩽0.1%.).
Gene-based CMC burden tests results
The p values in bold are significant after correcting for the number of genes in the according category. Non-syno stands for non-synonymous. The four categories resulting from the filtering of MAF (minor allele frequency) and functional prediction included: category A: non-synonymous variants with MAF ⩽1% and predicted as potentially damaging; category B: non-synonymous variants with MAF ⩽1%; category C: non-synonymous variants with MAF ⩽0.1% and predicted as potentially damaging; category D: non-synonymous variants with MAF ⩽0.1%
In gene set-based CMC burden test, among 3367 curated gene sets from MSigDB (Liberzon et al., Reference Liberzon, Birger, Thorvaldsdottir, Ghandi, Mesirov and Tamayo2015), one pathway, Neuroactive Ligand Receptor Interactive (NLRI, M13380; p = 1.29×10−5) showed significantly higher burden in cases with MDD (online Supplementary Fig. S3). NLRI comprises 272 genes, where seven stop-gain and one start-loss variants in different genes were found in patients. All variants were successfully validated by Sanger sequencing as shown in online Supplementary Table S2.
Imaging–genetic analyses with pICA
Imaging–genetic correlation
Ten genetic components (VAR1-10) derived from 2460 variants and eight imaging components (GMV1-8) (online Supplementary Fig. S4) were extracted in pICA. One pair of components (VAR9-GMV4) displayed significant correlation (r = 0.38, p = 9.92×10−6, Bonferroni corrected) and passed 10-fold cross-validation. We observed high levels of similarity for VAR9 and GVM4 in each run, and also stable correlations (online Supplementary Table S3). In a 1000-run permutation test, the absolute values of the top VAR-GMV correlations ranged from 0.14 to 0.45 with a median of 0.29, yielding a significant p value of 0.03 (online Supplementary Fig. S5). After adjusting for covariates, the correlation remained significant (r = 0.34, p = 4.06×10−4, Bonferroni corrected). Patients with MDD tended to show stronger VAR-GMV correlation than controls as depicted in Fig. 2, although statistical significance was not achieved in Fisher's r-to-z transformation. We then further examined the partial correlation within each group and observed correlation coefficients of 0.37 (p = 1.61×10−3, uncorrected) and 0.27 (p = 17.04×10−3, uncorrected) for cases and controls, respectively.
Scatter plots of loading coefficients (LCs) of GMV4 and VAR9 components (MDD, major depressive disorder; HC, healthy control).
The spatial map of the GMV4 component with absolute value of Z-score >2.0 and volume >1.0 cm3 is illustrated in Fig. 3. The identified brain regions included bilateral precuneus, posterior cingulate, middle/inferior/superior temporal gyrus, inferior parietal lobule, postcentral gyrus, middle frontal gyrus, supramarginal gyrus, and insula (online Supplementary Table S4). All the representative regions had positive Z-score and patients with MDD had lower LCs of GMV4 compared with controls (p = 0.01, Fig. 2), suggesting that cases had reduced GMV in these regions. Based on the absolute value of the component weight, we selected 79 genes with variants ranked in the top 5% in VAR9 for subsequent pathway enrichment analysis (online Supplementary Table S5 and Fig. S6). The results highlighted two significant pathways: Singling by G-protein coupled receptors (GPCRs, M746, FDR q = 3.23×10−4) and Alzheimer Disease Up (M12921, FDR q = 6.12×10−4) (online Supplementary Table S6).
Brain regions in GMV4. Only regions with absolute value of Z-score >2.0 and volume >1.0 cm3 are displayed. Z-score is represented by color, with yellow meaning higher positive value. No region with >1.0 cm3 volume had a negative Z-score with absolute value larger than 2.0.
Post hoc analyses
Correlation analysis found significance between LCs of GMV4 and three cognitive measures, i.e. full IQ (r = 0.27, p = 0.002, uncorrected), verbal IQ (r = 0.26, p = 0.004, uncorrected), and performance IQ (r = 0.22, p = 0.013, uncorrected). No cognitive measure was significantly correlated with LCs of VAR9. No significant association between LCs and HAMD-17 score was observed.
In the mediation analysis, indirect effects were significant for full IQ (OR 0.86, 95% CI 0.73–0.97, p = 0.02) and verbal IQ (OR 0.85, 95% CI 0.71–0.97, p = 0.01), but not performance IQ (online Supplementary Fig. S7a). Standardized regression coefficients for the relationship between GMV and MDD as mediated by IQ are illustrated in online Supplementary Fig. S7b. Direct effects were not significant which suggested complete mediation by full IQ and verbal IQ. However, the OR estimates of direct effects deviated more from the null value (i.e. 1) than the estimates for indirect effects, and the lack of statistical significance may be due to inadequate sample size resulting in wide CI. Thus we concluded that full IQ and verbal IQ at least partially mediated the effect of GMV on MDD.
Discussion
In this study, we integrated rare non-synonymous/damaging variants from WES into sets to explore the association of genes/pathways and MDD. We found that the burden of rare variants in MDD is largely concentrated in two genes (CSMD1, CNTNAP5) and in one pathway (NLRI). Furthermore, we showed for the first time that reduced GMV in a number of brain regions in patients with MDD (precuneus, posterior cingulate, temporal gyrus, inferior parietal lobule, postcentral gyrus, middle frontal gyrus, supramarginal gyrus, and insula) were significantly associated with a genetic component. Genes underlying these GMV reductions are enriched in two pathways (Singling by GPCRs and Alzheimer Disease Up) hitting genes involved in neuroinflammation and synaptic plasticity. Furthermore, the identified reduced GMV was significantly associated with poorer cognitive functioning.
Genes and gene set with higher burden of rare variants in patients with MDD
CSMD1 protein was found to be predominantly expressed in brain tissue according to Genotype-Tissue Expression (GTEx, online Supplementary Fig. S8). It can inhibit the deposition of complement component 3 (C3) and act as an important regulator of complement activation in the developing central nervous system (CNS) (Kraus et al., Reference Kraus, Elliott, Chute, Horan, Pfenninger, Sanford, Foster, Scully, Welcher and Holers2006). C3 has an essential role in both the peripheral and CNS immune system (Schafer et al., Reference Schafer, Lehrman, Kautzman, Koyama, Mardinly, Yamasaki, Ransohoff, Greenberg, Barres and Stevens2012). In addition, the classical complement cascade is involved in the plasticity of the CNS (Stevens et al., Reference Stevens, Allen, Vazquez, Howell, Christopherson, Nouri, Micheva, Mehalow, Huberman, Stafford, Sher, Litke, Lambris, Smith, John and Barres2007; Boulanger, Reference Boulanger2009), particularly, disrupting complement receptor 3 (CR3)/C3 signaling in microglia resulting in sustained deficits in synaptic connectivity (Schafer et al., Reference Schafer, Lehrman, Kautzman, Koyama, Mardinly, Yamasaki, Ransohoff, Greenberg, Barres and Stevens2012). Among current hypotheses to explain the etiology of MDD, neuroinflammation (Kim and Won, Reference Kim and Won2017) and synaptic plasticity (Marsden, Reference Marsden2013; Wainwright and Galea, Reference Wainwright and Galea2013) are now widely accepted due to substantial evidence. Patients with MDD show increased levels of proinflammatory cytokines and depressive symptoms may be induced in healthy participants by cytokines (Berk et al., Reference Berk, Williams, Jacka, O'NEIL, Pasco, Moylan, Allen, Stuart, Hayley, Byrne and Maes2013). Furthermore, anti-inflammatory drugs exhibit antidepressant effects in patients with MDD (Kohler et al., Reference Kohler, Benros, Nordentoft, Farkouh, Iyengar, Mors and Krogh2014). In addition, reduced GMV, lowered levels of neuropil, lowered expression of synaptic proteins and growth factors suggest abnormal synaptic plasticity of MDD (Marsden, Reference Marsden2013; Wainwright and Galea, Reference Wainwright and Galea2013). The most recent PGC study reported that major histocompatibility complex (MHC) region, which is involved in schizophrenia and other mental disorders (Schizophrenia Working Group of the Psychiatric Genomics Consortium, 2014), also plays a role in MDD (Major Depressive Disorder Working Group of the Psychiatric Genomics Consortium, 2018). In schizophrenia, reduced numbers of synapses in the brains of schizophrenic patients, which may result from excessive complement activity caused by structurally diverse alleles of gene C4, was reported (Sekar et al., Reference Sekar, Bialas, de Rivera, Davis, Hammond, Kamitaki, Tooley, Presumey, Baum, Van Doren, Genovese, Rose, Handsaker, Daly, Carroll, Stevens and Mccarroll2016). In addition, an enrichment of significant loci in genes involved in cytokine and immune response was also found in MDD (Major Depressive Disorder Working Group of the Psychiatric Genomics Consortium, 2018). Our findings of CSMD1 and previous findings of MHC both indicate the importance of neuroinflammation and abnormal synaptic plasticity in MDD. In addition to marginal association of CSMD1 and MDD in previous GWASs [rs 2 616 976, p = 3.60×10−6 (Shi et al., Reference Shi, Potash, Knowles, Weissman, Coryell, Scheftner, Lawson, Depaulo, Gejman, Sanders, Johnson, Adams, Chaudhury, Jancic, Evgrafov, Zvinyatskovskiy, Ertman, Gladis, Neimanas, Goodell, Hale, Ney, Verma, Mirel, Holmans and Levinson2011); rs 117 511 508, p = 1.44×10−6 (Hyde et al., Reference Hyde, Nagle, Tian, Chen, Paciga, Wendland, Tung, Hinds, Perlis and Winslow2016); rs 116 959 294, p = 6.88×10−6 (Hyde et al., Reference Hyde, Nagle, Tian, Chen, Paciga, Wendland, Tung, Hinds, Perlis and Winslow2016)], linkage analysis also found that 8p23, where CSMD1 is located, had linkage with MDD (Bulayeva et al., Reference Bulayeva, Lencz, Glatt, Takumi, Gurgenova, Kawakami and Bulayev2012). Knockout mouse models have shown anxiety-like and depression-like behaviors in CSMD1 mutants (Steen et al., Reference Steen, Nepal, Ersland, Holdhus, Naevdal, Ratvik, Skrede and Havik2013). Taken together, the above evidence strongly suggests that CSMD1 is not only involved in MDD, but is also a core gene for MDD according to the newly developed omnigenic model (Boyle et al., Reference Boyle, Li and Pritchard2017).
CNTNAP5 encodes one of the five protein members of neurexin family (CNTNAP1-5), which act as key regulators of synapse formations and neurotransmission (Anderson et al., Reference Anderson, Galfin, Xu, Aoto, Malenka and Sudhof2012; Karayannis et al., Reference Karayannis, Au, Patel, Kruglikov, Markx, Delorme, Heron, Salomon, Glessner, Restituito, Gordon, Rodriguez-Murillo, Roy, Gogos, Rudy, Rice, Karayiorgou, Hakonarson, Keren, Huguet, Bourgeron, Hoeffer, Tsien, Peles and Fishell2014). CNTNAP5 is mainly expressed in brain tissue based on GTEx (online Supplementary Fig. S8). In addition to the evidence showing the participation of CNTNAP5 in schizophrenia (McClay et al., Reference Mcclay, Adkins, Aberg, Stroup, Perkins, Vladimirov, Lieberman, Sullivan and van den Oord2011; Levinson et al., Reference Levinson, Shi, Wang, Oh, Riley, Pulver, Wildenauer, Laurent, Mowry, Gejman, Owen, Kendler, Nestadt, Schwab, Mallet, Nertney, Sanders, Williams, Wormley, Lasseter, Albus, Godard-Bauche, Alexander, Duan, O'DONOVAN, Walsh, O'NEILL, Papadimitriou, Dikeos, Maier, Lerer, Campion, Cohen, Jay, Fanous, Eichhammer, Silverman, Norton, Zhang, Hakonarson, Gao, Citri, Hansen, Ripke, Dudbridge and Holmans2012) and bipolar disorder (Djurovic et al., Reference Djurovic, Gustafsson, Mattingsdal, Athanasiu, Bjella, Tesli, Agartz, Lorentzen, Melle, Morken and Andreassen2010), another member of this family, CNTNAP2, is found to be associated with MDD (Ji et al., Reference Ji, Li, Pan, Tao, Ju, Wen, Fu, An, Zhao, Wang, He, Feng, Yi and Shi2013). Moreover, several members of this family are involved in synaptic plasticity. For example, knockdown of CNTNAP2 leads to a global decline in excitatory and inhibitory synapse numbers and a decrease in synaptic transmission (Anderson et al., Reference Anderson, Galfin, Xu, Aoto, Malenka and Sudhof2012). CNTNAP4 is localized presynaptically and its loss leads to a reduction in the output of cortical PV-positive GABAergic basket cells (Karayannis et al., Reference Karayannis, Au, Patel, Kruglikov, Markx, Delorme, Heron, Salomon, Glessner, Restituito, Gordon, Rodriguez-Murillo, Roy, Gogos, Rudy, Rice, Karayiorgou, Hakonarson, Keren, Huguet, Bourgeron, Hoeffer, Tsien, Peles and Fishell2014). Given shared features among CNTNAP1-5, it is possible that the CNTNAP5 is also heavily implicated in synaptic plasticity.
TMEM190 encodes a single-pass transmembrane protein with trefoil domain. Although no involvement in depression has been reported yet, TMEM190 has been found to express in MDD-related brain regions such as hippocampus and hypothalamus (online Supplementary Fig. S8). Copy number gain in Chr19q13.33–13.43, where TMEM190 is located in, is associated with intellectual and developmental disabilities (Kaminsky et al., Reference Kaminsky, Kaul, Paschall, Church, Bunke, Kunig, Moreno-De-luca, Moreno-De-luca, Mulle, Warren, Richard, Compton, Fuller, Gliem, Huang, Collinson, Beal, Ackley, Pickering, Golden, Aston, Whitby, Shetty, Rossi, Rudd, South, Brothman, Sanger, Iyer, Crolla, Thorland, Aradhya, Ledbetter and Martin2011) (refer to online Supplementary materials for more details).The pathway with a higher burden of disruptive variants, NLRI is involved in neurodevelopment and has potential participation in many neuropsychiatric diseases, including schizophrenia (Putnam et al., Reference Putnam, Sun and Zhao2011; Adkins et al., Reference Adkins, Khachane, Mcclay, Aberg, Bukszar, Sullivan and van den Oord2012), ADHD (Alemany et al., Reference Alemany, Ribases, Vilor-Tejedor, Bustamante, Sanchez-Mora, Bosch, Richarte, Cormand, Casas, Ramos-Quiroga and Sunyer2015; Khadka et al., Reference Khadka, Pearlson, Calhoun, Liu, Gelernter, Bessette and Stevens2016) and Parkinson's disease (Kong et al., Reference Kong, Liang, Liu, Zhang, Wan, Gan and Yuan2015). Among genes with disruptive variants in NLRI, we found two (F2RL1, PTGER4) were explicitly implicated in inflammatory processes and two (TACR1, NPBWR2) related to inflammatory processes. The activation of PAR2 encoded by F2RL1 in glia suppresses the astrocytic anti-inflammatory processes (Afkhami-Goli et al., Reference Afkhami-Goli, Noorbakhsh, Keller, Vergnolle, Westaway, Jhamandas, Andrade-Gordon, Hollenberg, Arab, Dyck and Power2007). PTGER4 protein can suppress the proinflammatory gene in the model of innate immunity (Shi et al., Reference Shi, Johansson, Woodling, Wang, Montine and Andreasson2010). Moreover, PEG2, the ligand of PTGER4 has increased circulating levels in patients with MDD (Zorrilla et al., Reference Zorrilla, Luborsky, Mckay, Rosenthal, Houldin, Tax, Mccorkle, Seligman and Schmidt2001) and antidepressant medication can reduce the PEG level (Bah et al., Reference Bah, Benderdour, Kaloustian, Karam, Rousseau and Godbout2011). Gene TACR1 or NK-1R encodes the receptor of substance P, which is a neuropeptide with multiple functions, including initiating expression of cytokines in CNS (Rameshwar et al., Reference Rameshwar, Gascon and Ganea1992; Palma and Manzini, Reference Palma and Manzini1998) and participating in pain perception (De Felipe et al., Reference De Felipe, Herrero, O'BRIEN, Palmer, Doyle, Smith, Laird, Belmonte, Cervero and Hunt1998). Multiple lines of evidence indicate that substance P and its receptor are involved in MDD, including elevation of substance P in cerebrospinal fluid (Geracioti et al., Reference Geracioti, Carpenter, Owens, Baker, Ekhator, Horn, Strawn, Sanacora, Kinkead, Price and Nemeroff2006), decreased binding to NK-1 receptors across all cortical layers (Stockmeier et al., Reference Stockmeier, Shi, Konick, Overholser, Jurjus, Meltzer, Friedman, Blier and Rajkowska2002) in patients with MDD, and antidepressive effects of antagonist of substance P (Yang et al., Reference Yang, Yu, Jin and Zhao2014). In addition, increased levels of substance P in patients with MDD may relate to the pathophysiological process of unexplained painful physical symptoms reported by patients with MDD (Jaracz et al., Reference Jaracz, Gattner, Jaracz and Gorna2016). Gene NPBWR2 encodes the receptors of neuropeptide B and W, which have an analgesic effect in inflammatory pain (Yamamoto et al., Reference Yamamoto, Saito, Shono and Tanabe2005). Neuropeptide B gene knockout mice has also exhibited hyperalgesia to inflammatory pain (Kelly et al., Reference Kelly, Beuckmann, Williams, Sinton, Motoike, Richardson, Hammer, Garry and Yanagisawa2005).
Our results did not replicate the previous two WES studies in MDD patients (Tammiste et al., Reference Tammiste, Jiang, Fischer, Magi, Krjutskov, Pettai, Esko, Li, Tansey, Carroll, Uher, Mcguffin, Vosa, Tsernikova, Saria, Ng, Eller, Vasar, Nutt, Maron, Wang and Metspalu2013; Tombacz et al., Reference Tombacz, Maroti, Kalmar, Csabai, Balazs, Takahashi, Palkovits, Snyder and Boldogkoi2017) or in two population cohorts (Amin et al., Reference Amin, Belonogova, Jovanova, Brouwer, van Rooij, van den Hout, Svishcheva, Kraaij, Zorkoltseva, Kirichenko, Hofman, Uitterlinden, van Ijcken, Tiemeier, Axenovich and van Duijn2017a, Reference Amin, Jovanova, Adams, Dehghan, Kavousi, Vernooij, Peeters, de Vrij, van der Lee, van Rooij, van Leeuwen, Chaker, Demirkan, Hofman, Brouwer, Kraaij, Willems, Hankemeier, van Ijcken, Uitterlinden, Niessen, Franco, Kushner, Ikram, Tiemeier and van Duijn2017b). This may be due to different characteristics of samples. The two previous WES studies had small sample sizes. One focused on pharmacological responses to an antidepressant, the other on suicide victims with MDD. The former used Mini International Neuropsychiatric Interview (Sheehan et al., Reference Sheehan, Lecrubier, Sheehan, Amorim, Janavs, Weiller, Hergueta, Baker and Dunbar1998) and a ‘moderate or above’ rating on the Montgomery–Asberg Depression Rating Scale (Montgomery and Asberg, Reference Montgomery and Asberg1979), which is comparable in quality of case definition to our study but included a very small sample. The second used case files to diagnose MDD in deceased subjects. In the population cohort studies, a self-report measure of depressive symptoms in the prior week was used as a continuous score to screen a population for ‘possible depression’ usually defined at a cutoff of 16 out of a maximum score of 60.
GMV reduction in MDD and its underlying genetic components
Several brain regions in GMV4 are related to the default mode network (DMN) which is a brain network of interacting regions activated during resting consciousness. Precuneus and posterior cingulate cortex (PCC) are core nodes of DMN (Buckner et al., Reference Buckner, Andrews-Hanna and Schacter2008). Reduced precuneus volume (Shen et al., Reference Shen, Cheng, Yang, Dai, Ye, Liu, Lu, Li, Liu, Lu, Sun and Xu2016), thinner cortical gray matter in PCC (Schmaal et al., Reference Schmaal, Hibar, Samann, Hall, Baune, Jahanshad, Cheung, van Erp, Bos, Ikram, Vernooij, Niessen, Tiemeier, Hofman, Wittfeld, Grabe, Janowitz, Bulow, Selonke, Volzke, Grotegerd, Dannlowski, Arolt, Opel, Heindel, Kugel, Hoehn, Czisch, Couvy-Duchesne, Renteria, Strike, Wright, Mills, de Zubicaray, Mcmahon, Medland, Martin, Gillespie, Goya-Maldonado, Gruber, Kramer, Hatton, Lagopoulos, Hickie, Frodl, Carballedo, Frey, van Velzen, Penninx, van Tol, van der Wee, Davey, Harrison, Mwangi, Cao, Soares, Veer, Walter, Schoepf, Zurowski, Konrad, Schramm, Normann, Schnell, Sacchet, Gotlib, Macqueen, Godlewska, Nickson, Mcintosh, Papmeyer, Whalley, Hall, Sussmann, Li, Walter, Aftanas, Brack, Bokhan, Thompson and Veltman2016), and different correlation between activity in precuneus and non-DMN regions (Peng et al., Reference Peng, Liddle, Iwabuchi, Zhang, Wu, Liu, Jiang, Xu, Liddle, Palaniyappan and Fang2015) have been reported in MDD patients. In a recent study of 9707 participants by UK Biobank, several SNPs were significantly associated with precuneus volume, which is strong evidence that the volume of the precuneus is affected by genes (Elliott et al., Reference Elliott, Sharp, Alfaro-Almagro, Shi, Miller, Douaud, Marchini and Smith2018). Another region in GMV4, insula that encompasses emotional processing, has clusters that are negatively correlated with the PCC (Uddin et al., Reference Uddin, Kelly, Biswal, Castellanos and Milham2009). Decreased functional connectivity is also found between insula and the left superior temporal gyrus in drug-naive MDD patients (Guo et al., Reference Guo, Liu, Xiao, Zhang, Liu, Yu, Zhang and Zhao2015). The temporal lobe is the largest element of GMV4. Involved in memory and perception (Baxter, Reference Baxter2009; Lee and Rudebeck, Reference Lee and Rudebeck2010), temporal lobe is involved in the neurobiology of mood disorders (Cusi et al., Reference Cusi, Nazarov, Holshausen, Macqueen and Mckinnon2012). In patients with MDD, reduced temporal gyrus volume (Du MY et al., Reference Du, Wu, Yue, Li, Liao, Kuang, Huang, Chan, Mechelli and Gong2012; Shen et al., Reference Shen, Cheng, Yang, Dai, Ye, Liu, Lu, Li, Liu, Lu, Sun and Xu2016; Zhang et al., Reference Zhang, Li, Wu, Chen, Hu, Chen, Zhu, Jia and Gong2016a) as well as oligodendroglial abnormalities in the temporal cortex (Aston et al., Reference Aston, Jiang and Sokolov2005) has been reported. In addition, a meta-analysis identified reduced right middle frontal gyrus volume in MDD compared with controls (Lai, Reference Lai2013).
In the present study, we found significant correlations between LCs of GMV4 and IQ, supporting the biological associations between gene, brain volume, and intelligence. We did not observe significant association of the LCs with HAMD-17 scores which are only available for patients. However, mediation analysis showed that full IQ, and more specifically verbal IQ, partially mediated the effect of GMV4 on MDD. Brain regions in GMV4 have been reported to be involved in cognition, especially verbal intelligence-related cognition. Verbal creativity has been found to be negatively correlated with functional connectivity within the medial prefrontal cortex of DMN (Zhu et al., Reference Zhu, Chen, Xia, Beaty, Yang, Tian, Sun, Cao, Zhang, Chen and Qiu2017). Internet use has been associated with decrease of verbal intelligence and smaller longitudinal increase in GMV of widespread brain areas, including insula and temporal area. Lifelong association between ventral temporal volume and verbal knowledge has been reported (Hoffman et al., Reference Hoffman, Cox, Dykiert, Munoz, Valdes, Bastin, Wardlaw and Deary2017). Temporal lobe resection for patients with temporal lobe epilepsy increases the risk of verbal memory decline (Wagner et al., Reference Wagner, Uherek, Horstmann, Kadish, Wisniewski, Mayer, Buschmann, Metternich, Zentner and Schulze-Bonhage2013).
Genes with the greatest weights in VAR9 are enriched in Singling by GPCRs pathway, which has been linked to GMV reduction of schizophrenia in a previous pICA study (Jagannathan et al., Reference Jagannathan, Calhoun, Gelernter, Stevens, Liu, Bolognani, Windemuth, Ruano, Assaf and Pearlson2010). Approximately 90% GPCRs are located in the human brain and regulate neurotransmitters (Vassilatis et al., Reference Vassilatis, Hohmann, Zeng, Li, Ranchalis, Mortrud, Brown, Rodriguez, Weller, Wright, Bergmann and Gaitanaris2003). Among the overlapping genes (online Supplementary Table S6), several (CCR2, DRD5, GRM1, and ITPR2) are worth noting. In a mouse model of depression, microglia and brain macrophages with increased CCR2 protein exhibited anti-inflammatory function (Chabry et al., Reference Chabry, Nicolas, Cazareth, Murris, Guyon, Glaichenhaus, Heurteaux and Petit-Paitel2015). In mice treated with lipopolysaccharide, the percentage of CCR2-expressing microglia cells increased significantly, which caused a reduction of 5-HT release in the amygdala and hippocampus (Cazareth et al., Reference Cazareth, Guyon, Heurteaux, Chabry and Petit-Paitel2014). DRD5 is important given the monoamine hypothesis of MDD, which hypotheses depression is caused by underactivity in brain monoamines, such as dopamine. Nominal association of SNP in DRD5 and corresponding depressive scores were observed (Nyman et al., Reference Nyman, Sulkava, Soronen, Miettunen, Loukola, Leppa, Joukamaa, Maki, Jarvelin, Freimer, Peltonen, Veijola and Paunio2011). The glutamate system is disrupted in MDD and higher expression levels of GRM1 protein were detected in the female patients with MDD in a large cohort of postmortem subjects (Gray et al., Reference Gray, Hyde, Deep-Soboslay, Kleinman and Sodhi2015). ITPR2 protein participates in neuronal calcium shuttling and ITPR2 knockout mice were observed to exhibit depressive-like behaviors (Cao et al., Reference Cao, Li, Wang, Wu, Hu, Zhang, Fang, Zhang, Li, Xiong, Yan, Gao, Liu, Li, Sun, Zeng, Zhu and Gao2013).
Another pathway enriched in VAR9 is Alzheimer Disease Up, composed of genes upregulated in hippocampus samples from patients with Alzheimer's disease (AD) (Blalock et al., Reference Blalock, Geddes, Chen, Porter, Markesbery and Landfield2004). According to a recent review, MDD is a risk factor for AD (Herbert and Lucassen, Reference Herbert and Lucassen2016). Prominent participation of neuroinflammation in both diseases is generally acknowledged. Moreover, it is postulated that progress from depression to dementia could result from activation of macrophages in the blood and microglia in the brain (Leonard, Reference Leonard2007). Previous studies have found shared genetic risks between MDD and AD, e.g. variants in gene CR1 are associated with both MDD and AD (Hamilton et al., Reference Hamilton, Evans, Macintyre, Deary, Dominiczak, Smith, Morris, Porteous and Thomson2012). Among 15 overlapping genes, three (C7, NFASC, and GPR1) are involved in neuroinflammation and one (SEMA4G) participates in synaptic plasticity. C7 plays a vital role in complement activation. NFASC encodes neurofascin, a member of the immunoglobulin cell adhesion molecule superfamily. GPR1 encodes chemerin receptors and chemerin acts as a potent chemoattractant for macrophages (Hart and Greaves, Reference Hart and Greaves2010). Additionally, protein semaphorin 4g partakes in neurotransmission and synaptic plasticity by interacting with post-synaptic density proteins (Burkhardt et al., Reference Burkhardt, Muller, Badde, Garner, Gundelfinger and Puschel2005). Chronic treatment with vortioxetine can regulate the expression of SEMA4G in the brain. More importantly, semaphorin 4g is an axonal guidance molecule and many identified genes associated with imaging phenotypes have a role in exon guidance (Elliott et al., Reference Elliott, Sharp, Alfaro-Almagro, Shi, Miller, Douaud, Marchini and Smith2018).
Both synaptic plasticity and neuroinflammation could lead to changes in GMV. For example, positive experience and training increase regional volume in humans through influencing synaptic plasticity and, vice versa, negative experiences may decrease volume (Anderson, Reference Anderson2011). In a postmortem study, schizophrenia patients of ‘high inflammation’ status had reduced cortical GMV compared with cases of ‘low inflammation’ status (Zhang et al., Reference Zhang, Catts, Sheedy, Mccrossin, Kril and Shannon2016b).
Limitations
The current study has several limitations. Firstly, variants were selected from the association results of WES data with MDD for the analysis of pICA, resulting in non-exomic or non-associated variants not being included in our analysis. However, we effectively focused our core analyses on potential disease-related genes and pathways using this approach (Pearlson et al., Reference Pearlson, Liu and Calhoun2015). It is also argued that the filtering strategy of an arbitrary p value threshold in association test could bias the significance level of group difference in the genetic component, rather than the significance of identified genetics–imaging correlations (Chen et al., Reference Chen, Calhoun, Pearlson, Ehrlich, Turner, Ho, Wassink, Michael and Liu2012b; Gupta et al., Reference Gupta, Chen, Liu, Damaraju, Wright, Perrone-Bizzozero, Pearlson, Luo, Michael, Turner and Calhoun2015). In addition, as a preliminary WES study of MDD, the sample size is relatively small compared with previous studies of schizophrenia and autism. Independent replication is important for verifying the reliability of the results in the present study, while the expensive price of WES did not allow us have another available sample to do further validating test. We call for further replication studies with independent and larger sample sets.
Conclusion
In conclusion, we have identified three genes and a pathway embracing higher burden of rare variants in patients with MDD and revealed a genetic component underlying the GMV reductions in MDD. These two results are in line with the important neuroinflammation and synaptic plasticity hypotheses of the disorder. More importantly, cognition mediated the effect of GMV changes on MDD. This indicates that the identified genetic component may act on MDD via GMV alteration and through cognitive impairment. GMV, as an intermediate phenotype, can play an important role in decoding the genetic mechanism of MDD.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/S0033291719000072.
Acknowledgements
Samples in this study were recruited by psychiatrists in Prof. Tao Li's laboratory. The analyses were performed with the help of senior students in Prof. Pak C. Sham's laboratory. Many thanks to them.
Financial support
This work was supported by the National Nature Science Foundation of China Key Project (T.L., grant numbers 91332205, 81130024, 81630030); National Key Technology R & D Program of the Ministry of Science and Technology of China (T.L., grant number 2016YFC0904300); National Natural Science Foundation of China/Research Grants Council of Hong Kong Joint Research Scheme (T.L., grant number 8141101084); Natural Science Foundation of China (W.D. grant number 8157051859); Sichuan Science & Technology Department (Q.W., grant number 2015JY0173); and 1.3.5 Project for disciplines of excellence, West China Hospital of Sichuan University (T.L., grant numbers ZY2016103, ZY2016203). This research received no specific grant from any funding agency, commercial or not-for-profit sectors.
Conflict of interest
None.
Ethical standards
The authors assert that all procedures contributing to this work comply with the ethical standards of the relevant national and institutional committees on human experimentation and with the Helsinki Declaration of 1975, as revised in 2008.
