Genetic and epigenetic regulation of gene expression in fetal and adult human livers

1The liver plays a central role in the maintenance of homeostasis and health in general. Given the substantial inter-individual variation seen in metabolism, regulation of nutrients, protein synthesis, and detoxification of xenobiotics. It is essential to have a better understanding on inter-individual variation of gene expression, methylation and genetic effects specific to liver, and on different conditions, e.g. developmental stages. These variations can affect the liver’s metabolic properties, leading to high levels of metabolites, either in the forms of lipids, proteins or xenobiotics, which can result in serious diseases or toxic side-effects. For example, several single nucleotide polymorphisms (SNPs) associated with liver function and related diseases have been identified through genome-wide association (GWA) studies. We and others have studied how these SNPs affect liver gene expression levels by mapping expression quantitative trait loci (eQTL), and several genetic variants that regulate genes involved in the absorption, distribution, metabolism and excretion of drugs (ADME genes) have also been identified.

2Apart from genetic variation, epigenetic mechanisms (DNA methylation and histone modifications) also play an important role in regulating tissue-specific gene expression. In particular, such mechanisms can influence the expression of hepatic ADME genes. For example, the methylation status of a CpG island in exon 2 of CYP1A2 was shown to correlate with interindividual differences in the expression of this gene in human livers. Given that CYP1A2 is an important drug-metabolizing enzyme, those factors that influence its epigenetic state may also contribute to the individual drug response. Interestingly, the epigenetic state of ADME genes, at least in rodent livers, can change in response to xenobiotic exposure, thus opening the perspective for epigenetics-mediated drug-drug interactions. More examples on epigenetic regulation of ADME genes have been reviewed by Kacevska et al.. However, the majority of such data come from studies of epigenetic alterations observed either in tumors, or in cell lines treated with DNA demethylating agents. So far it is not clear, to which extent such cancer-related or experimentally induced epigenetic alterations correspond to the natural epigenetic variability in human livers. Hence, it is essential to include epigenetic variation when studying the regulation of hepatic gene expression, to further explain the causes of differences in drug response and the etiology of diseases associated with liver function.

3Here we present a comprehensive survey of the methylome and transcriptome of the human liver (Figure 1A). First, we addressed the regulation of gene expression in the developing human liver by comparing genome-wide expression and methylation levels in 96 adult and 14 fetal livers from the Karolinska Liver Bank. Then we used genetic, epigenetic and gene expression data from the adults, along with an extra cohort of 85 Dutch adult liver samples to investigate the regulation of gene expression in the human liver. Finally, we explored the tissue specificity of the identified associations between SNPs, methylation and expression in other tissues from the Dutch adult samples.

4We first compared the epigenomes of 8- to 21-week-old fetal livers with adult livers. We assessed the methylation levels of 366,074 variable CpG sites and found 28,917 CpG sites (annotated to 12,619 unique genes) that showed a significant difference (absolute mean beta value difference > 0.2, FDR < 0.05) between fetal and adult liver tissue (see Supplementary Online Methods in the Additional files 1 and 2). Although the number of hypomethylated CpG sites in fetal liver (53.4%) was similar to the number of hypermethylated sites (46.6%) in this cohort, we observed an age-specific association between the genomic location of CpG sites and whether they were hypo- or hypermethylated (chi-squared test p-value < 2.2 × 10-16). In fetal livers, the majority (86%) of the differentially methylated CpG sites that are located within CpG islands (CGI) were hypomethylated, whereas this was not the case for CpG sites outside CGIs, where roughly 50% of the CpG sites were either hypo- or hypermethylated in fetal livers (Figure 1B). This is particularly interesting because in both adult and fetal livers, close to 80% of the CpG sites within CGI are not methylated, with > 95% overlap between the two age groups. Accordingly, the CpG sites within CGIs that were hypomethylated in the fetal livers mostly had intermediate methylation levels in the adult liver samples (Additional file 3).

