Super-silencers are crucial for development and carcinogenesis in B cells

Abstract

The strength of the repressive histone H3K27me3 signal varies across silencers. Focusing on regions with unusually strong signals-super-silencers-we show that B-cell super-silencers are initially linked to gene upregulation in development, with target genes highly expressed in stem cells. About 13% of B-cell super-silencers convert to super-enhancers in B-cell lymphoma; 22% of these recur in over half of patients. Genes like BCL6 and BACH2 tied to these conversions are downregulated faster by JQ1, a super-enhancer-disrupting anti-cancer agent. Super-silencers are enriched for B-cell cancer-associated variants-both somatic and germline-and translocation breakpoints, exceeding levels in other regulatory elements like CTCF binding sites. Over 80% of B-cell lymphoma t(3;14)(q27;q32) translocations fuse BCL6 super-silencers with enhancer-rich regions. Super-silencer repression depends on CpG content: CpG-rich elements block promoter-enhancer contacts; CpG-poor - inhibit looping. These findings highlight super-silencers' key role in B-cell regulation and suggest their alteration may be a primary factor of B-cell carcinogenesis.

Fig. 1. Identification of SSs in GM12878 cells.

a Distribution of H3K27me3 ChIP-seq intensity in the silencers (blue) and of H3K27ac in the enhancers (pink). Intensity levels were linearly normalized to have a maximum of 1 for display purposes. b Results of ROSE with the intensities of H3K27me3 in GM12878 silencers. c Genomic distribution and overlap with CGIs of different enhancer and silencer types. d Average methylation levels of CpG sites in different enhancer and silencer types. Expression of genes (e) most proximal to and (f) having contacts with different enhancer/silencer types. Tissue-specificity of genes (as measured by tau) (g) most proximal to and (h) having contacts with different enhancer/silencer types. In (e–h), “All” represents all genes profiled in RNA-seq data. “SS comp” and “SE comp” denote SS and SE components, respectively. P values ($**P<{10}^{-10}$) were determined against “All” (i.e., all genes assayed in RNA-seq experiments, see Supplementary Information) using two-sided Wilcoxon rank sum tests without adjustment (n = 14,388). The SS components (represented by dark blue) are frequently located close to or target the genes lowly expressed in GM12878 and with a high tissue specificity. i Average ATAC-STARR-seq scores of different enhancer and silencer types. P values ($**P<{10}^{-10}$) were determined against “TE” using two-sided Wilcoxon rank sum tests without adjustment (n = 645 for “SS comp”; 3848 for “TS”; 3726 for “SE comp” and 3848 for “TE”). j Activity similarity among SS components within the same SS regions across cell types. k TFBS similarity between SS components located within the same SS regions across cell types. In ( j, k), “SS comp” and “SE comp” represent SS and SE components, respectively. In the left panel of ( j, k), “BK” is the TSs randomly selected to match the distribution of distances between SS components, and “random” represents the randomly selected silencers. Similarly, “BK” and “random” in the ( j, k) right panel are SE-distribution-matching and randomly selected TEs, respectively. In ( j, k), $**P<{10}^{-10}$ were determined against “BK” by two-sided Wilcoxon rank sum tests without adjustment. In the left panels, n = 13,671 for “SS comp”, 15,028 for “BK”, and 29,170 for “random” in left panels, while n = 35,518 for “SE comp”, 9728 for “BK”, 24,562 for “random” in right panels. In (e–k), the center line (in red) in a box shows the median; the box bounds represent the lower and upper quartiles; the whiskers extend to the minima and maxima points up to a maximum of 1.5× the interquartile range; and dots denote outliers.

Fig. 2. Experimental validation of SSs.

a Elements tested in the experiments, accompanied by the genomic and epigenetic profiles of the regions hosting these elements in normal and cancer genomes. The numbers in the track of “DLBCL SE” and “CLL SE” are the numbers of patients having the corresponding SEs out of all patients. b Luciferase activities of the tested SSs in GM12878 and K562. Data are presented as the median ± SD (n  =  9 for controls, and n = 3 for SSctrl1 … SSctrl3 and SS1… SS4, i.e., the sequences located within the SS). All these samples are biological replicates. P values were determined against control elements by two-sided t tests without adjustment. Source data are provided in Supplementary Data 3.

Fig. 3. Functional analysis on GM12878 SSs.

