Genetic Differences in the Immediate Transcriptome Response to Stress Predict Risk-Related Brain Function and Psychiatric Disorders

Abstract

Depression risk is exacerbated by genetic factors and stress exposure; however, the biological mechanisms through which these factors interact to confer depression risk are poorly understood. One putative biological mechanism implicates variability in the ability of cortisol, released in response to stress, to trigger a cascade of adaptive genomic and non-genomic processes through glucocorticoid receptor (GR) activation. Here, we demonstrate that common genetic variants in long-range enhancer elements modulate the immediate transcriptional response to GR activation in human blood cells. These functional genetic variants increase risk for depression and co-heritable psychiatric disorders. Moreover, these risk variants are associated with inappropriate amygdala reactivity, a transdiagnostic psychiatric endophenotype and an important stress hormone response trigger. Network modeling and animal experiments suggest that these genetic differences in GR-induced transcriptional activation may mediate the risk for depression and other psychiatric disorders by altering a network of functionally related stress-sensitive genes in blood and brain.

Figure 1

Summary Figure Illustrating the Sequence of Experiments and Analyses Applied in This Study The main hypothesis tested in this study is that common genetic variants that alter the short-term transcriptional response to GR activation also alter the risk for stress-related psychiatric disorders and related neural endophenotypes.

Figure 2

GR-Response-Modulating cis-eQTLs (A) Study design for GR-stimulated gene expression in whole blood of 160 male individuals from the Max Planck Institute of Psychiatry cohort. (B) Circularized Manhattan plot displaying cis-associations for GR-response eQTL bins (n = 320) and their respective significance (−log₁₀ Q values). Displayed from the outer to the inner circle are the number of chromosomes, the ideograms for the human karyotype (hg18), genes nearby eSNPs, and Manhattan plots for the eQTL bins that survived correction for multiple testing. (C and D) Boxplots of human gene expression values for ADORA3, which is an example of a significant GR-response eQTL. Expression levels are stratified based on the eSNP genotypes for ADORA3. Baseline (6 p.m.) measures are displayed in blue and GR-stimulated measures (9 p.m.) in red. Microarrays data are displayed in (C) and their qPCR validation in (D). Q value in (C) is derived from GR-response cis-eQTL analysis and the p value in (D) from the qPCR linear regression model.

Figure 3

GR-Response eSNPs Are Enriched in Enhancer Regions in Multiple Tissues Bar graph illustrating the enrichment of GR-response eSNPs for enhancers in multiple tissues from the Roadmap Epigenome Project, including brain tissue. The x axis shows the fold enrichment and the y axis all brain enhancers all well as the mean fold enrichment among all hematopoietic cells (see Figure S2) and brain enhancers. The fold enrichment for GR-response eSNPs is illustrated in red and for the permuted baseline eSNPs in blue. Only the GR-response eSNP enrichment, which passes a Bonferroni corrected significance threshold (corrected for the number of all tested tissues or cells, n = 62) is illustrated. ^∗ p ≤ 0.05, obtained by binomial enrichment test and Bonferroni correction, error bars ± SD.

Figure 4