5To explore the functions of the genes that were differentially methylated in fetal liver compared to adult liver, we used the GREAT pathway tool. The CpG sites that were hypomethylated in the adult livers and hypermethylated in fetal liver were strongly enriched for metabolic pathways, such as the steroid metabolic process the regulation of lipid metabolic processes, regulation of generation of precursor metabolites and energy, and regulation of glycolysis (all with p-values < 1.15 × 10-44) (Table 1A). However, the genes that were associated with hypomethylated CpG sites in the fetal samples were strongly enriched for pathways of insulin receptor signaling, regulation of glycogen synthase activity, differentiation processes, and developmental functions (Table 1B).

6Comparison of gene expression levels between the fetal and adult liver samples yielded 3,284 differentially expressed probes (absolute log2-fold change > 1.0, FDR < 0.05, Additional file 4). Pathway analysis, using Gene Network, confirmed that 1,396 genes with higher expression in the adult livers were strongly enriched for metabolic functions like monocarboxylic acid, steroid and bile acid metabolic processes, as well as the response to xenobiotic process (Table 2A). In contrast, 1,277 genes that were highly expressed in fetal tissue were associated with regulating organelle organization, chromosome organization, and tetrapyrrole (e.g. hemoglobin) biosynthetic processes (Table 2B). These observations are in line with the fetal development, which is characterized by tissue differentiation and growth and by the fact that the liver is predominantly a hematopoietic organ during this period.

7We found 1,655 genes that showed both differential expression and differential methylation in adult vs. fetal livers (Additional file 5). More specifically, 657 genes were linked to probes with higher expression levels in adults, and 1,000 genes linked to probes that were more highly expressed in fetal livers (with an overlap of two genes). As expected, these genes are even more significantly enriched for developmental stage-specific functions, such as drug response for the adult cohort (p-value 4.0 × 10-131) and liver development for the fetal cohort (p-value 6.0 × 10-90). In the majority of the genes with more than one detection probe, the differences in expression levels were very similar between fetal and adult livers. However, in two genes (TGM2 and INS-IGF2), one of the probes was more highly expressed in fetal livers, while the other probe reflected higher expression in adult livers. The location of the differentially methylated CpG sites differed significantly in relation to CGIs, depending on the expression and methylation differences between fetal and adult livers (chi-squared test p-value < 2.2 × 10-16, Figure 1C): for genes with a lower expression in fetal livers, the hypomethylated CpG sites more often map within CpG islands, shores and shelves, while the hypermethylated CpG sites map further away from CGI regions.

8Regions within 2 kb of the transcription start site (TSS) of the 1,655 genes are enriched for binding sequences of transcription factors essential for the development or function of the liver, specifically HNF4A (adjusted p-value = 2 × 10-73) and HNF1A (adj. p-value = 6 × 10-38); hematopoietic transcription factors GATA1 (adj. p-value = 8 × 10-36), STAT5A (adj. p-value = 2 × 10-43), and STAT5B (adj. p-value = 1 × 10-49); and YY1 (adj. p-value = 2 × 10-36), which plays a fundamental role in embryogenesis and differentiation. We therefore investigated the expression levels of these transcription factors and observed that transcripts for the HNF1A and HNF4A genes were more highly expressed in adult livers, and GATA1, STAT5A, STAT5B and YY1 were all more highly expressed in fetal livers (Figure 2, Additional file 6).

9Table 3 lists the 20 genes with the largest differences in expression and methylation, clearly illustrating the fetal-specific expression of genes involved in differentiation and hematopoiesis (e.g. DLK1, HBZ, HBM, AHSP, EPB42 and NFE2), and the adult-specific expression of genes involved in drug metabolism, catabolism and other biosynthesis processes. CYP2E1 and CYP2C8 are the cytochrome P450 (CYP) genes; these show the most significant difference in expression levels between fetal and adult liver, with an approximately 7-fold higher expression level in adult liver.

