4See: A Flexible Browser to Explore 4C Data

Abstract

It is established that transcription of many metazoan genes is regulated by distal regulatory sequences beyond the promoter. Enhancers have been identified at up to megabase distances from their regulated genes, and/or proximal to or within the introns of unregulated genes. The unambiguous identification of the target genes of newly identified regulatory elements can thus be challenging. Well-studied enhancers have been found to come into direct physical proximity with regulated genes, presumably by the formation of chromatin loops. Chromosome conformation capture (3C) derivatives that assess the frequency of proximity between different genetic elements is thus a popular method for exploring gene regulation by distal regulatory elements. For studies of chromatin loops and promoter-enhancer communication, 4C (circular chromosome conformation capture) is one of the methods of choice, optimizing cost (required sequencing depth), throughput, and resolution. For ease of visual inspection of 4C data we present 4See, a versatile and user-friendly browser. 4See allows 4C profiles from the same bait to be flexibly plotted together, allowing biological replicates to either be compared, or pooled for comparisons between different cell types or experimental conditions. 4C profiles can be integrated with gene tracks, linear epigenomic profiles, and annotated regions of interest, such as called significant interactions, allowing rapid data exploration with limited computational resources or bioinformatics expertise.

Keywords: 4C; biological replicates; browser; chromatin loops; epigenomics; quantile normalization.

Figure 1

Overview of the 4C method and analysis. (A) The 4C method entails chromatin fixation, restriction digestion and re-ligation to generate hybrid sequences between fragments that were physically proximal during fixation. The DNA is purified, digested with a secondary restriction enzyme and re-ligated under dilute conditions to generate DNA circles. Chimeric products linked to a specific bait fragment of interest (orange) are amplified by inverse PCR with bait-specific primers (orange arrows) flanked by Illumina sequencing adapters (black overhangs). The PCR products are then directly loaded onto Illumina flow-cells for high-throughput sequencing. (B) Pre-processing steps before 4See; tools denoted in bold accompany this manuscript. The fastq sequences are first trimmed to remove bait restriction fragment sequence (orange), leaving just the prey DNA sequence (green) for mapping to the reference genome with Bowtie. The mapped genomic coordinates are converted to restriction fragment space by a custom perl script, which counts the total number of reads mapping to each fragment on the same chromosome as the bait. This “cis” file can then be directly input into 4See.

Figure 2

Overview of graphing options from existing 4C analysis methods. All methods have been applied to two ES replicates and one NPC 4C data set for the Sox2 SCR bait (see also Figure 3 ). (A) 4Cseqpipe (van de Werken et al., 2012) results shown independently for one ES replicate and the NPC data set. Running median scores are plotted as a line graph (5 kb resolution), with domainograms plotted underneath as a heat map for median scores at steadily increasing window sizes. Positions of the CTCF site within the SCR and the Sox2 gene are indicated by arrows. Note that the independent normalization means that the SCR-Sox2 interaction differences between the two cell types is not evident, compared to other methods. (B) r3Cseq (Thongjuea et al., 2013) results for the combined three data sets, showing panels, from top to bottom: positions of Refseq genes; restriction fragment coverage; averaged profile for the “experiment” (ES) condition, with called interactions at different confidence levels highlighted; profile for the “control” (NPC) condition, with called interactions at different confidence levels highlighted; plot of the log2-ratio of experiment vs control 4C signal. The position of the SCR bait is indicated by a red dashed line. Note that r3Cseq appears to call a very large number of interactions. (C) fourSig (Williams et al., 2014) results shown independently for one ES and one NPC replicate, showing smoothed 4C plots (21-fragment windows). Called interactions (“Categories” 1, 2 or 3, for different confidence levels) would be highlighted on the plot, but none were called by fourSig. (D) FourCSeq (Klein et al., 2015) results shown independently for one ES and one NPC replicate, showing normalized and processed 4C signal plots (gray line graphs and black points), alongside the positions of known genes. The green line indicates the centralized 4C value, and the dashed blue lines indicate the threshold values for a z-score difference > 2. Significant interactions would be highlighted, but were not detected by FourCSeq. Note that differences between the two cell types is not so evident, compared to most other methods. (E) peakC (Geeven et al., 2018) results shown for the combined analysis of the two ES replicates, independently of the NP replicate, giving smoothed 4C plots (21-fragment windows) as histograms. The red regions indicate called interactions. (F) 4C-ker (Raviram et al., 2016) results shown for the combined analysis, with line plots of combined ES (red) and NPC (blue) 4C signal. Note that many interactions were called by 4C-ker within the plotted window for both ES and NP, but that the documentation did not provide a means to plot them alongside the shown line plots.