Long-Range Chromatin Interaction of GR-Response eQTLs (A) Long-range chromatin interaction as exemplified by the eSNP region containing the NRTN locus (chr10: 5,690,000–5,840,000; hg19) was confirmed by 3C in five lymphoblastoid cell lines (LCLs) each, homozygous for the two opposite SNP alleles, both in the presence and absence of dexamethasone. A SNP in the NRTN locus (rs1379868) affects the differentially regulated gene expression of LONP1 in human whole blood cells (based on GR-response eQTL analysis). Baseline (6 p.m.) measures are displayed in blue and GR-stimulated measures (9 p.m.) in red. (B) SNP effect on GR-dependent gene transcription was validated by qPCR in the LCLs used for the 3C assay. (C) Characterization of the eSNP locus. Top panel, ideogram for chromosome 19 (p13.3). A red box isolates the region shown (enlarged) in the bottom panel. Bottom panel: 3C-primers (green track) were designed at the LONP1 TSS (C1, anchor) and multiple regions (P1–P6) in and around the eSNP bin. The eSNP bin includes a GR binding site in blood cells (pink track). ChIA-PET tags from the leukemia cell line (brown and green tracks) validate a direct chromatin interaction between the NRTN eSNP locus and the regulated gene LONP1. The paired ChIA-PET tags coincide with DNaseI hypersensitivity sites in the leukemia cell line (red track) and blood cells (yellow track). (D) Chromatin conformation capture interaction data. A 3C physical interaction between the LONP1 TSS and eSNP bin (P4), emphasized by a gray box, was found in the 3C libraries made from LCLs (p = 3.35 × 10⁻²³, χ² = 115.15) with a stronger interaction following stimulation with the GR-agonist (p = 0.06, χ² = 3.35). Q values in (A) are derived from GR-response cis-eQTL analysis, and p values in (B) and (D) are derived from linear mixed model; error bars ± SD.

Figure 5

GR-Response eSNPs Are Enriched among Variants Associated with MDD (A) The dotted red line shows the enriched number of GR-response eSNPs that overlap with SNPs in our meta-analysis for MDD (= MDD-related GR eSNPs; 8,864 cases and 8,982 controls). The distribution of the observed overlap for sets of 1,000 random SNPs (gray) and 1,000 random baseline eSNPs (blue) are represented as histograms (null distributions). Both permuted data sets never reached the same overlap with MDD-associated SNPs as the GR-response eSNPs. (B) The distribution of the MDD-related GR eSNP genetic risk profile scores (GRPSs) for an independent sample of MDD cases (n = 1,005 cases; red) and controls (n = 478; gray) are represented as density plots. Individuals with MDD display higher GRPSs (p = 0.00017). p value by logistic regression model.

Figure 6

Functional Annotation of Transcripts Regulated by MDD-Related GR Risk Variants (A) Gene network produced using GeneMANIA. The network consists of 43 genes (circles) connected by 164 interactions (edges). Genes that are within a black filled circle indicate our MDD-related GR transcripts (n = 24), while those within a gray filled circle indicate additional genes (n = 20). The interactions found between these genes, which were more enriched than expected, are shown (co-expression: purple lines, shared protein domains: yellow lines, and co-localization: blue lines). (B) Heatmap of gene expression changes (log₂) between stress versus vehicle groups of mouse in brain and blood (n = 17 mice, left panel) as well as between baseline and GR-stimulation in human blood cells (blue, middle panel) and in mouse brain and blood (n = 22 mice, right panel). Investigated tissues are labeled within the bottom row of the heatmap (prefrontal cortex [PFC], hippocampus [HC], and amygdala [AM]). p values were computed by using linear regression model, and significance is indicated by a black box (FDR ≤ 0.1, dotted box p ≤ 0.05).

Figure 7

GR-Response MDD-Related eSNP GRPS Correlate with Overgeneralized Amygdala Reactivity (A) Statistical parametric map illustrating left centromedial amygdala reactivity to facial expressions with an “Angry & Fearful > Neutral” contrast in the entire sample (15 contiguous voxels; max voxel MNI coordinate, x = −24, y = −10, z = −14, t = 4.35, p = 7.76 × 10⁻⁶). (B) Higher MDD-related GR eSNP genetic risk profile scores (GRPSs) in the European-American subsample of the DNS cohort (n = 306) predicted amygdala reactivity to threat-related facial expressions in comparison to neutral facial expressions. (C and D) Post hoc analyses revealed that GRPSs did not predict amygdala reactivity to threat-related expressions (C), but that higher GRPSs predicted elevated amygdala reactivity to neutral facial expressions (D) in comparison to non-face control stimuli. The 95% confidence interval is displayed as gray shaded band in (B)–(D).