10We next assessed whether DNA methylation is correlated to gene expression levels in the adult samples. We combined data from the Karolinska Liver Bank and Dutch liver samples (total number of samples with expression and methylation data = 158) and compared expression probes with CpG sites that map within 250 kb of these probes. We did not include the fetal samples due to the large developmental differences reported above, and we estimated that the fetal samples would not add any considerable statistical power for the analyses. We identified 3,238 significant methylation-expression associations (eQTMs, Additional file 7), comprising 1,988 unique expression probes (in 1,798 genes) and 2,980 CpG sites (reflecting 2,057 unique genes), with a permutation p-value < 0.05. As expected, there are more eQTMs with a negative correlation between expression levels and CpG methylation levels (58.4%), irrespective of the CpG site location in relation to CpG islands. Furthermore, for CpG sites with strong correlation between expression and methylation levels, and/or within 50 kb of the expression probes, we observed an overrepresentation of negative correlations (chi-squared test p-value < 2.2 × 10-16, Figure 3).

11We next explored the effects of genetic variation on liver gene expression levels. Expression quantitative trait locus (eQTL) mapping in the adult livers (meta-analysis of the two cohorts, combined number of samples with expression and genotype data = 171) yielded a total of 47,168 significant SNP-probe pair correlations (FDR < 0.05), representing 751 unique genes (Additional file 8). The eQTL probes are significantly enriched for liver-specific genes (area under the curve (AUC) 0.67, p-value 4 × 10-57, as reported by Gene Network) and are strongly enriched for genes encoding drug-metabolizing enzymes (p-value 2.0 × 10-19). We compared our results with reported liver cis-eQTLs and observed that we could replicate 667 reported eQTL genes, however we also identified 84 new eQTL genes (Additional file 8).

12We investigated the effects of SNPs on CpG methylation (meQTL) in adult liver samples (meta-analysis, combined number of samples with methylation and genotype data =161). In total we found significant cis-meQTL for 28,447 unique methylation probes (FDR < 0.05, mapping to 12,054 unique genes), reflecting 1,477,126 different SNP-CpG site combinations. In contrast to the eQTL, we did not observe any enrichment of liver functions for these 12,054 meQTL associated genes. Looking further into the SNPs affecting DNA methylation and gene expression, we identified 215 unique genes and 10,432 unique SNPs associated with both an eQTL and meQTL, resulting in a total of 30,644 overlapping QTL effects. Interestingly, for most of the 215 genes (69.3%) influenced by both an eQTL and meQTL we observed an opposite effect direction, i.e. the same genotype was associated with higher methylation levels and lower expression levels, or vice versa (Additional file 9). This effect is strongest in the CpG islands and CpG island shores, where it occurs in more than 75% of the cases (Additional file 10).

13Once we had identified eQTL and eQTMs, we ascertained to what extent SNPs and DNA methylation could jointly explain the variation in liver gene expression levels. We selected 293 expression probes (reflecting 274 unique genes) that had both a significant cis-eQTL and significant eQTM effect. We then tested four different linear models (see Supplementary Online Methods in the Additional file 1) to assess the proportion of variation in gene expression that could be explained. For 83% of these 293 expression probes, most of the expression variation was explained by a SNP (Additional file 11), whereas for the remaining 17% the expression variation was most strongly explained by a specific CpG site. For the latter cases, we observed that these expression-associated CpG sites were likely to have a meQTL effect (chi-squared p-value = 0.035). As expected, when we combined the SNP genotype and CpG site methylation levels, we could explain more of the expression variation than by using either SNP or methylation levels alone. Given the correlations between genotypes and methylation levels, we also estimated the unique contributions of the two on gene expression levels (Additional file 12). Overall, SNP genotypes uniquely explain a greater proportion of the variation in gene expression (median 0.1, standard deviation 0.122) than methylation levels (median 0.029, standard deviation 0.049). The SNPs and CpG sites with particularly high correlations with the expression levels were generally closer to the transcription start site of the corresponding genes (Additional file 13).

