Skip main navigation

Epigenetic Analyses of Human Left Atrial Tissue Identifies Gene Networks Underlying Atrial Fibrillation

Originally publishedhttps://doi.org/10.1161/CIRCGEN.120.003085Circulation: Genomic and Precision Medicine. 2020;13:e003085

Abstract

Background:

Atrial fibrillation (AF) often arises from structural abnormalities in the left atria (LA). Annotation of the noncoding genome in human LA is limited, as are effects on gene expression and chromatin architecture. Many AF-associated genetic variants reside in noncoding regions; this knowledge gap impairs efforts to understand the molecular mechanisms of AF and cardiac conduction phenotypes.

Methods:

We generated a model of the LA noncoding genome by profiling 7 histone post-translational modifications (active: H3K4me3, H3K4me2, H3K4me1, H3K27ac, H3K36me3; repressive: H3K27me3, H3K9me3), CTCF binding, and gene expression in samples from 5 individuals without structural heart disease or AF. We used MACS2 to identify peak regions (P<0.01), applied a Markov model to classify regulatory elements, and annotated this model with matched gene expression data. We intersected chromatin states with expression quantitative trait locus, DNA methylation, and HiC chromatin interaction data from LA and left ventricle. Finally, we integrated genome-wide association data for AF and electrocardiographic traits to link disease-related variants to genes.

Results:

Our model identified 21 epigenetic states, encompassing regulatory motifs, such as promoters, enhancers, and repressed regions. Genes were regulated by proximal chromatin states; repressive states were associated with a significant reduction in gene expression (P<2×10−16). Chromatin states were differentially methylated, promoters were less methylated than repressed regions (P<2×10−16). We identified over 15 000 LA-specific enhancers, defined by homeobox family motifs, and annotated several cardiovascular disease susceptibility loci. Intersecting AF and PR genome-wide association studies loci with long-range chromatin conformation data identified a gene interaction network dominated by NKX2-5, TBX3, ZFHX3, and SYNPO2L.

Conclusions:

Profiling the noncoding genome provides new insights into the gene expression and chromatin regulation in human LA tissue. These findings enabled identification of a gene network underlying AF; our experimental and analytic approach can be extended to identify molecular mechanisms for other cardiac diseases and traits.

Introduction

Over the past 15 years, genome-wide association studies (GWAS) have identified thousands of genetic variants associated with diseases and traits. Despite this success, a major challenge has been the difficulty in linking disease-associated variants to specific genes and mechanisms. The vast majority of GWAS variants reside outside of the protein-coding region of the genome and are presumed to modulate the function of enhancers or repressors that, in turn, regulate the expression of a proximal gene or genes. Thus, understanding of the mechanism of GWAS variants will require careful annotation of the noncoding genome in a disease-relevant tissue.

In addition to a myriad studies by independent laboratories, large consortia efforts, such as ENCODE (Encyclopedia of DNA Elements)1 and the Roadmap Epigenomics project,2 have gone to great lengths to annotate the noncoding genome. Initially, these chromatin regulatory elements were identified by hand, generally by perturbing enhancers in mouse or drosophila, and examining developmental changes.3,4 Currently, these elements can be identified en masse using chromatin immunoprecipitation sequencing (ChIP-seq) and algorithms to identify colocalization of histone post-translational modifications or transcription factors. These annotations rely upon profiling DNA accessibility, as well as multiple histone post-translation modifications, after which one can identify defined chromatin regulatory motifs, such as enhancers, promoters, bivalent, and silenced regions.5

The ENCODE project focused on profiling cultured cell lines and used methods such as DNase Hypersensitive regions to profile open chromatin regions that were likely to be accessible to transcriptional machinery. The Roadmap Epigenomics project profiled primary human tissue and cell lines from over 100 distinct lineages and has produced chromatin state maps in 3 chambers of the heart, including the left ventricle (LV), right ventricle, and right atrium, but not in the left atrium (LA). The LA plays a particularly important role in diseases, such as atrial fibrillation (AF) and cardiac conduction disorders.6 As many features of chromatin topology are unique to a given cell or tissue type, it is important that the tissue studied be relevant to the disease or phenotype in question to implicate mechanisms of genetic risk.7 It was the goal of this study to fill this important knowledge gap.

To this end, in the present work, we have created a model of LA chromatin states and integrated this data with (1) available chromatin state models from the other chambers of the heart, (2) LA transcription as assessed by RNA sequencing, (3) DNA methylation, (4) long-range chromatin interaction data from primary human LA and LV derived from HiC,8 and (5) GWAS data for AF and electrocardiographic traits. Using this approach, we have defined a gene regulatory network underlying AF and PR interval susceptibility loci.

Methods

Data Access