a Biological processes associated with SS components and TSs. These are the evaluation results of GREAT in which two-sided binomial tests were performed. P values presented were the results without adjustment. b Expression of genes associated with different REs across the cell types relevant to B-cell differentiation and the Psoas muscle (used as reference). P values above the boxplot $\left(**P<{10}^{-10}\right)$ were determined against all genes assayed in RNA-seq experiments (n = 14,388, see Supplementary Information) using two-sided Wilcoxon rank sum tests without adjustment. c Top-ranked biological processes enriched by the genes in which the promoters interact with SS components. These are the evaluation results of DAVID in which Fisher’s exact tests were used. P values presented were the results without adjustment. Analysis results for TSs are presented in Supplementary Fig. 11. d Expression of genes linked with different REs by Hi-C contacts. P values $\left(**P<{10}^{-10}\right)$ were determined using two-sided Wilcoxon rank sum tests without adjustment (n = 14,388 of all genes assayed in RNA-seq experiments, see Supplementary Information). In (b, d), the center line (in red) in a box shows the median; the box bounds represent the lower and upper quartiles; the whiskers extend to the minima and maxima points up to a maximum of 1.5× the interquartile range; and dots denote outliers. Throughout the figure, “SS comp” and “SE comp” represent SS and SE components, respectively.

Fig. 4. Disease association of SSs.

Enrichment of GWAS SNPs associated with (a) all traits and (b) traits relevant to immunity, B-cell lymphoma, and CLL. “WG” is the whole human genome, and “BK” is the background sequences, i.e., randomly selected DNase-seq peaks in other cell types. “SS comp” and “SE comp” represent SS and SE components, respectively, throughout this figure. P values $\left(**P<{10}^{-10}\right)$ were determined against “WG” using two-sided binomial tests without adjustment. c Enrichment fold of GWAS SNPs associated with individual traits. P values $\left(**P<{10}^{-10}\right)$ were determined for significant enrichment of the traits having an enrichment fold of > 2 in a subgroup of traits as compared to all GWAS traits, using two-sided binomial tests without adjustment (n = 2722 for all, 137 for immunity-associated, 8 for B-cell-lymphoma-associated, and 3 for CLL-associated GWAS traits). The center line (in red) in a box shows the median; the box bounds represent the lower and upper quartiles; the whiskers extend to the minima and maxima points up to a maximum of 1.5× the interquartile range; and dots denote outliers. d Fraction of replicated GWAS SNPs in different enhancer/silencer types. P values $\left(**P<{10}^{-10}\right)$ were determined against “SS comp” (representing SS components) using one-sided binomial tests without adjustment.

Fig. 5. Enrichment of B-cell-cancer variants in GM12878 SSs.

a Fractions of B-cell-cancer SNVs in different RE types. Fractions of recurrent (b) B-cell-cancer SNVs and (c) non-B-cell cancer SNVs in different enhancer and silencer types. d Fractions of B-cell-cancer TLBPs in and by different enhancer and silencer types. “IN” and “BY” represent the RE sequences and the genomic regions proximal to enhancers or silencers within $\pm 50\text{kb}$. In (a–d), “SS comp” and “SE comp” represent SS and SE components, respectively, while “BK” is the background sequences, i.e., randomly selected DNase-seq peaks in other cell types. In (a, b, d), P values ($**P<{10}^{-10}$) were determined (against “BK”, if not specified) using two-sided binomial tests without adjustment (n = 81,565,198 cancer SNVs and n = 1916 B-cell-cancer TLBPs documented in the ICGC data portal). e Distribution of target enhancers or silencers of B-cell-cancer translocations (TLs). “BK” is randomly selected DNase-seq peaks in other cell types. f Genomic and epigenetic profiles in the BCL6 SS (left panel) and IGHM SE (right panel) regions in normal and cancer genomes. “SS comp” and “SE comp” represent SS and SE components, respectively. These two regions are merged by the DLBCL translocation t(3;14)(q27;q32). The distribution of non-B-cell-cancer mutations (i.e., SNVs and TLBPs) is not presented here since close-to-zero non-B-cell-cancer mutations are reported in these regions. Another DLBCL translocation t(3;14)(q27;p14) joins this BCL6 SS region with the RHOH SE region (Supplementary Fig. 10). The numbers in the track of “DLBCL SE” are the numbers of patients having an SE and all patients.

Fig. 6. Conversion of GM12878 SSs to SEs during carcinogenesis.