14The contributions of SNPs and DNA methylation levels to the proportion of variation explained in gene expression levels are illustrated in Additional file 14 and Table 4 by 16 unique ADME genes that had both significant eQTL and eQTMs. The ADME gene list was extracted from http://www.pharmaADME.org. For the GSTT1, GSTM1, UGT1A1, GST01 and PON1 genes, DNA methylation explains a larger proportion of the variation in gene expression levels compared to SNP genotypes. Overall, we found that adding more CpGs to the model, which were all associated to the selected expression probes of the same gene, did not significantly increase the power to explain more of the variation in gene expression (p-value < 2.2 × 10-16). In addition to the ADME genes, we also investigated the role of SNPs and CpG site methylation in the regulation of genes associated with diseases and liver function by querying all SNPs from the GWAS catalog (http://www.genome.gov/gwastudies/) in our list of identified eQTLs. We identified cis-acting SNPs and DNA methylation differences that were associated with the expression of 47 genes previously identified in different GWA studies with complex traits, including enzyme and metabolite levels as well as cardiovascular and inflammatory bowel diseases (Additional file 15).

15Since we had also generated methylation and expression data for three other tissues (muscle, subcutaneous- (SAT) and visceral adipose tissue (VAT)) from the same individuals in the Dutch cohort, we could assess the tissue-specificity of the detected eQTL, meQTL and eQTM effects. We had previously compared liver eQTL with other tissues for only a limited number of samples, so we re-did this analysis with the new adult liver samples from the Karolinska Liver Bank. For liver eQTL, approximately 40-50% of the effects found in one tissue could also be significantly detected in another tissue (Figure 4A). We identified only a few opposite allelic effects (<1%) between the tissues (Additional file 16A and B), suggesting that if a SNP affects expression in multiple tissues, the allelic direction is mostly identical. The eQTL effects (n = 32,863) that were only present in liver and not in the other three tissues were related to genes strongly specific to liver function (p-value 5 × 10-53) and metabolic and catabolic processes (p-values < 5 × 10-20).

16Contrary to the strong tissue-specificity of eQTL, meQTL were much more stable across the different tissues. On average, 70% of the meQTL are shared between at least two tissues, with over 98% of their effects having the same allelic direction (Figure 4B, Additional file 16C & D). As we had observed for eQTL, there were also a few significant meQTL that showed an opposite allelic direction between liver and the other three tissues (Additional file 17). The CpG sites of the meQTL with opposite effects were more often located outside the gene bodies (p-value 1.53 × 10-11), but when they were in gene bodies, they were in exons rather than introns (p-value 1.7 × 10-90).

17As expected, we observed very strong tissue-specificity for the identified eQTMs. Only up to 4% of the eQTMs found in one tissue were also detectable with an identical effect direction in another tissue (Figure 4C, Additional file 16E).

18We hypothesized that the expression of genes with important functions in the adult liver should be under strict genetic and epigenetic control. We thus focused on the set of probes with significantly higher and lower expression levels in adult liver compared to fetal liver (from the previous section “The transcriptome of the developing liver”), and formed a matched set of probes that were not differentially expressed between the two groups but displayed similar median expression levels and standard deviations in the adult liver samples. We observed that the expression of these adult liver-specific probes are much more likely to be affected by SNPs than the matched set of probes (1.43 times more than expected, chi-squared test p-value = 8.8 × 10-7). Furthermore, we observed that these probes were 1.24-fold enriched for liver-specific eQTL probes compared to a matched set of probes with eQTL in multiple tissues (chi-squared test p-value = 8.847 × 10-6). Vice versa, probes with lower expression levels in adult liver compared to fetal liver did not differ from the matched set of probes in terms of having eQTL and liver-specific eQTL effects. Furthermore, we did not observe any enrichment of meQTL in the adult liver-specific methylation probes.

19Previous studies on the regulation of gene expression in human liver have only accounted for the effect of genetic variation in adult samples. In this study, we investigated the developmental regulation of gene expression in human livers by comparing the expression and methylation levels of genes in adult and fetal livers. In addition, we used both genetic variants and DNA methylation differences in order to explain the variability in transcript levels observed in adult livers. Comparison of the fetal and adult liver methylomes and transcriptomes revealed that hypomethylated CpG sites and up-regulated genes were closely related to the tissue-specific functions: with fetal livers enriched for developmental and hematopoietic functions, while catabolic and metabolic processes were more prominent in adult livers. This has been described in the transcriptome of fetal livers at different stages of development in mice.