Figure 3

4See provides flexibility in handling multiple replicates and/or experimental conditions. (A) The 4See dialog box for conditions settings automatically opens when cis files are first loaded, in this instance two ES replicates and one NPC 4C data set for the Sox2 SCR bait. The two ES replicates have been assigned different integers to be treated independently, and the NPC data set has been omitted by assigning it 0. (B) The resultant 4See plot from the conditions set in (A), whereby the two ES 4C replicates are quantile normalized and plotted separately, one in black, the other gray. The plot has normalized 4C signal as the y-axis and genomic coordinate of the interacting fragment as the x-axis. The position of the SCR bait is denoted by a black vertical line, and gene position (blue arrows) and the ES H3K27ac ChIP-seq profile (black) is shown underneath the 4C plot. The profiles are highly consistent between replicates, with a strong interaction peak centered on the Sox2 gene; note that both the gene and enhancer have a strong enrichment for H3K27ac. (C) As for (A), but in this case the two ES replicates are given the same value to be averaged together, and the NPC data set is included as a different integer to the ES data sets. (D) As for (B), but with the settings conditions of (C), and the redundant gene and H3K27ac tracks omitted. The averaged ES 4C plot is given in black and the NPC 4C plot in red, showing a strong perturbation of the SCR-Sox2 interaction on differentiation.

Figure 4

4See provides flexibility in running mean window sizes. (A) The main control panel of 4See, including options to set the x-axis (“start coordinate,” “end coordinate,” and “plot window”) and y-axis (“max y plot”) plot limits, to choose a bait name for the plot title (“bait name”), and to set the running mean window size (“smooth window”) in number of restriction fragments. (B) 4See plots for a 4C data set in DP thymocytes with the Satb1 gene promoter as bait. Two instances of the same x- and y-axis limits are shown, with a running mean window of 21 (left) or 55 (right) fragments. The position of the Satb1 bait is denoted by a black vertical line, and gene position (blue arrows) and the DP H3K27ac ChIP-seq profile (black) is shown underneath the 4C plots. Pink rectangles denote regions called as interacting by peakC for the equivalent window size as the plot. For the long-range interaction, the smaller window size appears to have more spurious called interactions, less evidently linked to H3K27ac peaks; the link is better seen with greater smoothing from a larger window size. (C) As for (B), but with the Ikzf1 gene promoter as bait. In this case, the smaller window size (17 fragments) seems to give better resolution of specific interactions with distinct putative enhancers, which are merged into one at larger window sizes.

Figure 5

4See provides flexibility in handling epigenomic track scales. (A) The same 4C profile as Figure 4B (55-fragment window size; black) is plotted alongside the 4C profile for ES cells, where the locus is silent, and thymocyte enhancer interactions are not evident. (B) The 4See dialog box for managing epigenomic tracks defines how the different tracks are scaled. In this case, the ChIP-seq tracks for H3K27ac in ES and DP cells are treated independently, so the autoscaling of the ES track creates some spurious peaks from noise above background. (C) As for (B), but this time the two H3K27ac tracks have been set the same integer, making the lower ES signal much more visually apparent in the plot.

Figure 6

4See exploration can uncover previously overlooked features of chromatin topology. (A) Reproduced Supplemental Fig 7A from Narendra et al. (2015), used with permission. 4C profiles from the Hoxa5 bait are shown as domainogram heat maps for wild-type ES cells (top), as well as lines that have had deletions of one CTCF site (middle; site of deletion denoted by red asterisk) or two (bottom; additional deletion site denoted by green asterisk). The CTCF ChIP-seq profile is shown above the 4C sets. No chromatin topology differences are apparent between these cell lines, and the original study concluded that spatial phenotypes only occurred on cell differentiation (Narendra et al., 2015). (B) The same data, processed and plotted using 4See. The position of the Hoxa5 bait is denoted by a black vertical line, and gene position (blue arrows) and the ES CTCF ChIP-seq profile (black) is shown underneath the 4C plots. Red and green asterisks denote the positions of the single and double CTCF site deletions, as for (A). In this plot, CTCF-dependent restriction of interactions between Hoxa5 and more caudal regions (e.g. the gene body of Hoxa10) seems apparent.