Interpreting Coronary Artery Disease Risk Through Gene-Environment Interactions in Gene Regulation

Abstract

GWAS and eQTL studies identified thousands of genetic variants associated with complex traits and gene expression. Despite the important role of environmental exposures in complex traits, only a limited number of environmental factors were measured in these studies. Measuring molecular phenotypes in tightly controlled cellular environments provides a more tractable setting to study gene-environment interactions in the absence of other confounding variables. We performed RNA-seq and ATAC-seq in endothelial cells exposed to retinoic acid, dexamethasone, caffeine, and selenium to model genetic and environmental effects on gene regulation in the vascular endothelium-a common site of pathology in cardiovascular disease. We found that genes near regions of differentially accessible chromatin were more likely to be differentially expressed [OR = (3.41, 6.52), [Formula: see text]]. Furthermore, we confirmed that environment-specific changes in transcription factor binding are a key mechanism for cellular response to environmental stimuli. Single nucleotide polymorphisms (SNPs) in these transcription response factor footprints for dexamethasone, caffeine, and retinoic acid were enriched in GTEx eQTLs from artery tissues, indicating that these environmental conditions are latently present in GTEx samples. Additionally, SNPs in footprints for response factors in caffeine are enriched in colocalized eQTLs for coronary artery disease (CAD), suggesting a role for caffeine in CAD risk. By combining GWAS, eQTLs, and response genes, we annotated environmental components that can increase or decrease disease risk through changes in gene expression in 43 genes. Interestingly, each treatment may amplify or buffer genetic risk for CAD, depending on the particular SNP or gene considered.

Keywords: TWAS; gene expression response; gene regulation; genetic variation; transcription factor.

Figure 1

Changes in gene expression and chromatin accessibility in HUVECs treated with dexamethasone, retinoic acid, caffeine, and selenium. (A) Diagram of artery layers and study design. HUVECs were treated separately with four compounds and vehicle controls prior to RNA- and ATAC-seq (Materials and Methods). (B) QQ-plot of P-values from DESeq2 analysis of gene expression in response to retinoic acid. The purple dots represent genes within 50 kb of a DAR, while the gray dots represent genes not within 50 kb of a DAR. (C) Enrichment of DE genes within 50 kb of a DAR. Odds ratio was estimated using Fisher’s exact test, and for all treatments P$<2.2\times {10}^{-16}$.

Figure 2

Shared DE genes and DARs across treatments. (A) DE gene sharing. Vertical bars indicate the number of genes differentially expressed in each group of treatments. Inset heatmap depicts the Spearman correlation of the logFC of all DE genes between treatments. (B) DAR sharing. Vertical bars indicate the number of DARs in each group of treatments. Inset heatmap depicts the Spearman correlation of the logFC of all DARs between treatments.

Figure 3

Identification of response factors. Topmost enriched transcription factor footprints in DARs in each condition. For selenium, the top factors are enriched in regions of closing chromatin. For all other treatments, the top factors are enriched in regions of opening chromatin.

Figure 4

Latent environments in GTEx and CAD risk. (A) Enrichment of GTEx eQTLs in SNPs within response factor footprints for each treatment. (B) Enrichment of coronary artery disease risk loci in SNPs within response factor footprints for each treatment. (C) For the genes shown in (D), number of concordant (treatment changes gene expression in the same direction as TWAS risk) and discordant (vice versa) gene/treatment pairs. (D) Comparison of logFC of gene expression vs. CAD TWAS z-score. A positive logFC indicates that the treatment increases gene expression. A positive TWAS value indicates that increased expression of the gene is associated with CAD risk.

Figure 5

reQTLs. (A) Boxplot depicting logFC for selenium reQTL for MAT2A (B) Same reQTL as in (A), but with trendlines representing gene expression in the treatment and control conditions. (C) All reQTLs. Colors indicate gene expression in the treatment condition (blue, dexamethasone; red, caffeine; orange, selenium; purple, caffeine; black, control).

Figure 6

reQTLs modulate environmental effects on complex traits. (A) Table relating TWAS, gene expression, and reQTL effects. (B) Risk associated with treatment vs. genetic modulation by reQTLs. The risk associated with treatment represents the TWAS effect multiplied by the sign of the logFC of the treatment. Positive values indicate that the treatment amplifies CAD risk. Genetic modulation represents the reQTL effect multiplied by the sign of the logFC of the treatment. Positive values indicate that the risk allele amplifies the treatment effect.