20As the differences in methylation between fetal and adult livers were very large, when attempting to characterize the effects of variable methylation on gene expression levels in adults, we performed the eQTM analysis using a panel of only adult liver samples. Similarly to Gutierrez-Arcelus et al., we observed both positive and negative correlations between DNA methylation and gene expression across the samples, with similar distributions across different genomic regions. Bell et al. have also observed a modest but significant excess of negative correlations between DNA methylation and variation in gene expression levels across individuals. It has been reported that the role of DNA methylation appears to depend on the genomic context: for example, CpG sites located near the genes and/or with a stronger correlation between the methylation and expression were more likely to display a negative correlation. Interestingly, CpG sites downstream of the expression probes displayed less negative correlations than those upstream of the probes, indicating that methylation in gene bodies is associated with active gene expression, as known from the early days of DNA methylation research. This paradox – in which methylation in the promoter is negatively correlated with the expression, whereas methylation in the gene body is positively correlated with expression – can be explained by the fact that, in mammals, DNA methylation silences the initiation of transcription, but not transcription elongation.

21Our eQTL mapping in adult livers revealed 751 unique genes, which were strongly liver-specific and enriched for drug metabolizing functions. Of these, 84 genes were new associations, while others have already been reported. The new associations are probably due to the larger number of samples and imputation of SNPs not present on previously used genotyping arrays, using data from the 1000 Genomes project. While we observed liver-specific associations with eQTL, the meQTL were not enriched for liver-specific functions. Furthermore, when we analyzed the SNPs that had significant effects on both methylation and expression, in most of the genes the same SNP allele had an opposite effect on gene expression compared to the methylation level, and this effect was most evident in CpG islands and shores (Additional file 10). These results show that, although there are many associations between SNPs and methylation levels, the relationships between them are not clear and do not reflect tissue-specific functions.

22Inter-individual variability in ADME gene expression has been shown to affect drug efficacy, toxicity, and susceptibility to environmental toxins. When we focused on the expression of ADME genes, we observed very strong cis-acting SNP and DNA methylation effects for 16 genes (Table 4), including members of the glutathione S-transferases (GSTs) family of phase II ADME isozymes: GSTA4, GSTM1, GSTM3, GSTO2 and GSTT1; solute carrier transporters SLC19A1 and SLC22A18, responsible for the transmembrane transfer of multiple drugs and endogenous compounds; and FMO4, GPX7, PON1 and UGT1A1. GSTT1 is involved in the conjugation of a variety of compounds, while GSTM1 functions in the detoxification of exogenous/endogenous toxins. The effects of epigenetic modifications on the expression of these genes have been reported in blood and brain tissues. In our study, we observed that both SNPs and DNA methylation contribute to the variability of the expression of these genes. For example, the SNP rs2739330, downstream of the GSTT1 gene and upstream of the DDT gene, has been reported to be associated with gamma-glutamyl transferase levels in plasma. This SNP, together with methylation levels of a nearby CpG site cg05380919, explains 78% of the variability in the expression of GSTT1, possibly with a stronger contribution from the methylation levels of the CpG site. Similarly, for GSTM1 the strongest SNP only explains 11% of the variation in its expression, while methylation levels of the CpG site cg18938907 has a much stronger association with the expression of the gene, and may be responsible for up to 55% of the variation (Additional file 11). The CpG site falls within a CpG island that spans the promoter and a portion of the gene’s first intron. Interestingly, they are located near the transcription factor binding site of TBP, which has been shown to bind to the promoter of GSTM1 in HepG2 cells, according to ENCODE ChIP-Seq data.

23A substantial portion of the overall phenotypic variance in hepatic enzyme PON1 activity between individuals remains unexplained. Besides a variety of non-genetic factors, numerous transcription factors and miRNA regulation, various functional PON1 polymorphisms have been shown to influence serum PON1 levels and activity. The SNP rs705379 has been shown to be associated with approximately 50% mean reductions in serum PON1 protein levels as well as transcript levels. In our study, it was interesting to see that this SNP was associated with increased methylation of nine CpG sites in its vicinity and with a lower expression of PON1.

