Integration of Hi-C and ChIP-seq data reveals distinct types of chromatin linkages

Abstract

We have analyzed publicly available K562 Hi-C data, which enable genome-wide unbiased capturing of chromatin interactions, using a Mixture Poisson Regression Model and a power-law decay background to define a highly specific set of interacting genomic regions. We integrated multiple ENCODE Consortium resources with the Hi-C data, using DNase-seq data and ChIP-seq data for 45 transcription factors and 9 histone modifications. We classified 12 different sets (clusters) of interacting loci that can be distinguished by their chromatin modifications and which can be categorized into two types of chromatin linkages. The different clusters of loci display very different relationships with transcription factor-binding sites. As expected, many of the transcription factors show binding patterns specific to clusters composed of interacting loci that encompass promoters or enhancers. However, cluster 9, which is distinguished by marks of open chromatin but not by active enhancer or promoter marks, was not bound by most transcription factors but was highly enriched for three transcription factors (GATA1, GATA2 and c-Jun) and three chromatin modifiers (BRG1, INI1 and SIRT6). To investigate the impact of chromatin organization on gene regulation, we performed ribonucleicacid-seq analyses before and after knockdown of GATA1 or GATA2. We found that knockdown of the GATA factors not only alters the expression of genes having a nearby bound GATA but also affects expression of genes in interacting loci. Our work, in combination with previous studies linking regulation by GATA factors with c-Jun and BRG1, provides genome-wide evidence that Hi-C data identify sets of biologically relevant interacting loci.

Figure 1.

Flow chart of data processing. The sequential analytical steps of this study rely on several types of experimental input. The process is as follows: (i) analyze Hi-C data using a MPRM and a power-law decay background; (ii) group interacting loci using hierarchical clustering based on the loci’s epigenetic status as determined by ChIP-seq analysis of modified histones and regions of open chromatin; (iii) apply the Apriori algorithm to identify transcription factor (TF) association networks using ChIP-seq analysis of transcription factors and (iv) validate the effect of TF-induced DNA loops on gene regulation through comparison of gene expression profiling using RNA-seq before and after knockdown of a single transcription factor.

Figure 2.

Hi-C analysis and genomic interactions in K562 cells. (A) Hi-C data analysis. Several steps were taken to select real interactions from the initial set of hybrid fragments. First, self-ligation was filtered based on its special properties (Supplementary Methods). Second, a MPRM was used to eliminate random loops. Next, the proximate ligation threshold was determined (Supplementary Figure 1C), and those interactions that pass the threshold were further filtered by a power-law decay background model. (B) Circos plot of the identified genome-wide inter-chromosomal interactions. The outer layer ring represents the human genome, with each chromosome labelled a different color. Each link between two genomic loci denotes one chromosomal interaction. (C) Circos plot of intra-chromosomal interactions for two representative chromosomes showing that thousands of interactions form within a relatively short distance on individual chromosomes; in fact, intra-chromosomal interactions (with the majority of the two loci being within 1 million bp) represent 95% of all interactions.

Figure 3.

Clustering interacting loci. (A) Interacting loci were divided into 12 groups using hierarchical clustering based on the epigenetic status of loci 1. In many cases, the epigenetic status of interacting loci 2 shows a similar pattern to that of the partner loci 1. (B) Comparison of the distribution of the average intensities for each mark in each cluster. The plots of broad peak marks, including H3K9me3, H3K27me3, H4K20me1 and H3K36me3, are centered at the position of interacting loci, and the plots of sharp peak marks, including open chromatin, H3K4me1/2/3, H3K9ac and H3K27ac, are centered at the middle point of each peak. The X-axis is the relative distance to the central position of the analyzed mark, and the Y-axis is the tag density of that mark at each position.

Figure 4.

Relative distance of the interacting loci to a transcription start site. The X-axis is the relative distance between the nearest transcription start site (TSS) and the interacting loci, and the Y-axis is the occurrence frequency (i.e. the frequency that an interacting loci and the nearest TSS have a distance of x) in each cluster; 20 kb up- and downstream of the nearest TSS was analyzed, using a bin size of 2 kb.

Figure 5.

Interacting loci are bound by a network of transcription factors. (A) Hierarchical clustering suggests the existence of at least three discriminating groups of factors, namely, insulator-binding factors such as CTCF and RAD21 (bound to regions of open chromatin), promoter-binding factors such as c-Myc and E2F4 (bound to regions marked by open chromatin, H3K4 methylation and H3K9/27 acetylation) and non-promoter-binding factors such as GATA1 and GATA2 (bound to regions marked by open chromatin and H3K4 methylation). Promoter-binding TFs are mainly present in cluster 7 and 8, whereas non-promoter-binding factors showed less frequency in cluster 7 and 8 but higher frequency in cluster 9. This indicates that different transcription factors preferentially bind to regions with distinct epigenetic status. The left color bar corresponds to the color scale for the epigenetic marks, and the right color bar corresponds to the color scale for the TF-binding sites. (B) TF interaction network revealed by the Apriori algorithm. The color of each node represents the binding preference of the TF showed in Figure S9.

Figure 6.

Gene expression analyses reveal two types of chromatin linkages. Gene expression heatmaps show the expression of genes in loci 1 and 2 versus a random set of genes in each cluster. This analysis revealed two types of chromatin linkages. Type I: genes associated with both interacting loci in each pair have a higher gene expression level than a random set of genes; the interacting genes are likely co-activated (permissive chromatin linkages). Type II: genes associated with both interacting loci in each pair have lower expression levels than a random set of genes (non-permissive chromatin linkages).

Figure 7.

Gene expression changes induced by GATA1 and GATA2 knockdown. RNA-seq was performed using control cells and cells treated with siRNAs to either GATA1 or GATA2, and genes whose expression was altered upon knockdown were identified. By intersection of the list of deregulated genes with genes having nearby bound GATA factors, direct target genes were identified (blue pie segments). A second set of direct target genes were identified by intersection of the list of deregulated genes with genes linked to a bound GATA factor via chromatin looping (maroon pie segments). All other genes were classified as indirectly regulated genes (green pie segments). The colors in the pie chart correspond to these different categories of GATA-regulated genes, as labelled in left panel. For each category, the percentage of genes that were upregulated versus downregulated in the knockdown cells is shown by the purple and orange closed pie graphs. (G: GATA factor; Pol: RNA Polymerase II, TF: a transcription factor whose expression is regulated by a GATA factor).

Figure 8.

Schematic model of GATA-regulated chromatin linkages. Dynamic chromatin interactions form globule structures, which may function to initiate or stabilize three-dimensional gene regulatory structures. Multiple CREs, including enhancers and promoters, are bridged by different groups of transcription factors (TFs) and mediators under different conditions. In the example shown, a pair of interacting loci from cluster 9 is represented. A loop is formed between a promoter and a distal region via interactions of GATA1, BRG1, INI1 and SIRT6 (a histone deacetylase), all bound to a region having H3K4me1 but no marks of an active enhancer or promoter; the nearby gene is repressed. On loss of GATA1 (via reduction of levels of the protein by treatment with siRNAs or on normal physiological changes or owing to displacement of GATA1 by another DNA-binding factor), a different set of enhancer-binding factors are recruited to the distal open chromatin region (which now gains the H3K27Ac mark), the loop changes from a repressive to an activating structure, active histone modifications are placed on the promoter region, RNA Polymerase II is recruited and begins transcription, and the transcribed region gains H3K36me3. We note that GATA1 can also activate transcription, and thus other functions of GATA1-mediated loops can be envisioned.