a Enrichment of SS components in the loci of DLBCL-essential and DLBCL suppressor genes. “SS comp” and “SE comp” represent SS and SE components, respectively, while “BK” is DNase-seq peaks randomly selected from other cell types, throughout this figure. $**P<{10}^{-10}$ were determined against “BK” using two-sided binomial tests without adjustment. b Enrichment of enhancers and silencers in each DLBCL-essential and DLBCL suppressor gene locus. c Fractions of GM12878 enhancers and silencers coinciding with DLBCL SE regions. $**P<{10}^{-10}$ were determined against “BK” using two-sided binomial tests without adjustment. d Numbers of DLBCL patients in which SEs were detected. SS-to-SE and TS-to-SE represent the DLBCL SEs converted from GM12878 SSs and TSs, respectively. P values ($**P<{10}^{-10}$) were determined (against all DLBCL SEs if not specified) using a one-sided binomial test without adjustment (n = 2892). e Average binding intensities of BRD4 (as measured by TF ChIP-seq signals) in DLBCL SEs and BK sequences in DLBCL Ly4 cells. These intensity levels were detected at 24 h after JQ1 and DMSO treatments. The P value was determined against Ly1-cell SEs carrying BRD4 ChIP-seq signals using a two-sided t-test without adjustment (n = 127). f Fractions of genes significantly upregulated and downregulated at different time points following JQ1 treatment. The genes are categorized based on their association with different DLBCL SE groups. P values ($**P<{10}^{-10}$) were determined against all tested genes using two-sided binomial tests without adjustment (n = 28,869). g Expression levels of BCL2 and BACH2 at different time points following JQ1 treatment in DLBCL cells. P values between DMSO and JQ1 were determined using two-sided t tests (n = 29 biological samples). h Numbers of CLL patients in which SEs were detected. SS-to-SE and TS-to-SE represent CLL SEs converted from GM12878 SSs and TSs, respectively. P values ($**P<{10}^{-10}$) were determined (against all CLL SEs if not specified) using a one-sided binomial test without adjustment (n = 1010). In (d, e, g, h), the center line of a box shows the median; the bounds of a box represent the lower and upper quartiles; the whiskers extend to the minima and maxima points up to a maximum of 1.5× the interquartile ranges.

Fig. 7. Enrichment of SSs in TAD-shores.

a Distribution of REs in different TAD sections (boundaries, shores, and centers). The pi chart presents the distribution of SSs in TADs. “SS comp” and “SE comp” represent SS and SE components, respectively, throughout this figure. b Enrichment of silencers and enhancers in SS-counterpart-TAD-shores. The schematic shows two neighboring TADs. The Shore2 is the counterpart-TAD-shore of the Shore1, and vice versa. Each row in this heatmap represent a TAD boundary. Enrichment fold is calculated through comparing with the whole genome. c Enrichment of enhancers in the SS-counterpart-TAD-shores. The values above bars are the significance level as compared to all TAD-shores. CTCF-TAD-shores have no significant preference for enhancers. Data are presented as the median±SEM. P values ($**P<{10}^{-10}$) were determined against all TAD boundaries in GM12878 (n = 6240 as documented in the Peakachu) using two-sided Wilcoxon rank sum tests without adjustment. d Enrichment of disease-associated mutations in different RE types across TAD sections. Pie charts are the distribution of mutations in SS components. “BK” is DNase-seq peaks randomly selected in other cell types. P values ($**P<{10}^{-10}$) were determined against “BK” using one-sided binomial tests without adjustment.

Fig. 8. Repression models of SSs.

a Enrichment of chromatin contacts to different enhancer and silencer types (as illustrated in bar plots) and contact enrichment between two enhancer and silencer types (in the bubble plots). Throughout the figure, “SS comp” and “SE comp” represent SS and SE components, respectively. P values ($**P<{10}^{-10}$) were determined (against “BK”, i.e., the background sequences, which were randomly selected DNase-seq peaks in other cell types, if not specified) using one-sided binomial tests without adjustment. “PR” represents all protein-coding promoters. b Two major repression models utilized by silencers. c Contact enrichment in the proximities of different enhancer and silencer types. Data are presented as the median±SEM. P values ($**P<{10}^{-10}$) above bars were determined against all chromatin contact anchors in GM12878 (as documented in the Peakachu) using two-sided t tests without adjustment (n = 22,817). d Contact enrichment of promoters, which are categorized based on their enhancer or silencer contacts. Data are presented as the median±SEM. P values ($**P<{10}^{-10}$) above bars were determined against all promoters having Hi-C contact in GM12878 (as reported in the Peakachu) using two-tailed t tests without adjustment (n = 6951). e Numbers of ChIP-seq TFBSs located within silencers and enhancers. Data are presented as the median±SEM. $**P<{10}^{-10}$ were determined by two-sided t tests without adjustment (n = 2,755,484 ChIP-seq peaks for 168 TFs or co-transcription factors, see Supplementary Information). f Signatures of ChIP-seq TFBSs and epigenetic marks in CGI and non-CGI silencers. All enhancers were used as the background in the top four rows. CGI enhancers were used as background for the bottom two rows. High similarity between the top two and bottom two rows suggests that CGI is one of the factors determining TFBS features of CGI silencers.