24Glutathione peroxidases (GPX) constitute a major antioxidative damage enzyme family and are thus important in cancer therapy. Not only genetic but also epigenetic mechanisms of gene regulation have been proposed for GPX7, while recently a CpG island was identified as a key player in regulating GPX7 expression. In total, we identified 87 SNPs in the GPX7 gene affecting the methylation of nine CpG sites (meQTL), with six of the sites being directly implicated in quantitative differences in the gene expression. In addition, we replicated the two eQTL reported in human liver samples for GPX7, and for the first time identified their association with differences in methylation of CpG sites, which are further correlated with changes in GPX7 expression levels. One of the expression-associated SNPs discovered in this study, rs11810754, appears to explain most of the variation in the expression levels of the gene (48%), while the CpG site with the strongest correlation with expression (cg11953272) did not add any extra information to the variability of the expression of the gene (only 1%, Table 4).

25Three other tissues (muscle, SAT, VAT) were used to assess the tissue-specificity of both eQTL and meQTL effects. Our eQTL results showed that over half of the associations of SNPs with expression in one tissue could not be detected in another tissue (with identical eQTL allelic effect directions). This is similar to other studies. In contrast, SNP-methylation correlations were much less tissue-specific than SNP-expression correlation: approximately 70% of the meQTL were also identified in any of the other tissues. To the best of our knowledge, this has not been described before, but could be due to sequence-dependent DNA methylation or the fact that genetic variation in a similarly methylated region can affect the entire region (given that we have excluded the direct effects of SNPs on methylation probes). On the other hand, we observed that DNA methylation associated with expression levels (i.e. eQTMs) are highly tissue-specific, in accordance with the fact that DNA methylation plays an important role in regulating tissue-specific gene expression. Thus, conclusions drawn from eQTL or eQTM data in one tissue cannot be extrapolated to other tissues, whereas the effect of SNPs on methylation is more likely to be detectable in an alternative tissue, for example DNA in blood, which is more readily accessible.

26The greatest limitation of our study was the use of microarrays instead of massively parallel sequencing. Despite stringent filtering and remapping of expression and methylation probe sequences, we cannot rule out all technical artefacts inherent to microarray studies. Another drawback of microarrays is also their lower coverage of the genome, with the expression arrays only covering a few exons per gene, and the methylation array containing approximately 1% of the CpG sites in the human genome. Future studies using RNA and genome sequencing should be able to generate a more complete picture of the factors involved in the regulation of gene expression in human liver or other tissues.

27By performing a genome-wide survey of genomic and epigenomic variation and their associations with gene expression in fetal and adult human liver, we have generated a comprehensive resource for the analysis of factors involved in the regulation of hepatic gene expression. The investigation of fetal livers allowed us to explore the developmental changes in the hepatic methylome and transcriptome. Although the role of DNA methylation in different regions of the genome is still unclear, our results elucidate the coordinated effects of SNPs and methylation, as well as the tissue specificity of their effects on gene expression. This strengthens the hypothesis that knowledge of inter-individual variability, driven by genetic polymorphisms and DNA methylation marks and their interaction, is crucial for understanding the causes of differences in drug response and the etiology of diseases associated with liver function.

28The materials and methods of this study are described in detail in the Supplementary Online Methods in the Additional file 1. Briefly, our study was performed on two different cohorts, 14 fetal and 96 adult liver samples from the Karolinska Liver Bank cohort, and 85 adult samples from the Dutch tissue cohort MORE (BBMRI obesity cohort). For both datasets, the number of samples for which there is full expression, methylation and genotype data is not 100%. We therefore report the number of samples per specific analysis. DNA from the samples were genotyped using HumanOmni BeadChips (Illumina), according to the manufacturer’s instructions. We imputed both datasets using the GIANT release from the 1000 Genomes project, resulting in 5,763,069 unique SNPs, which were used in all our downstream analyses. Gene expression data was generated using HumanHT-12 BeadChips (Illumina), according to the standard protocol. Bisulfite-converted DNA samples were hybridized to Infinium HumanMethylation450 BeadChips (Illumina), following the Illumina Infinium HD Methylation protocol.

29The data sets supporting the results of this article are available in the GEO repository, Dutch BBMRI more expression data: GSE22070; methylation data: GSE61454; Karolinska expression and methylation data: GSE61279.