The LA tissue ChIP-seq data sets and chromatin state models generated in this study will be made freely available upon publication of the article and can be found at the Broad Cardiovascular Disease Knowledge Portal (www.broadcvdi.org). The raw fastq and aligned bam files for the LA tissue ChIP-seq data sets will be available at dbGaP (https://www.ncbi.nlm.nih.gov/gap/), accession number phs001539.v1.p1 upon publication.

Human Tissue Samples

Adult human myocardial samples of European ancestry were collected from deceased organ donors by the Myocardial Applied Genetics Network (www.med.upenn.edu/magnet). For all donors, clinical examination and medical history displayed no indications of structural heart disease. Hearts were transported to the lab in ice-cold cardioplegia solution until cryopreservation (always <4 hours). Written informed consent for research use of donated tissue was obtained from next of kin in call cases. Research use of tissues was approved by the relevant institutional review boards at the Gift-of-Life Donor Program, the University of Pennsylvania, Massachusetts General Hospital and the Broad Institute.

Results

The overall goal of this study was to profile a series of histone post-translational modifications in human LA tissue and define regulatory motifs to better understand the underlying functionality of noncoding DNA regions. An overview of the workflow of this study is detailed in Figure 1.

Figure 1.

Figure 1. Identifying enhancers specific to the left atria (LA) using human LA tissue. Overview of the study design. We profiled histone modifications in human LA tissues and used a HMM (hidden Markov model) method5 to produce a model of chromatin states. This model of distinct chromatin states was then integrated with gene expression profiles, and DNA methylation to validate the functionality of chromatin states by expression and methylation levels. To better understand the specific biology of enhancing elements in the LA, we used enhancer lists from the Roadmap Epigenomics Project2 to remove any identified enhancers that overlap with those from left ventricle (LV), right atrium (RA), or right ventricle (RV) tissues. Finally on these LA-specific enhancers, we performed de novo motif identification using HOMER (Hypergeometric Optimization of Motif Enrichment),9 identified gene-enhancer connections using HiC performed in primary LA, analyzed these genes for enriched pathways, and prioritized effector genes at cardiovascular disease (CVD) genome-wide association studies (GWAS) loci.

Chromatin States Identified by the Model

We began by performing ChIP-seq using tissue samples from the lateral wall of the LA obtained from 5 individuals aged from 51 to 62 years old with no evidence of cardiac pathology. The baseline characteristics of the donors are displayed in the Table. To identify regulatory regions of diverse functionality in LA tissue samples, we modified existing chromatin immunoprecipitation ChIP protocols,10,11 and profiled 5 active histone modifications (H3K4me1, H3K4me2, H3K4me3, H3K27ac, and H3K36me3) and 2 repressive histone modifications (H3K9me3 and H3K27me3), as well as the multifunctional insulator binding protein CTCF (CCCTC-binding factor). Overall, a median of 107 million paired-end reads was generated per sequenced IP experiment, with an average alignment rate of 97.5%. The average properly paired rate was high, at 96.4%, and sequenced libraries were complex, with an average duplication rate of only 3.2%. After calling significant peaks using MACS2 (model-based analysis of ChIP-seq), the sum total of hg38 genome coverage from all experiments was 5.99% (this is at the target level, please see Methods for details). Similarly, RNA-sequencing was of high quality, the average number of paired-end reads per library was 33.6 million, with an average alignment rate of 88.22%. To model these data, we used the fold enrichment over background for all ChIP data to identify peak regions for each experiment. These peak regions were ingested into ChromHMM,5 and we identified 21 distinct chromatin states, each a unique combination of histone post-translational modifications. We annotate the putative functionality of these chromatin states and the relative frequency of each mark at each state in Figure 2.

Table. Clinical Characteristics of the LA Tissue Donors

SampleSexEthnicityAge, yWeight, kgHeight, cmHW, gLV mass, gLVEF, %AFVT/VFDiabetesHTN
P01202MWhite5689172431NA65NoNoNoNo
P01221FWhite5259156300NA75NoNoNoNo
P01294MWhite5310619154023363NoNoNoNo
P01600FWhite5168162.621313450NoNoNoYes
P01623MWhite6211518853231850NoNoNoYes

AF indicates atrial fibrillation; HTN, hypertension; HW, heart weight; LA, left atria; LVEF, left ventricular ejection fraction; NA, not available; VF, ventricular fibrillation; and VT, ventricular tachycardia.

Figure 2.

Figure 2. A 21-state model of chromatin in left atrial (LA) tissue —chromatin states and characteristics. The above table describes the output of the 21-state model of chromatin in human LA developed using ChromHMM5 for this study. The heatmap from column 3 (H3K27me3) to column 10 (H3K27ac) defines the percentage of the time that the described chromatin mark (in the column header) is present for a given chromatin state (rows). The broad functional categories in column 11 (functional category) are further described in column 12 (putative function). Numbers are provided as percentages. TS indicates tissue specific.

As expected, transcribed genes displayed active modifications, and transcriptionally silent genes displayed repressive modifications at their promoters (Figure 3). An exemplar region is illustrated in Figure 3A that distinguishes active chromatin over the promoter for CALM1 (calmodulin),, a calcium sequestration and signaling gene highly expressed in the heart. In contrast, there is silent chromatin over the promoter for NEUROG2, a gene with expression confined to neuronal tissue (Figure 3B). When correlating peak data from all experiments using the Pearson R, experiments tended to cluster by histone mark, rather than by sample identity (Figure 3C), illustrating that the biological state of the chromatin was the primary driver of experiment clustering. When directly clustering all peak experiments in the heatmap illustrated in Figure 3D, this same trend was evident. The peaks from CTCF, in particular, illustrate that high Z-scores across this set of experiments tend to co-occur; this is also visually prominent for H3K4me3 experiments. The cooccupancy of these marks reflects biology, with active chromatin marks, such as H3K4me3, H3K4me2, and H3K27ac colocalizing, whereas H3K9me3 occupancy is distinct, given that repressive and active marks very rarely occupy the same space on the genome in the same cell type.

Figure 3.

Figure 3. Epigenetic profiles in left atria (LA) samples.A, Active chromatin (see box) over the promoter for CALM1, which is highly expressed in LA tissue. Arrow indicates direction of transcription. B, Silent chromatin over the promoter for NEUROG2, which is expressed only in neural tissue, and not in heart tissue. Arrow indicates direction of transcription. C, A clustered correlation heatmap of all chromatin profiles generated in this study. The pairwise correlation coefficients across the 18 309 genomic loci showing measurable signal in at least 5 experiments are shown, with clusters of chromatin marks indicated by the text. D, Heatmap of normalized ChIP-seq signals in the same 18 309 genomic loci shown in C. Both genomic loci (vertical) and tissue experiments (horizontal) were hierarchically clustered.

Chromatin States Affiliated With Transcription and Methylation

A view of the ChIP-seq reads surrounding promoter, enhancer, bivalent, and polycomb silenced chromatin states illustrated that the identified chromatin states faithfully recapitulate the underlying combinations of histone marks (Figure 4A). We affiliated chromatin states in each sample with matched, normalized RNA-seq data. Here, the presence of a repressive mark was associated with a significant decrease in expression when compared with all other chromatin states (t test, P<2×10−16). Reciprocally, promoter states were also associated with a significant increase in median expression, when compared with any repressed state (t test, P values range from 7.7×10−6 to <2×10−16, promoters are chromatin states 5 to 8, and 12 to 15, Figure 4B, Figure I in the Data Supplement).

Figure 4.

Figure 4. Association of chromatin state with transcriptional output and methylation status in the left atrium (LA).A, Epigenetic signal in a 20 Kb region surrounding 4 different chromatin states in the sample LA1623. The data were sorted by H3K4me3 signal for promoters, H3K4me1 signal for enhancers, and H3K27me3 signal for bivalent regions and polycomb-repressed regions. B, Average expression of genes closest to each of the 4 states displayed in A. Normalized log2 transformed TPM counts across all tissues were derived from DEseq2.12C, Average methylation across each of the states in A. Fractional methylation levels were derived from whole-genome bisulfite sequencing (WGBS) data in heart (average of left ventricle [LV], right ventricle [RV], right atria [RA]) from the Roadmap Epigenomics project2 and intersected with our data using bedtools intersect.13

We used available whole-genome bisulfite sequencing from LV, right ventricle, and right atrium tissue2 to approximate methylation levels corresponding to each chromatin state. The fractional methylation data plotted in Figure 4C and Figure II in the Data Supplement represents an average over all 3 primary tissues. The canonical promoter state (state 8, H3K4me3, and H3K27ac) had very low levels of DNA methylation, while canonical enhancers (state 10) were highly variable in their levels of methylation, bivalent regions (state 4) tended to be highly methylated, and polycomb regions (state 1) less methylated. Although genes proximal to enhancer regions are more highly expressed on average than genes proximal to only promoter regions (without an enhancer region present), the enhancing regions themselves are more highly methylated (t test, P<2×10−16). A possible explanation for this is that whole-genome bisulfite sequencing is unable to distinguish between methylation and 5-hydroxymethylcytosine; methylation of cytosines is generally considered to be silencing, while 5-hydroxymethylcytosine specifically binds enhancers and gene bodies in embryonic stem cells.14 Consistent with previous reports, CTCF binding sites were associated with high levels of methylation (Figure II in the Data Supplement, state 16).15

LA Enhancers: Strength, Pathway Enrichment, and Comparisons With Other Data Sets

The active enhancer state in our model is defined by H3K27ac and H3K4me1 colocalization4 and corresponds to state 10 of the model. However, there are other states which may have enhancing activity as well, such as colocalizations of H3K27ac, H3K4me2, and CTCF (state 9, 11, 15, and 19), and our analysis of enhancers includes all of these enhancing states. In total, we identified 67 183 individual enhancing elements across all samples. Since some of these enhancers overlapped or were bookended across multiple samples, merging all enhancer regions resulted in 37 782 distinct enhancing regions. When comparing the expression of genes proximal to enhancer regions (the enhancer is wholly contained within the gene body or proximal to but outside of the gene body) with genes proximal to promoter regions without an enhancer near the promoter or gene body, enhancing regions drive gene expression higher (P<2×10−16) than promoter regions (Figure 4B). Similarly, to ENCODE and the Roadmap Epigenomics Project, enhancers had a wide range of widths, from 200 bp (the minimum bin size in ChromHMM) to over 30 Kb (Figure IIIA in the Data Supplement).

To understand how enhancer strength contributes to level of expression of a connected gene, we scored enhancers by affiliating them with fold enrichment scores derived from the primary peak calls in MACS2 (ChromHMM does not consider peak strength in chromatin state calls, see Introduction and Methods in the Data Supplement). We scored enhancer states (states 9, 10, 11, 15, and 19, Figure 2) for each sample using CTCF, H3K27ac, and H3K4me2 peak score data. After scoring and merging bookended enhancers as above, we were able to score 31 948 peaks across all samples (Figure IIIB in the Data Supplement). Enhancer scores were moderately correlated with enhancer width (Pearson R: 0.53, Figure IIIC in the Data Supplement). To directly understand the effects of enhancer width and score on gene expression, we associated enhancers with HiC loops derived from LA and LV tissue and identified 14 210 total enhancer-gene interactions. Using a normalized average expression value in LA for each gene linked to an enhancer, we were able to affiliate enhancers with a connected gene and gene expression value.

Enhancer width did not have a strong effect on the expression of the looped gene (Pearson R: 0.05, Figure IIID in the Data Supplement). To identify pathways enriched in genes linked to the strongest LA enhancers, we identified LA and LV concordant genes from the top 10% of identified enhancers and performed pathway analysis on this list of 222 protein-coding genes using the GO Biological Process database via PantherDB16 (Figure IIIE in the Data Supplement).

Over 15 000 Enhancers Are Unique to the LA

Having defined the total set of 37 782 LA enhancers, we compared these with other publicly available enhancer data sets. After comparing all LA enhancers with the enhancers identified in Roadmap2 enhancers, we found that 15 545 were unique to LA and not present in the other 3 chambers of the heart. Nearly 2000 of these enhancers intersected with the cardiac enhancers identified by Dickel et al17 (Figure 5A). We then sought to define the characteristics of these enhancers through motif identification, enhancer-gene connections using HiC, and annotation of cardiovascular disease GWAS loci.

Figure 5.

Figure 5. Identification and characterization of left atria (LA) specific enhancer states.A, Number of overlapping and unique LA-specific enhancers in the present study when compared with the Roadmap Epigenomics project and cardiac enhancers from Dickel et al.17B, Motifs identified from the top 10% of LA-specific enhancers, using HOMER.9 Repressors were all polycomb regions (state 1 in the 21-state model described in Table 2). C, Enriched GO Biological Process (BP) pathways from the top 25% of LA enhancers and connected genes using HiC. D, The percentage of cardiovascular disease (CVD) loci which overlap an LA enhancer (shared with other adult tissues from the Roadmap data set), or an LA-specific enhancer, identified only in our data set. The percentage of loci annotated by an LA-specific enhancer for each trait is as follows: atrial fibrillation (AF): 12.7, blood pressure (BP): 10.7, coronary artery disease (CAD): 11.5, PR: 11.4, QT: 20. E, Interconnectivity between PR interval and AF susceptibility loci (including proxy SNPs with an R2>0.8), and heart HiC loops. The genes are connected to the susceptibility loci by any loop present in the HiC heart data set. The size of the nodes is scaled by the number of connected edges. The weight of the edges reflects the combined score metric reported from StringDb,18 and this information is available in full in Table II in the Data Supplement. A view of this network containing only nodes at the StringDb defined medium (>0.4) and high (>0.7) confidence combined scores is shown in Figure V in the Data Supplement. All genes identified are expressed in the LA, although SCN10A, GJC2, and GPR144 have very few reads across the gene body. Node coloring is a proxy for gene function, per the legend at the bottom. % Enh: percentage of the time the motif is found in enhancers; and % Repr: percentage of the time the motif is identified in polycomb-repressed regions.

LA-Specific Enhancers Are Enriched for Homeobox Motifs

We used HOMER9 to perform motif analysis de novo of the top 10% (1554) strongest LA-specific enhancers, sorted by score. To perform background correction, we used the total set of all polycomb-repressed regions (state 1, please see Table 2) as these elements are likely to be both functional yet not involved in enhancing activity. Motif analysis of these enhancers identified several homeobox TFs (transcription factors) involved in heart biology, including SMAD3, NKX2-3, and the PITX (paired like homeodomain) family of motifs (Figure 5B, Table I in the Data Supplement).

We then examined enhancer-gene connections using heart HiC data from the top 25% of scored LA-specific enhancers and used Enrichr19 to identify GO Biological Process terms. Interestingly, we observed signatures of cardiac conduction biology, IL (interleukin)-15 signaling, and actin crosslinking (Figure 5C). Several enriched pathways reflect specific metabolic processes, such as succinyl-CoA and citrate metabolism. LA-specific enhancers are connected to genes involved with cardiac muscle cell membrane potential, as well as the cellular response to IL-15, and the triglyceride biosynthetic process pathway.

Chromatin States Proximal to Cardiovascular and Electrocardiographic Trait Loci

We next sought to understand the contribution of LA enhancers to genetic susceptibility for cardiovascular diseases. We annotated GWAS loci from AF, blood pressure, coronary artery disease, PR interval, and QT interval with LA enhancers and highlighted which enhancers are specific to the LA (Figure 5D). For these cardiovascular diseases and traits, the inclusion of LA-specific enhancers increased the number of loci which reside within an enhancing region identified in the LA. This allows us to assign a functionality to the locus, by demonstrating that at least one small nucleotide polymorphism (SNP) in this region is located within an enhancer identified in the LA. Using the comprehensive 21-state model most SNPs were located in quiescent regions (state 20; AF: 84.4; blood pressure: 88.8; coronary artery disease: 85.6; PR interval: 88.1; and QT interval: 88.8). When considering those SNPs located within actively regulated regions (Figure IVA in the Data Supplement), the most common region of intersection is state 9, which consists of solely H3K27ac (an enhancing state), followed by the cell type-specific promoter state (state 7) and strong promoters (state 8). States 1, 2, 3, and 4 are silenced and make up a small minority of the SNP/state combinations, intersecting with each trait <1% of the time. Overall, when SNPs are located within a ChromHMM annotated region in LA, they are located in regions conferring active chromatin status, particularly one of the 5 identified enhancer states.

To identify a LA-specific gene interaction network, we intersected AF and PR SNP loci plus proxies with heart HiC loop data.20 Connections within 85 genes were determined using the StringDB18 database, and 81 genes form the interaction network illustrated in Figure 5E. This network is highly interconnected, with 10 genes having at least 10 connections to neighboring genes. The interconnectivity of the network is greater than expected and highly significant (P=1.19×10−07; hypergeometric test). The genes NKX2-5, TBX3, SCN5A, OBSCN, and XPO1 are prominent connection hubs (each having at least 14 connections to first-degree neighbors). Similar to the overall and the LA-specific enhancer data sets above, these interconnected genes are also significantly enriched (false discovery rate<0.05) for many gene ontology terms involved in cardiovascular function, such as cardiac ventricle development (GO:0003231), muscle contraction (GO:0006936), regulation of myotube differentiation (GO:0010830), and His-Purkinje system cell differentiation (GO:0060932). The full output from StringDb18 is provided in Table II in the Data Supplement, whereas Figure V in the Data Supplement illustrates the AF and PR network containing only nodes connected by edges at the StringDb defined cutoffs of medium (>0.4) and high (>0.7) confidence scores.

Discussion

While the nucleotide sequence of the genome has been extensively studied, and many GWAS variants have been identified, most of these disease-related variants map to noncoding regions of the genome. Mechanistically, these variants are assumed to act by modulating enhancer, repressor, or transcription factor binding activity, which in turn alters the expression of topologically proximal transcripts. Understanding how these variants function will require annotation of the noncoding genome, and this optimally should be in the most disease-relevant tissue, as chromatin topology is often unique for a given cell type.7

Integration of multiple data sets can improve identification of disease-associated candidate genes near (or connected to) a GWAS locus. The advent of chromatic conformation capture-based technologies such as 4C-seq allowed for interrogation of the topology of the genome in the nucleus, for example, long-range physical contacts originating from any single SNP locus.21 Concurrently, genotyping or sequencing of many tissue samples coupled with RNA-seq identified expression quantitative trait loci, which elicit changes in gene transcription.22 More recent approaches such as HiC have allowed interrogation of physical long-range contacts of DNA in an all by all matrix of the genome, while approaches such as STARR-seq assess the functional/transcriptional characteristics of specific susceptibility loci.23 Prioritization of effector genes near GWAS loci is a difficult problem to solve, but integrating chromatin functionality, gene expression, and chromatin conformation data allows for ranking of genes, near a GWAS locus when there are multiple genes, or can confirm a long-range SNP/gene interaction.

Here, we have generated chromatin state maps using 7 histone post-translational modifications, plus the insulator binding protein CTCF in 5 primary human LA tissue samples. Our model of chromatin states in the LA consists of 21 distinct states, which recapitulate known chromatin regulatory motifs like enhancers, promoters, and silenced regions. Unique to this model, we have profiled the tissue-type specific histone modification H3K4me2,24 which allows some identification of regulatory regions specific to the LA. Coupled with matched RNA-sequencing data, we were able to assess the amount of transcription originating from each state. We used these data, together with available data from primary heart tissue on DNA methylation, chromatin topology, and GWAS data for AF and electrocardiographic traits to build a model of the gene regulatory network underlying AF, involving enhancers as well as tissue-type specific promoters. Importantly, LA tissue is heterogeneous, all ChIP signals are an average of the cell types present in each tissue sample; in the LA, around 80% of cells are cardiomyocytes.25

Outside of the ENCODE1 and Roadmap Epigenomics2 projects, there is a relative dearth of epigenomic data in primary heart tissue—the vast majority of studies that consider the epigenome use either mouse tissue or cultured cardiomyocytes.17,26–32 These data can be quite helpful but have inherent limitations. For example, the biology of the cardiac conduction system in mouse and human is obviously different, with a 10-fold difference in the average heart rate. Another issue is the changes that occur in cultured cell lines versus primary tissue. Primary tissue provides a direct measurement of the state of a given (heterogeneous) cell population.33 Moreover, efforts to profile a single histone post-translational modification (such as H3K27ac) to identify regulatory elements in the heart are unable to reconcile the combinatorial nature of histone modification—active histone modifications can occur alone, but in combination may define many types of regulatory elements, not distinguishable using a single mark.

Here, we have used chromatin state models in human LA to define a model of chromatin states and annotated that model with external data that is also derived from human LA tissue. Our goal was to identify enhancers and tissue-specific promoters and determine how these elements can inform an understanding of cardiovascular biology with respect to disease. Using AF as an example, we made 2 primary observations.

First, we identified 15 545 enhancers specific to the LA chamber of the heart. These enhancers are enriched for motifs involved in SMAD/TGF-β (transforming growth factor) signaling, zinc fingers, and homeobox genes, including the PITX family, as well as the NKX2 (NK2 Homeobox) family. Using HiC from the LA and LV to connect LA-specific enhancers to genes illustrates enrichment in pathways for cardiac muscle cell membrane potential, and the cellular response to IL-15, which is implicated in cardiomyocyte survival under stress conditions.34 These identified pathways from healthy tissue recapitulate pathways which are differentially expressed in AF, such as TGF-β signaling and extracellular matrix formation.35 These LA-specific enhancers also annotate a substantial percentage of cardiovascular disease susceptibility loci, implying a gain in understanding of these loci using LA-specific enhancers.

Second, we used the loci identified in genetic studies of AF and the electrocardiographic PR interval to select a subset of arrhythmia specific genes connected via HiC loops. We identified a highly interconnected gene interaction network that underlies conduction of electrical signals through the LA, and thus AF. This interaction network identified many known genes involved in the pathology of AF and illustrates the potential interaction partners of these genes. Familiar members, such as NKX2-5,36SYNPO2L,37 and SCN5A,38 are intimately involved in the processes of heart development, AF risk, and cardiac conduction, respectively. We also identify several genes previously identified as differentially expressed in AF, such as GPR22, COG5, and MYH7.39,40 Using HiC to identify gene linked to susceptibility loci allows for expansion of the network into novel genes, as well as illustrating the prevalence of sarcomere, and genome regulatory genes in the governance of the cardiac conduction system. This network will be of assistance when building cellular and animal models of AF, to understand what nontarget genes may be affected in a mutant or knockout model.

Our study was subject to a number of potential limitations. First, epigenetic profiles were performed in bulk LA tissue and thus represent an admixture of cell types. Although this approach prevents cell-specific analyses, it was in keeping with the approach used by Roadmap Epigenomics Consortium and enabled comparisons to this data set. It will be interesting to compare the bulk epigenetic maps to tissue-specific epigenetic analyses as these methods develop. Second, there is heterogeneity in the clinical and environmental status of the samples used for our LA studies. Given our limited sample size, we were unable to resolve any specific epigenetic or transcriptomic differences across important potential epigenetic modifiers, such as age or hypertension. In ongoing work, we are increasing the number of healthy controls and future epigenetic studies in diseases such as AF cases versus controls will allow us to evaluate interindividual variation. Despite these potential confounders, it is reassuring that our results were qualitatively similar to data from the Roadmap Epigenetics Consortium. Finally, the GWAS associated genes that we have identified for AF and other traits will require future validation in cell or animal model systems.

In conclusion, we have developed a high-resolution epigenetic map of the primary adult human LA that enabled functional annotation of variants in electrocardiographic traits and cardiovascular disease and discovery of novel enhancer-gene and enhancer-promoter associations which are representative of pathways in cardiovascular health and disease. The LA-specific enhancers we identified substantially improved the annotation for 5 cardiovascular disease loci.

Nonstandard Abbreviations and Acronyms

AF

atrial fibrillation

CTCF

CCCTC-binding factor

ChIP-seq

chromatin immunoprecipitation sequencing

GWAS

genome-wide association studies

LA

left atrium

LV

left ventricle

SNP

small nucleotide polymorphism

Footnotes

*Drs Tucker and Ellinor are joint senior authors.

The Data Supplement is available at https://www.ahajournals.org/doi/suppl/10.1161/CIRCGEN.120.003085.

For Sources of Funding and Disclosures, see page 597.

Correspondence to: Patrick T. Ellinor, MD, PhD, Broad Institute of MIT and Harvard, 75 Ames St, Cambridge, MA 02124. Email

References

  • 1. ENCODE Project Consortium. An integrated encyclopedia of DNA elements in the human genome.Nature. 2012; 489:57–74. doi: 10.1038/nature11247CrossrefMedlineGoogle Scholar
  • 2. Kundaje A, Meuleman W, Ernst J, Bilenky M, Yen A, Heravi-Moussavi A, Kheradpour P, Zhang Z, Wang J, Ziller MJ, et al. Integrative analysis of 111 reference human epigenomes.Nature. 2015; 518:317–330. doi: 10.1038/nature14248CrossrefMedlineGoogle Scholar
  • 3. Hoffmann AD, Yang XH, Burnicka-Turek O, Bosman JD, Ren X, Steimle JD, Vokes SA, McMahon AP, Kalinichenko VV, Moskowitz IP. Foxf genes integrate tbx5 and hedgehog pathways in the second heart field for cardiac septation.PLoS Genet. 2014; 10:e1004604. doi: 10.1371/journal.pgen.1004604CrossrefMedlineGoogle Scholar
  • 4. Rada-Iglesias A, Bajpai R, Swigut T, Brugmann SA, Flynn RA, Wysocka J. A unique chromatin signature uncovers early developmental enhancers in humans.Nature. 2011; 470:279–283. doi: 10.1038/nature09692CrossrefMedlineGoogle Scholar
  • 5. Ernst J, Kellis M. ChromHMM: automating chromatin-state discovery and characterization.Nat Methods. 2012; 9:215–216. doi: 10.1038/nmeth.1906CrossrefMedlineGoogle Scholar
  • 6. Nattel S. New ideas about atrial fibrillation 50 years on.Nature. 2002; 415:219–226. doi: 10.1038/415219aCrossrefMedlineGoogle Scholar
  • 7. de Laat W, Duboule D. Topology of mammalian developmental enhancers and their regulatory landscapes.Nature. 2013; 502:499–506. doi: 10.1038/nature12753CrossrefMedlineGoogle Scholar
  • 8. van Berkum NL, Lieberman-Aiden E, Williams L, Imakaev M, Gnirke A, Mirny LA, Dekker J, Lander ES. Hi-C: a method to study the three-dimensional architecture of genomes.J Vis Exp. 2010; 39:1869. https://pubmed.ncbi.nlm.nih.gov/20461051/. Accessed December 8, 2020. doi: 10.3791/1869Google Scholar
  • 9. Heinz S, Benner C, Spann N, Bertolino E, Lin YC, Laslo P, Cheng JX, Murre C, Singh H, Glass CK. Simple combinations of lineage-determining transcription factors prime cis-regulatory elements required for macrophage and B cell identities.Mol Cell. 2010; 38:576–589. doi: 10.1016/j.molcel.2010.05.004CrossrefMedlineGoogle Scholar
  • 10. Landt SG, Marinov GK, Kundaje A, Kheradpour P, Pauli F, Batzoglou S, Bernstein BE, Bickel P, Brown JB, Cayting P, et al. ChIP-seq guidelines and practices of the ENCODE and modENCODE consortia.Genome Res. 2012; 22:1813–1831. doi: 10.1101/gr.136184.111CrossrefMedlineGoogle Scholar
  • 11. Hall AW, Battenhouse AM, Shivram H, Morris AR, Cowperthwaite MC, Shpak M, Iyer VR. Bivalent chromatin domains in glioblastoma reveal a subtype-specific signature of glioma stem cells.Cancer Res. 2018; 78:2463–2474. doi: 10.1158/0008-5472.CAN-17-1724CrossrefMedlineGoogle Scholar
  • 12. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2.Genome Biol. 2014; 15:550. doi: 10.1186/s13059-014-0550-8CrossrefMedlineGoogle Scholar
  • 13. Quinlan AR, Hall IM. BEDTools: a flexible suite of utilities for comparing genomic features.Bioinformatics. 2010; 26:841–842. doi: 10.1093/bioinformatics/btq033CrossrefMedlineGoogle Scholar
  • 14. Stroud H, Feng S, Morey Kinney S, Pradhan S, Jacobsen SE. 5-Hydroxymethylcytosine is associated with enhancers and gene bodies in human embryonic stem cells.Genome Biol. 2011; 12:R54. doi: 10.1186/gb-2011-12-6-r54CrossrefMedlineGoogle Scholar
  • 15. Wang H, Maurano MT, Qu H, Varley KE, Gertz J, Pauli F, Lee K, Canfield T, Weaver M, Sandstrom R, et al. Widespread plasticity in CTCF occupancy linked to DNA methylation.Genome Res. 2012; 22:1680–1688. doi: 10.1101/gr.136101.111CrossrefMedlineGoogle Scholar
  • 16. Mi H, Muruganujan A, Casagrande JT, Thomas PD. Large-scale gene function analysis with the PANTHER classification system.Nat Protoc. 2013; 8:1551–1566. doi: 10.1038/nprot.2013.092CrossrefMedlineGoogle Scholar
  • 17. Dickel DE, Barozzi I, Zhu Y, Fukuda-Yuzawa Y, Osterwalder M, Mannion BJ, May D, Spurrell CH, Plajzer-Frick I, Pickle CS, et al. Genome-wide compendium and functional assessment of in vivo heart enhancers.Nat Commun. 2016; 7:12923. doi: 10.1038/ncomms12923CrossrefMedlineGoogle Scholar
  • 18. Szklarczyk D, Gable AL, Lyon D, Junge A, Wyder S, Huerta-Cepas J, Simonovic M, Doncheva NT, Morris JH, Bork P, et al. STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets.Nucleic Acids Res. 2019; 47(D1):D607–D613. doi: 10.1093/nar/gky1131CrossrefMedlineGoogle Scholar
  • 19. Kuleshov MV, Jones MR, Rouillard AD, Fernandez NF, Duan Q, Wang Z, Koplev S, Jenkins SL, Jagodnik KM, Lachmann A, et al. Enrichr: a comprehensive gene set enrichment analysis web server 2016 update.Nucleic Acids Res. 2016; 44(W1):W90–W97. doi: 10.1093/nar/gkw377CrossrefMedlineGoogle Scholar
  • 20. Bianchi V, Geeven G, Tucker N, Hilvering CRE, Hall AW, Roselli C, Hill M, Martin JF, Margulies KB, Ellinor PT, et al. Detailed regulatory interaction map of the human heart facilitates gene discovery for cardiovascular disease.SSRN Electron J. 2019; 18:146. doi: 10.1186/s13059-017-1279-yGoogle Scholar
  • 21. Allahyar A, Vermeulen C, Bouwman BAM, Krijger PHL, Verstegen MJAM, Geeven G, van Kranenburg M, Pieterse M, Straver R, Haarhuis JHI, et al. Enhancer hubs and loop collisions identified from single-allele topologies.Nat Genet. 2018; 50:1151–1160. doi: 10.1038/s41588-018-0161-5CrossrefMedlineGoogle Scholar
  • 22. Aguet F, Brown AA, Castel SE, Davis JR, He Y, Jo B, Mohammadi P, Park YS, Parsana P, Segrè AV, et al. Genetic effects on gene expression across human tissues.Nature. 2017; 550:204–213. doi: 10.1038/nature24277CrossrefMedlineGoogle Scholar
  • 23. van Ouwerkerk AF, Hall AW, Kadow ZA, Lazarevic S, Reyat JS, Tucker NR, Nadadur RD, Bosada FM, Bianchi V, Ellinor PT, et al. Epigenetic and transcriptional networks underlying atrial fibrillation.Circ Res. 2020; 127:34–50. doi: 10.1161/CIRCRESAHA.120.316574LinkGoogle Scholar
  • 24. Pekowska A, Benoukraf T, Ferrier P, Spicuglia S. A unique H3K4me2 profile marks tissue-specific gene regulation.Genome Res. 2010; 20:1493–1502. doi: 10.1101/gr.109389.110CrossrefMedlineGoogle Scholar
  • 25. Tucker NR, Chaffin M, Bedi KC, Papangeli I, Akkad AD, Arduini A, Hayat S, Eraslan G, Muus C, Bhattacharyya RP, et al; Human Cell Atlas Lung Biological Network; Human Cell Atlas Lung Biological Network Consortium Members. Myocyte-specific upregulation of ACE2 in cardiovascular disease: implications for SARS-CoV-2-mediated myocarditis.Circulation. 2020; 142:708–710. doi: 10.1161/CIRCULATIONAHA.120.047911LinkGoogle Scholar
  • 26. Boogerd CJ, Aneas I, Sakabe N, Dirschinger RJ, Cheng QJ, Zhou B, Chen J, Nobrega MA, Evans SM. Probing chromatin landscape reveals roles of endocardial TBX20 in septation.J Clin Invest. 2016; 126:3023–3035. doi: 10.1172/JCI85350CrossrefMedlineGoogle Scholar
  • 27. He A, Gu F, Hu Y, Ma Q, Ye LY, Akiyama JA, Visel A, Pennacchio LA, Pu WT. Dynamic GATA4 enhancers shape the chromatin landscape central to heart development and disease.Nat Commun. 2014; 5:4907. doi: 10.1038/ncomms5907CrossrefMedlineGoogle Scholar
  • 28. Papait R, Serio S, Pagiatakis C, Rusconi F, Carullo P, Mazzola M, Salvarani N, Miragoli M, Condorelli G. Histone methyltransferase G9a is required for cardiomyocyte homeostasis and hypertrophy.Circulation. 2017; 136:1233–1246. doi: 10.1161/CIRCULATIONAHA.117.028561LinkGoogle Scholar
  • 29. Gomez-Velazquez M, Badia-Careaga C, Lechuga-Vieco AV, Nieto-Arellano R, Tena JJ, Rollan I, Alvarez A, Torroja C, Caceres EF, Roy AR, et al. CTCF counter-regulates cardiomyocyte development and maturation programs in the embryonic heart.PLoS Genet. 2017; 13:e1006985. doi: 10.1371/journal.pgen.1006985CrossrefMedlineGoogle Scholar
  • 30. Preissl S, Weichenhan D, Gilsbach R, Hein L, Grüning BA, Gelb BD, Jacobs AR, Backofen R, Krane M, Braun C, et al. Distinct epigenetic programs regulate cardiac myocyte development and disease in the human heart in vivo.Nat Commun. 2018; 9:391. doi: 10.1038/s41467-017-02762-zCrossrefMedlineGoogle Scholar
  • 31. Ang SY, Uebersohn A, Spencer CI, Huang Y, Lee JE, Ge K, Bruneau BG. KMT2D regulates specific programs in heart development via histone H3 lysine 4 di-methylation.Development. 2016; 143:810–821. doi: 10.1242/dev.132688CrossrefMedlineGoogle Scholar
  • 32. Goldman JA, Kuzu G, Lee N, Karasik J, Gemberling M, Foglia MJ, Karra R, Dickson AL, Sun F, Tolstorukov MY, et al. Resolving heart regeneration by replacement histone profiling.Dev Cell. 2017; 40:392–404.e5. doi: 10.1016/j.devcel.2017.01.013CrossrefMedlineGoogle Scholar
  • 33. Schnabel M, Marlovits S, Eckhoff G, Fichtel I, Gotzen L, Vécsei V, Schlegel J. Dedifferentiation-associated changes in morphology and gene expression in primary human articular chondrocytes in cell culture.Osteoarthritis Cartilage. 2002; 10:62–70. doi: 10.1053/joca.2001.0482CrossrefMedlineGoogle Scholar
  • 34. Yeghiazarians Y, Honbo N, Imhof I, Woods B, Aguilera V, Ye J, Boyle AJ, Karliner JS. IL-15: a novel prosurvival signaling pathway in cardiomyocytes.J Cardiovasc Pharmacol. 2014; 63:406–411. doi: 10.1097/FJC.0000000000000061CrossrefMedlineGoogle Scholar
  • 35. Thomas AM, Cabrera CP, Finlay M, Lall K, Nobles M, Schilling RJ, Wood K, Mein CA, Barnes MR, Munroe PB, et al. Differentially expressed genes for atrial fibrillation identified by RNA sequencing from paired human left and right atrial appendages.Physiol Genomics. 2019; 51:323–332. doi: 10.1152/physiolgenomics.00012.2019CrossrefMedlineGoogle Scholar
  • 36. Zhang L, Nomura-Kitabayashi A, Sultana N, Cai W, Cai X, Moon AM, Cai CL. Mesodermal Nkx2.5 is necessary and sufficient for early second heart field development.Dev Biol. 2014; 390:68–79. doi: 10.1016/j.ydbio.2014.02.023CrossrefMedlineGoogle Scholar
  • 37. van Eldik W, den Adel B, Monshouwer-Kloots J, Salvatori D, Maas S, van der Made I, Creemers EE, Frank D, Frey N, Boontje N, et al. Z-disc protein CHAPb induces cardiomyopathy and contractile dysfunction in the postnatal heart.PLoS One. 2017; 12:e0189139. doi: 10.1371/journal.pone.0189139CrossrefMedlineGoogle Scholar
  • 38. Arnolds DE, Liu F, Fahrenbach JP, Kim GH, Schillinger KJ, Smemo S, McNally EM, Nobrega MA, Patel VV, Moskowitz IP. TBX5 drives Scn5a expression to regulate cardiac conduction system function.J Clin Invest. 2012; 122:2509–2518. doi: 10.1172/JCI62617CrossrefMedlineGoogle Scholar
  • 39. Nothjunge S, Nührenberg TG, Grüning BA, Doppler SA, Preissl S, Schwaderer M, Rommel C, Krane M, Hein L, Gilsbach R. DNA methylation signatures follow preformed chromatin compartments in cardiac myocytes.Nat Commun. 2017; 8:1667. doi: 10.1038/s41467-017-01724-9CrossrefMedlineGoogle Scholar
  • 40. Deshmukh A, Barnard J, Sun H, Newton D, Castel L, Pettersson G, Johnston D, Roselli E, Gillinov AM, McCurry K, et al. Left atrial transcriptional changes associated with atrial fibrillation susceptibility and persistence.Circ Arrhythm Electrophysiol. 2015; 8:32–41. doi: 10.1161/CIRCEP.114.001632LinkGoogle Scholar

eLetters(0)

eLetters should relate to an article recently published in the journal and are not a forum for providing unpublished data. Comments are reviewed for appropriate use of tone and language. Comments are not peer-reviewed. Acceptable comments are posted to the journal website only. Comments are not published in an issue and are not indexed in PubMed. Comments should be no longer than 500 words and will only be posted online. References are limited to 10. Authors of the article cited in the comment will be invited to reply, as appropriate.

Comments and feedback on AHA/ASA Scientific Statements and Guidelines should be directed to the AHA/ASA Manuscript Oversight Committee via its Correspondence page.