Figure 1. Gene expression changes in the cerebral cortex of persons with schizophrenia (SZ). A, Venn diagram depicting the overlap of genes differentially expressed between persons with SZ and controls in Brodmann areas (BAs) 21, 32, 38, and 46. Differential expression was assessed with unpaired 2-class test analysis (using the significance analysis of microarrays package) at a false discovery rate of less than 0.05 and a fold change of greater than 1.3. B, Differentially expressed genes between persons with SZ and controls in each BA and paired overlap analysis between the 4 BAs assessed with the hypergeometric distribution. C, Heat map of top 200 genes differentially expressed between SZ and control samples of the cerebral cortex in BAs 21, 32, and 38. Scaled expression values are color-coded. The dendrogram depicts hierarchical clustering based on the Euclidean distance and the complete linkage method. The bars show variables for each sample: diagnosis (Dx) (green, persons with SZ; red, controls), BA (red, BA 21; green, BA 32; yellow, BA 38), sex (red, male; green, female), age of death (AOD), and postmortem interval (PMI). The corresponding scale for quantitative variables is shown on the left.
Figure 2. Gene coexpression modules associated with schizophrenia (ie, visualization of the OLIG1 [A] and GABA [B] modules). The heat map of intramodular expression genes is shown in the top of the left side for each module, and the corresponding module eigengene values, across samples according to Brodmann area (BA), are shown on the middle left side (black, BA 21; red, BA 32; green, BA 38; blue, BA 46). Visualization of the OLIG1 and GABA modules is provided on the right side. The top 1000 connections are shown for each module on the top right side, according to a degree-sorted circle layout. The intramodular hubs are shown in larger size at the bottom, where the relevant gene ontology categories enriched in the OLIG1 and GABA modules are listed on the bottom left side.
Figure 3. Genome-wide association enrichment analysis of modules associated with schizophrenia (SZ). Several genome-wide association studies examining a variety of illnesses using several cohorts (eg, SZ from the Psychiatric Genomics Consortium [PGC_SZ],18 bipolar disorder from the PGC [PGC_BD],19 attention-deficit/hyperactivity disorder [ADHD],20 and major depression disorder [MDD]21) and 3 nonpsychiatric disease studies as negative controls in human immunodeficiency virus (HIV),22 psoriasis,23 and Parkinson disease (PD)24 were included. Single-nucleotide polymorphisms located within the transcript boundaries and an additional 20 kilobase on the 5′ end and the 3′ end of the gene were included in the gene set analysis method.15 The analyses presented here are based on k = 5. For multiple testing corrections, a Benjamini-Hochberg correction was applied, and α was set at 10−4.
Figure 4. Cortical interregional differences shown as diminished in schizophrenia (SZ). A, Venn diagram depicting the overlap of differentially expressed genes between persons with SZ and controls in all regions (611 genes) and genes differentially expressed among the 4 different Brodmann areas (BAs) in controls (681 genes) and persons with SZ (82 genes). The 42 genes (in red font) are the genes that overlap between differentially expressed genes between persons with SZ and controls in all regions and genes differentially expressed among the 4 different BAs in controls. B, The 42 genes are differentially expressed among BAs in control samples and show different expression levels in control samples compared with SZ samples. The horizontal bars show the −log P values for differential expression among all BAs in SZ and control groups. The majority (n = 34) of these genes (ie, those with an asterisk) overlap with the cUNK1 module, which is 1 of the 3 modules that show a different regional expression in controls. C, Visualization of the top 1000 connections for the cUNK1 module is shown on the left side, and a more detailed view of the hub cUNK1 genes is shown on the right. Three of the 42 genes (SCARNA9, SMC3, and TOP1) were intramodular hub genes in cUNK1.
Figure 5. A hypothetical model based on the current data according to the disconnectivity hypothesis in schizophrenia.
Roussos P, Katsel P, Davis KL, Siever LJ, Haroutunian V. A system-level transcriptomic analysis of schizophrenia using postmortem brain tissue samples. Arch Gen Psychiatry.. Published online August 6, 2012. doi:10.1001/archgenpsychiatry.2012.704..
eAppendix. Supplementary methods and results.
eFigure 1. Composite preservation statistics of SZ modules in control samples.
eFigure 2. Module eigengene-based regional differences among cases and controls.
eTable 1 Schizophrenia cases and controls.
eTable 2 TaqMan gene expression assays used in the study.
eTable 3 Characteristics of the sample groups used in this study..
eTable 4 MetaCore analysis of genes differentially expressed in cerebral cortex between SZ and controls.
eTable 5 Comparison between the coexpression network built on the entire data set, including both SZ and control samples, and the premade set of brain-derived enrichment lists in the WGCNA package.
eTable 6 Eigengene significance for all modules in the coexpression network build on the entire data set, including both SZ and control samples.
eTable 7 Comparison between the control coexpression network and the premade set of brain-derived enrichment lists in the WGCNA package.
eTable 8 Comparison between the SZ coexpression network and the premade set of brain-derived enrichment lists in the WGCNA package.
Roussos P, Katsel P, Davis KL, Siever LJ, Haroutunian V. A System-Level Transcriptomic Analysis of Schizophrenia Using Postmortem Brain Tissue Samples. Arch Gen Psychiatry. 2012;69(12):1205-1213. doi:10.1001/archgenpsychiatry.2012.704
Author Affiliations: Department of Psychiatry, The Mount Sinai School of Medicine, New York (Drs Roussos, Katsel, Davis, Siever, and Haroutunian), and James J. Peters Veterans Affairs Medical Center, Bronx (Drs Roussos, Siever, and Haroutunian), New York.
Context Schizophrenia is a common, highly heritable, neurodevelopmental mental illness, characterized by genetic heterogeneity.
Objective To identify abnormalities in the transcriptome organization among older persons with schizophrenia and controls.
Design Weighted gene coexpression network analysis based on microarray transcriptomic profiling.
Setting Research hospital.
Patients Postmortem brain tissue samples from 4 different cerebrocortical regions (the dorsolateral prefrontal cortex, the middle temporal area, the temporopolar area, and the anterior cingulate cortex) from 21 persons with schizophrenia and 19 controls.
Main Outcome Measures Results from gene expression microarray analysis, from analysis of coexpression networks, and from module eigengene, module preservation, and enrichment analysis of schizophrenia-related genetic variants.
Results The oligodendrocyte, microglial, mitochondrial, and neuronal (GABAergic and glutamatergic) modules were associated with disease status. Enrichment analysis of genome-wide association studies in schizophrenia and other illnesses demonstrated that the neuronal (GABAergic and glutamatergic) and oligodendrocyte modules are enriched for genetically associated variants, whereas the microglial and mitochondrial modules are not, providing independent support for more direct involvement of these gene expression networks in schizophrenia. Interregional coexpression network analysis showed that the gene expression patterns that typically differentiate the frontal, temporal, and cingulate cortices in controls diminish significantly in persons with schizophrenia.
Conclusions These results support the existence of convergent molecular abnormalities in schizophrenia, providing a molecular neuropathological basis for the disease.
Genetic and phenotypic heterogeneity, epistatic interactions, and environmental and epigenetic factors may limit insights into the pathophysiology of schizophrenia (SZ) that can be derived from studies of genetic variants and specific molecular pathways alone. Expression profiling studies in SZ conducted using microarray platforms suggest changes in functional gene groups such as oligodendrocyte- and myelin-related genes, in metabolism, in synaptic transmission, and, to a lesser extent, in genes associated with inflammatory processes.1 Agreement across different microarray studies about which individual transcripts are the most affected is limited, perhaps reflecting differences such as age, disease chronicity, disease severity, and diagnostic heterogeneity across the studied postmortem cohorts. Furthermore, the reported outcome of high-throughput microarray studies is usually a list of genes with the majority of them showing only modest changes between cases and controls and nominal statistical significance after multiple testing corrections. A more integrative tactic that organizes the transcriptome in persons with SZ and comparison controls into regional and functional networks through a systems biology approach may be more helpful.
Weighted gene coexpression network analysis (WGCNA) is a systems biology method that, based on gene coexpression patterns, identifies the higher-order relationships among genes, providing a holistic view of transcriptome organization.2,3 Highly correlated genes are assigned into clusters (modules), and WGCNA identifies the most highly connected genes (hubs) within each module. The expression levels of each module can be summarized by the first principal component (the module eigengene). Module eigengenes and the expression levels of intramodular hub genes can be used to assess whether modules are related to clinical phenotypes or other experimental variables, such as brain region and genetic variants. In addition to the heuristic value of unsupervised grouping of gene expression patterns, such a data reduction strategy has the potential to alleviate the multiple testing correction uncertainties that are encountered in standard gene-centric methodological approaches for microarray data analysis and, hence, decreases the potential for type I and II statistical errors.
It is now widely accepted that SZ is associated with a global gene expression disturbance across many cortical regions, including the prefrontal, temporal, and cingulate cortices with some brain regions, such as the cingulate and temporal cortices exhibiting greater transcriptional vulnerability to SZ than others.4 Here we report a WGCNA in SZ from a large gene expression microarray data set in 4 different brain regions. Discrete modules of coexpressed genes associated with SZ were identified, including modules enriched in mitochondria, microglia, GABAergic, glutamatergic, neuronal, and oligodendrocyte-related genes. Using multiple published genome-wide association studies (GWASs), we found that only the neuronal (GABAergic and glutamatergic) and oligodendrocyte modules are significantly enriched for SZ genetic risk variants, providing additional support for more direct involvement of these gene expression networks in SZ. Furthermore, regional analysis suggested that the gene expression patterns that typically discriminate between brain areas in controls are attenuated in SZ.
Brain tissue specimens were obtained from the Brain Bank of the Department of Psychiatry of the Mount Sinai School of Medicine (New York, New York)/James J. Peters Veterans Affairs Medical Center (Bronx, New York). We analyzed postmortem brain tissue samples from 21 persons with SZ and 19 controls (eTable 1). Gray matter samples from the dorsolateral prefrontal cortex (Brodmann area [BA] 46), the middle temporal area (BA 21), the temporopolar area (BA 38), and the anterior cingulate cortex (BA 32) were obtained using stringent diagnostic, neuropathological, and RNA integrity inclusion-exclusion criteria and tissue-handling procedures described in detail elsewhere5- 8 (eTable 1).
DNA was extracted from 50 mg of human superior temporal gyrus brain tissue. Genomic analysis with Affymetrix 6.0 whole-genome genotyping assays confirmed the absence of large-scale genomic structural defects such as aneuploidy or copy number variants that have been directly associated with SZ9 (eAppendix). Thus, by multiple measures, these brains did not exhibit signs of neuropathological or genetic abnormalities.
Total RNA was extracted from 50 mg of frozen tissue, as described in detail elsewhere.4,10,11 Microarray analysis was performed using Affymetrix Human Genome U133 Plus 2.0 arrays, according to the manufacturer's protocol and analyzed as described elsewhere4,10,11 (eAppendix). Similarly prepared aliquots from the superior temporal gyrus were used in quantitative polymerase chain reaction analysis in a larger and independent set of brain tissue specimens using TaqMan probes and primer sets (Applied Biosystems) (eTable 2 and eAppendix).
Differential expression was assessed using the significance analysis of microarrays (SAM) package12 with the following parameters: a significance threshold false discovery rate of less than 0.05 and fold changes of greater than 1.3. Overall or regional-specific differential gene expression between persons with SZ and controls was assessed by an unpaired 2-class (SAM) test. Differences in gene expression among the different cortical regions was assessed by a multiclass (SAM) test (all 4 regions) or by a paired 2-class (SAM) test (paired comparison among 2 regions). The homogeneity of variance of gene expression was assessed using the Levene test.
Weighted gene coexpression networks were built using the WGCNA package in R (eAppendix). Network construction was performed using the “blockwiseModules” function in the WGCNA package.3 The expression levels of each module were summarized by the first principal component (the module eigengene). Differences among the modules were assessed with analysis of variance, paired or independent t test, depending on the comparison. P values were calculated based on 1000 permutations. A Bonferroni correction was applied, and α was set to .01. Network preservation statistics were calculated using the “modulePreservation” function in the WGCNA package.13 The composite statistics Zsummary and medianRank, which is independent of the module size, were used because both have been found to perform well in distinguishing between preserved and nonpreserved modules.13 Both Zsummary and medianRank statistics summarize density-based (Zdensity; medianRankdensity) and connectivity-based (Zconnectivity; medianRankconnectivity) preservation statistics. Density-based preservation statistics can be used to determine whether module nodes remain highly connected in the test network, whereas connectivity-based preservation statistics test whether the connectivity pattern between nodes in the reference network is similar to that in the test network. For the cell-type marker determination analysis, the “userListEnrichment” function in the WGCNA package, which measures the enrichment between inputted and premade collections of brain-related lists, was used. The top 1000 reciprocal within-module connections were graphically depicted using Cytoscape.14
Functional enrichment was assessed using the GeneGo MetaCore software. The statistical significance threshold level for multiple testing was adjusted using Benjamini-Hochberg false discovery rate analysis (P < .05).
All gene set overlap analyses were performed with the cumulative probability function of the hypergeometric distribution using the phyper function in R. The population size was defined as the total number of genes with evidence of robust expression (N = 20 420).
Genome-wide association enrichment analysis was performed using the gene set analysis (GSA) method, as implemented in the GSA-SNP software (http://gsa.muldas.org/).15 Only single-nucleotide polymorphisms (SNPs) located within the transcript boundaries and an additional 20 kilobase on the 5′ end and the 3′ end of the gene were used. A Benjamini-Hochberg multiple testing correction was used with α set at 10−4.
Postmortem brain tissue samples from 21 persons with SZ and 19 controls were analyzed (eTable 1). After filtering for high-quality array data, we retained 107 samples (from 53 persons with SZ and 54 controls) for further analysis (eTable 3). Patients with SZ had significantly longer postmortem intervals (eTable 3), but RNA quality (3′ to 5′ hybridization ratios of glyceraldehyde 3-phosphate dehydrogenase) was equivalent between the 2 groups and high.
An unpaired 2-class SAM test between SZ and control samples across all cortical regions identified 611 genes showing significant expression changes. A separate unpaired 2-class SAM test between SZ and control samples for each brain region demonstrated that gene expression in BAs 21, 32, and 38 were significantly altered in SZ, whereas changes in BA 46 were minimal (Figure 1A). In addition, the significances of overlap analyses were higher for BAs 21, 32, and 38 but not for BA 46 (Figure 1B), indicating that gene expression changes in this cohort of elderly patients with SZ are less pronounced in BA 46, similar to previous findings from our group using a different microarray collection.4 The results of these analyses were consistent with other databases derived from different cohorts and analyzed with identical methods (eAppendix).
Hierarchical clustering analysis of the top 200 differentially expressed genes between SZ and control samples showed distinct clustering of the majority of SZ cortex samples. More tight clustering was observed when only BAs 21, 32, and 38 were included (Figure 1C). Clustering was independent of age, sex, BA, or postmortem interval. Of 611 genes that were dysregulated in all cortical regions, 432 were downregulated in SZ samples, and 179 were upregulated. Gene ontology analysis showed that the upregulated genes were enriched for genes in the cytoskeleton and for genes implicated in the initiation of translation and the cell cycle, as has been shown previously4,10,11 (eTable 4). The downregulated genes in the SZ samples were enriched for gene ontology categories related to immune response and negative regulation of cell proliferation and cell adhesion, including markers of oligodendrocyte differentiation and axon × glia interaction as identified in previous studies4,7,10,16,17 (eTable 4). Overall, gene ontology analysis in each separate brain region or in the top 200 differentially expressed genes showed enrichment for categories related to transcription/translation, signal transduction, the cell cycle, cell adhesion, the immune response, apoptosis, and the cytoskeleton (eTable 4).
To identify discrete groups of coexpressed genes showing transcriptional differences between SZ and control samples, we constructed a coexpression network using the entire data set, composed of both SZ and control samples (eTable 5). Five network modules whose module eigengenes were significantly associated with disease status (but not associated with any of the potential confounding variables) were identified (eTable 6). The top module (OLIG1) showed highly significant enrichment for oligodendrocyte markers (eTable 5), as well as for genes belonging to cytoskeleton rearrangement, axonal guidance, synaptogenesis, and cell adhesion gene ontology categories (Figure 2). The gene context of this module, which was downregulated in the SZ brain, showed significant similarity to genes identified in previous studies of patients with SZ,4,7,10,16,17 including MOG, MAG, CNTN2, ERBB3, and CNP that were additionally identified as hub genes. These changes were more consistent for BAs 32 and 38 (Figure 2). The second module of coexpressed genes highly related to SZ disease status, MIT1, was enriched for genes related to the mitochondria, the neuron, and the glutamatergic synapse (eTable 5) and was overexpressed in the SZ brain, indicating that genes in this module were upregulated in the SZ brain. Pathway analysis indicated functional enrichment of MIT1 with gene categories involved in synaptic vesicle exocytosis, synaptogenesis, and axonal guidance. The next module of upregulated coexpressed genes in SZ, GABA, was enriched with genes associated with GABAergic neurons and also with a gene cluster with unknown function (eTable 5). These changes were more consistent for BAs 21 and 38 (Figure 2). Functional analysis (MetaCore; GeneGo/Thomson Reuters) showed that these genes belong to nuclear receptors of transcriptional regulation, cholecystokinin signaling, and antiapoptosis mediated by the mitogen-activated protein kinase pathway (Figure 2). The next 2 significant modules, MG3 and MG5, were both downregulated in SZ samples and were enriched in genes related to microglia and astrocytes (eTable 5). MetaCore analysis showed enrichment for pathways related to the regulation of epithelial-to-mesenchymal transition, cytoskeleton remodeling, and immune response (MG3), as well as regulation of EIF4F activity and the IL-6 signaling pathway (MG5). Correlation analysis of modular membership with levels of gene significance for SZ samples showed that genes significantly associated with SZ are also hubs in GABA and OLIG1, but not in microglia or mitochondria modules (eAppendix).
We then applied module preservation statistics to determine which modules were least preserved among SZ samples relative to control samples. Network preservation statistics demonstrated that the GLU1 was the least preserved module (medianRank = 33; Zsummary = 5.2) (eFigure 1). The GLU1 module was enriched for markers related to the neuron and glutamatergic synapse (eTable 5), and MetaCore analysis showed functional enrichment for pathways related to protein folding and neurogenesis. The nonpreservation for GLU1 was applicable to both density-based (medianRankdensity = 33; Zdensity = 4.2) and connectivity-based (medianRankconnectivity = 32; Zconnectivity = 6.3) preservation statistics. Similar results were obtained when the WGCNA was repeated after expression traits were first adjusted for relevant covariates such as sex, age, and postmortem interval using multiple linear regression analysis (eAppendix). These findings are illness-specific findings and are not secondary to the nature of tissue collection and the procedural and selection standards of the current brain bank because an identical analysis using the same cohort of controls and a group of patients with Alzheimer disease that were banked with identical procedural standards showed no overlap among the GABA, OLIG1, or GLU1 modules with any of the modules that were significantly associated with Alzheimer disease (eAppendix). Seven genes were randomly selected from the WGCNA, and gene expression changes were confirmed with quantitative polymerase chain reaction analysis in an independent and larger cohort (eAppendix).
To examine whether the SZ-associated transcriptional differences observed are likely to be primary changes (vs collateral, bystander, or environmentally induced changes), coexpression modules were tested for significant enrichment with SZ genetic association signals. Several GWASs examining a variety of illnesses (SZ,18 bipolar disorder,19 attention-deficit/hyperactivity disorder,20 and major depression disorder21) and 3 nonpsychiatric disease studies as negative controls in human immunodeficiency virus,22 psoriasis,23 and Parkinson disease24 were included. The GABAergic, oligodendrocyte, and glutamatergic modules showed highly significant enrichment for association signals (P < 10 × 10−4) in the SZ GWAS, but none of the microglia (MG3 and MG5) or mitochondria (MIT1) modules showed such enrichment (Figure 3). Of the 2179 total genes included in the significantly enriched modules, only 9 (4 for SZ samples and 5 for control samples) overlapped with copy number variants, indicating that the differential expression levels are unlikely to be due to differences in copy number variants.
To assess whether there are global differences in the organization of the brain transcriptome between persons with SZ and controls, separate coexpression networks were constructed by applying WGCNA. Both the control and SZ brain networks showed high similarity with previously described human brain coexpression networks (eTable 7 and eTable 8), consistent with the existence of robust modules of coexpressed genes related to specific cell types and biological functions.
Regional differences among the modules were assessed with an analysis of variance comparison of the module eigengenes. Three of the control modules (cUNK1, cGLU1, and cNEU3) and 2 of the SZ modules (sUNK1 and sNEU3) related to neurodevelopmental, pyramidal, and glutamatergic genes showed significant differences among the 4 cortical regions after Bonferroni correction (all P < 10−5). The cUNK1 module overlapped significantly with the sUNK1 module (P < 6 × 10−176, determined by use of a hypergeometric test), and both modules overlapped with modules from the premade set of brain-derived enrichment genes that have been extensively described previously25,26 but have an unknown function (eTable 7 and eTable 8). Pathway analysis indicated that both modules are enriched in genes related to retinoic acid, Ran regulation, and Hedgehog signaling pathways. The cNEU3 module overlapped significantly with the sNEU3 module (P < 3 × 10−15, determined by use of a hypergeometric test), and both modules overlapped with premade modules enriched in genes related to glutamatergic and pyramidal neurons27 (eTable 7 and eTable 8). The cGLU1 module overlapped with premade modules enriched in genes related to glutamatergic neurons and neurogenesis25,26,28 and to the sGLU1 module (P < 10−10, determined by use of a hypergeometric test).
The glutamatergic/neurogenesis-related module showed differential expression by region in control samples but not in SZ samples. The regional expression patterns of the module eigengenes were significantly altered among the cUNK1 and sUNK1 modules, as well as among the cNEU3 and sNEU3 modules, whereas the sGLU1 module, which overlaps with the cGLU1 module, showed dampened regional differences (eFigure 2). Differences between groups, as well as between brain regions, were most pronounced in BA 32 (eFigure 2).
Based on the more prominent interregional coexpression network differences found in the control samples, SAM tests were used to compare gene expression levels that differentiated the 4 brain regions between SZ samples and control samples. Remarkably, whereas 681 genes were differentially expressed among BAs 21, 32, 38, and 46 in the control samples (a false discovery rate of <1%), only 82 genes were differentially expressed in the same regional comparisons among the SZ samples (Figure 4A). This was not simply an issue of statistical thresholds because relaxing the statistical criteria for differential expression to a false discovery rate of 5% identified 1178 differentially expressed genes in control samples and only 247 differentially expressed genes in SZ samples, confirming the large difference (>75%) observed in regional cortical differential gene expression between SZ samples and control samples. The homogeneity of gene expression variance across the SZ and control groups was analyzed using the Levene test. Sixty-one genes showed a significant difference in variance (P < .05) between SZ and control groups, indicating that increased variance was not the major factor responsible for the striking difference in regional gene expression between SZ samples and control samples.
Pathway analysis showed that the 681 genes that were differentially expressed between regions in control samples were enriched in genes related to neurogenesis/synaptogenesis and cell adhesion. Seventy of 82 differentially expressed between-region genes in SZ samples overlapped with the regionally differentially expressed genes in control samples (Figure 4A). Pathway analysis showed that unique genes in the control cluster that did not overlap with any of the SZ cluster genes were enriched in genes related to transcription (nuclear receptors of transcriptional regulation: unique control genes, P = 1.3 × 10−9; unique SZ genes, P = .99; common genes between control and SZ clusters, P = 3.3 × 10−3). Forty-two of the differentially expressed genes in the control samples, but none of the regionally different genes in the SZ samples, overlapped with the 611 genes that showed significant expression changes between SZ samples and control samples in all cortical regions (Figure 4A and B). These 42 probes were enriched for genes related to the cytoskeleton and the cell cycle. Interestingly, 34 of 42 genes were part of the unknown function module in control samples (cUNK1; P = 2.7 × 10−36, determined by use of a hypergeometric test), and 3 of the 42 genes (SCARNA9, SMC3, and TOP1) were intramodular hub genes (Figure 4C). In the SZ-associated network, only 21 of the 42 probes were localized to the corresponding module (sUNK1), with 2 of the genes (SCARNA9 and SMC3) playing a central role. Separate regional paired SAM tests and WGCNA revealed that the more prominent regional differences in control samples compared with SZ samples were unrelated to the inclusion of subjects with incomplete data for all regions of interest (eAppendix). Similar results were obtained when analysis was limited to subjects with microarray data for all 4 BAs (eAppendix).
Moving from gene lists to function and brain circuits remains a challenge in most microarray studies, especially when the focus is on gene-by-gene analysis of differential expression in single brain regions, which relies on a priori biological categorizations to provide notions of biological structure and regional significance to disease. Network analysis moves beyond single-gene investigations to provide a system-level understanding of the relationships between members of a network.29 Using a multiregional, system-level analysis of the aged SZ brain transcriptome, we found the existence of convergent molecular abnormalities in SZ, revealing a more comprehensive molecular neuropathological basis for the disease than the single brain region transcriptomic studies or GWASs of the past. Although the cortical transcriptome of SZ and control networks is organized into distinct modules of coexpressed genes that overlap significantly with previously described modules,25,26 significant SZ-associated differences were found, including (1) significant differences in the expression levels of distinguishing module eigengenes, which also show specificity enrichment for SZ genetic risk factors but not other neuropsychiatric diseases; (2) a lack of preservation of coexpression relationships in a module, suggesting that major perturbations in gene connectivity and density exist among persons with SZ and controls; and (3) regional differences in gene expression patterns that are found in control samples but are either altered or absent in SZ samples. Although the present findings complement earlier observations4,7,16,17 regarding the dysregulation of genes such as MOG, MAG, CNTN2, ERBB3, and CNP, they underscore the importance of these changes by identifying them as hub genes within the oligodendrocyte module and place them in the context of other systems (glutamate and GABA) that have long been associated with SZ.
With the exception of 2 previous studies,4,30 multiregional genome-wide studies of gene expression in SZ have not been extensively reported. As reported previously,4 the dorsolateral prefrontal cortex, which is a frequent target of SZ studies, showed the lowest number of altered transcripts. On the other hand, the cingulate and temporal cortices were among the most profoundly affected brain regions examined with respect to overall dysregulation of gene expression.4 These findings are consistent with multiple meta-analyses of neuroimaging data, in which the most significant alterations are observed in the insula, temporal gyrus, and the anterior cingulate gyrus/medial frontal cortex.31- 33 However, compared with the other regions, the matching of samples for age and pH was poorer for BA 46, which could add a confounder in our analysis. Furthermore, it remains unclear whether BA 46 demonstrates more significant gene expression differences at earlier stages in disease that are attenuated with aging and not captured in our population.
The described unbiased analyses have identified significant and unique features of the transcriptomic network in SZ, and they have, at the same time, confirmed previously implicated gene expression abnormalities in SZ (oligodendroglial, GABAergic, and glutamatergic dysfunction) while suggesting that other changes, such as those associated with mitochondria and microglia, may be less disease specific. However, classification of primary and secondary changes in complex brain disorders, especially when derived from cross-sectional analyses and a homogenized whole cortex, is somewhat artificial, and one cannot be certain which set of changes contribute to the symptoms directly. Therefore, these data should be interpreted cautiously, and they do not diminish the importance of immune system activation and metabolism in SZ as potential therapeutic targets. As is the case with most studies in SZ, the current observations could be due to changes induced in gene expression by prolonged exposure to antipsychotic medications or nicotine and/or as part of the disease process. Although this possibility cannot be excluded with certitude, it is difficult to attribute all of the observed differences, including the associations of the oligodendroglial, GABAergic, and glutamatergic modules with SZ-specific GWASs, to these confounders. Another limitation is the chronic hospitalization and advanced age of the included subjects. Although this strategy minimizes potentially confounding factors such as dependence on alcohol and illicit substances, ambiguous or violent circumstances of death, and neurological or neuropsychiatric comorbidities, it only provides insight for later stages of the disease. Finally, the fact that the group of persons with SZ included more males and younger subjects could not completely exclude the possible confounding effect of age and sex on the reported results.
Gene expression patterns that typically distinguish brain regions were found to be significantly attenuated in the SZ brain, suggesting abnormalities in cortical circuit patterning. A subset of genes differentially expressed regionally in controls overlapped with genes that were significantly dysregulated among persons with SZ and controls. The majority of them are implicated in cellular processes that are related to DNA repair, RNA processing, and alternative splicing and were organized in a distinct module that showed regional differences in controls. Three of these genes were also intramodular hubs (SCARNA9, SMC3, and TOP1). The TOP1 gene encodes an essential enzyme involved in resolving the torsional stress associated with DNA replication, transcription, and chromatin condensation.34SMC3 is a component of the multimeric Cohesin complex that holds together sister chromatids during mitosis, enabling proper chromosome segregation.35SCARNA9 function remains unknown, although it has been suggested to be part of the Cajal bodies that play an important role in RNA processing.36 Overall, these data suggest that genes that are significantly associated with SZ and show typical regional differences, many of which are related to gene transcription and translation, are organized in distinct modules and that their expression levels are attenuated in the neocortex of aged patients with SZ. Such changes in the regional pattern of this module might be a potential pathophysiological driver leading to dampened regional differences in the gene expression patterns of aged patients with SZ and disruption of the networks that they subserve. Interestingly, a similar analysis conducted for autism spectrum disorder reported a similarly attenuated pattern of gene expression among frontal and temporal regions,37 and this similarity argues against medication effects. Although our findings are directly applicable to aged subjects, we speculate that such regional differences will be evident at early stages of the disease, affecting the developmental trajectory of modules related to myelination and synaptic neurotransmission. Overall, these findings suggest that abnormalities in cortical gene expression patterning may be common in neurodevelopmental disorders, such as autism and SZ, and may point to a potential pathophysiological substrate that, when combined with other disturbances, leads to different disease phenotypes.
There is strong evidence supporting the notion that SZ is caused by the developmentally regulated disconnection of higher-order association areas, which has led to the formulation of the disconnectivity hypothesis of SZ.17,38 Disconnectivity could result at a microstructural level (microconnectivity), as a result of altered developmental local synaptic connectivity or as deficits in macroconnectivity between cortico-cortical and cortico-subcortical loops. At a microstructural level, the disconnection could be due to aberrant wiring of association fibers during brain development, aberrant interactions between different neuronal elements (eg, glutamate and GABA neurons), impairments in synaptic plasticity, or all of the above. Similarly, disruption in macroconnectivity could be caused by abnormalities in axon guidance during development or by deficits in myelination, myelin function, neuron-myelin interactions, or any combination of the above. The present study provides evidence supporting abnormalities in both microconnectivity and macroconnectivity that share a common genetic cause and might function as the drivers for the disease. At the same time, changes in molecules historically associated with immune system induction and energy production that are not directly related to the disconnectivity hypothesis do not show enrichment for SZ-related genetic risk factors. Importantly, disconnection among different regions might lead to decoupling and diminished regional differences in the expression patterns of genes (Figure 5).
In summary, we have used WGCNA to elucidate the organization of gene expression, in persons with SZ and matched healthy controls, into a hierarchical, scale-free network containing modules of highly coexpressed genes. These modules correspond to multiple aspects of cellular and molecular function that have been systematically validated using published large-scale proteomic and genomic data sets. These data strongly suggest that networks related to synaptic neurotransmission, oligodendrocytes, and neuronal genes are crucial for SZ pathophysiology. Future studies should investigate whether these alterations are observed in earlier stages of the disease, as well as the effect of intramodular, hub genes on specific network functions. The identification of such hub genes and networks has the potential to lead to novel targets for the development of new treatments.
Correspondence: Dr Haroutunian, Department of Psychiatry, The Mount Sinai School of Medicine, Room 4F-33, 130 W Kingsbridge Rd, Bronx, NY 10468 (firstname.lastname@example.org).
Submitted for Publication: January 27, 2012; final revision received April 29, 2012; accepted May 3, 2012.
Published Online: August 6, 2012. doi:10.1001 /archgenpsychiatry.2012.704
Author Contributions: Drs Roussos and Haroutunian had full access to all of the data in the study and take responsibility for the integrity of the data and the accuracy of the data analysis.
Financial Disclosure: None reported.
Funding/Support: This work was supported by National Institutes of Health grants MH066392 and MH064673 (to Dr Haroutunian) and by Veterans Administration Mental Illness Research, Education, and Clinical Centers and Merit awards (to Drs Siever and Haroutunian, respectively).