Characterizing the genetic basis of methylome diversity in histologically normal human lung tissue

Abstract

The genetic regulation of the human epigenome is not fully appreciated. Here we describe the effects of genetic variants on the DNA methylome in human lung based on methylation-quantitative trait loci (meQTL) analyses. We report 34,304 cis- and 585 trans-meQTLs, a genetic-epigenetic interaction of surprising magnitude, including a regulatory hotspot. These findings are replicated in both breast and kidney tissues and show distinct patterns: cis-meQTLs mostly localize to CpG sites outside of genes, promoters and CpG islands (CGIs), while trans-meQTLs are over-represented in promoter CGIs. meQTL SNPs are enriched in CTCF-binding sites, DNaseI hypersensitivity regions and histone marks. Importantly, four of the five established lung cancer risk loci in European ancestry are cis-meQTLs and, in aggregate, cis-meQTLs are enriched for lung cancer risk in a genome-wide analysis of 11,587 subjects. Thus, inherited genetic variation may affect lung carcinogenesis by regulating the human methylome.

Figure 1. cis-meQTL structural characteristics

(a) Distribution of CpG probes and corresponding cis-meQTL numbers and proportions in gene and non-gene regions. meQTLs were detected based on EAGLE lung normal tissue samples (n=210). (b) Distribution of CpG probes and corresponding cis-meQTL numbers and proportions in CpG islands (CGIs), shores (< 2kb from the boundary of CGI), shelves (2–4kb from the boundary of CGI) and the remaining region or “open sea”. The box plots show the distribution of the methylation levels in each genotype category with error bars representing the 25% and 75% quantiles. (c) The strongest cis-association is between SNP rs10090179 and CpG probe cg19504605. P=1.5×10⁻⁷³, t-test. The SNP explains 79.8% of the phenotypic variance. (d, e): The x-coordinate is the average standard deviation (SD) of methylation levels for CpG probes in each category. The y-coordinate is the proportion of CpG probes detected with cis-meQTLs. The proportion of methylation probes detected with cis-meQTLs varied across categories, ranging from 4.0% for CGIs in 1st Exons to 15.7% for south shores in non-gene regions.

Figure 2. trans-meQTL structural characteristics

(a) Circos plot for trans-meQTLs. The outer rim shows the log₁₀ P-values Manhattan plots of trans-meQTL associations. The innermost network depicts spokes between all trans-meQTL SNPs and their target CpG sites. The red spikes show a master regulatory SNP rs1293229 located at 16p11.2 associated with methylation of CpG sites located in CGIs annotated to five genes. (b) Proportion of CpG probes detected with cis-meQTLs and trans-meQTLs across gene regions. The asterisks “*,**,***” indicate t-test P<0.05, 0.01, and 0.0001 for the comparison between CGI and non-CGI regions. CGI regions are strongly enriched with trans-meQTLs, while non-CGI regions are enriched with cis-meQTLs. CpG-sites in 3’UTR regions show an opposite trend. (c) The association between a SNP denoted as G and a distal CpG-site B may be mediated through a proximal CpG-site A. (d) For each trans-association (G, B) pair, the dots show their marginal v.s. partial correlation coefficients upon conditioning on the proximal A CpG probes. Analysis was based on 210 samples. Reduction of correlation coefficients by conditioning on A suggests the magnitude of the mediation effect.

Figure 3. Chromatin marks are increasingly enriched on meQTL SNPs with larger effect sizes

(a) We split cis-meQTL SNPs into five categories according to the meQTL association strength (P>10⁻⁷, 10⁻⁷>P>10⁻¹⁰, 10⁻¹⁰>P>10⁻¹⁵, 10⁻¹⁵>P>10⁻²⁰, P<10⁻²⁰). A SNP is determined to be related with a regulatory region if the SNP or any LD-related SNP (r² ≥ 0.8) resides in the ChIP-Seq peaks of the regulatory regions. Regulatory elements include CTCF binding sites, DNaseI hypersensitive sites and histone marks from small airway epithelial cells (SAEC) from ENCODE and human alveolar epithelial cells (hAEC) from our laboratory. For each p-value category, we calculated the proportions of cis-meQTL SNPs related with regulatory regions. The figures show that the proportions of cis-meQTL SNPs related with regulatory regions increase with the significance of meQTL associations except for the repressive mark H3K27me3.

Figure 4. DNA methylation regional associations for lung cancer GWAS SNPs in subjects of European ancestry

(a, b, f and g) Symbols represent the association between established lung cancer GWAS genetic loci in four regions and methylation levels in nearby CpG probes. Y-coordinate, P-value for association; x-coordinate, genomic location. For each SNP, the red solid circle or square represents the methylation probe with the strongest association, whereas other methylation probes are colored on the basis of their correlation (measured as r²) to the most-associated probe. For the most-associated probes, the P-values in EAGLE discovery set (n=210) and TCGA lung replication data (n=65) are shown. SNP locations are marked by a blue triangle. (c–e and h–j) show the associations between genotypes and methylation levels of the most associated CpG probes. The box plots show the distribution of the methylation levels in each genotype category with error bars representing the 25% and 75% quantiles.

Figure 5. Enrichment of cis-meQTL SNPs for lung cancer risk

Analysis based on NCI lung cancer GWAS data (5,739 cases and 5,848 controls). P-values were produced based on 10,000 permutations. AD, SQ, and SC represent adenocarcinoma, squamous cell carcinoma and small cell carcinoma. (a) Enrichment was tested using all cis-meQTL SNPs after LD pruning. (b and c) Strong enrichments were observed for cis-meQTL SNP associated with CpG probes annotated to north shores (b) and gene body (c) regions for SQ. (d) The enrichment in (c) was driven by the cis-meQTLs SNPs impacting CpG probes in non-CpG islands. (e) The enrichment in (d) is driven by the SNPs (or their LD SNPs with r² > 0.95) overlapping with CTCF binding sites or H3K27me3 mark regions.