Divergence of transcriptional landscape occurs early in B cell activation
Epigenetics & Chromatin volume 8, Article number: 20 (2015)
Signaling via B cell receptor (BCR) and Toll-like receptors (TLRs) results in activation of B cells with distinct physiological outcomes, but transcriptional regulatory mechanisms that drive activation and distinguish these pathways remain unknown.
Two hours after ligand exposure RNA-seq, ChIP-seq and computational methods reveal that BCR- or TLR-mediated activation of primary resting B cells proceeds via a large set of shared and a smaller subset of distinct signal-selective transcriptional responses. BCR stimulation resulted in increased global recruitment of RNA Pol II to promoters that appear to transit slowly to downstream regions. Conversely, lipopolysaccharide (LPS) stimulation involved an enhanced RNA Pol II transition from initiating to elongating mode accompanied by greater H3K4me3 activation markings compared to BCR stimulation. These rapidly diverging transcriptomic landscapes also show distinct repressing (H3K27me3) histone signatures, mutually exclusive transcription factor binding in promoters, and unique miRNA profiles.
Upon examination of genome-wide transcription and regulatory elements, we conclude that the B cell commitment to different activation states occurs much earlier than previously thought and involves a multi-faceted receptor-specific transcriptional landscape.
B cell activation, the transition from a naïve to an effector state, is important due to its essential role in immunity. Deregulated activation can have disastrous effects resulting in immune disorders and several B cell malignancies, some of which resemble activated B cell phenotypes [1, 2]. Mature resting splenic B cells maintain a quiescent G0 state with limited proliferative output . Upon encountering antigen, these cells become activated, leading to plasma cell differentiation and participation in immune responses. Activation of B cells can occur through either surface B cell receptor (BCR) [4–6] or various pathogen-associated molecular patterns such as bacterial lipopolysaccharide (LPS), which is mediated by Toll-like receptor (TLR) signaling and NF-κB pathways . Exposure of B cells to LPS via TLR4 can promote plasma cell differentiation [7, 8]. A properly regulated LPS activation appears critical as patients with deficient TLR signaling molecules, exhibit autoimmunity .
Splenic B cell differentiation can begin as early as 4 h and fully develop by 48–72 h . While much is known about signaling cascades during B cell activation at early and late time points [5, 11, 12], transcriptional changes during these times are still being addressed . In particular, a high-resolution picture reflecting the immediate transcriptional and epigenetic changes during early B cell activation, before mature B cells proceed toward proliferation and functional immune responses occur, is not available. Regardless of cell type, initial signaling events lead to rapid induction of primary response genes (PRGs) whose products initiate secondary waves of transcription resulting in egress from the G0 state and subsequently in proliferation and effector function . Regulatory mechanisms for these rapid responses, release of preformed “paused” transcription complexes, RNA polymerase II (Pol II) recruitment via transcription co-factors, and promoter accessibility/repression through histone modifications, are particularly well established [14, 15]. But how these mechanisms operate during BCR and LPS signaling in early activation of B cells is not known.
Activation of resting splenic B cells ex vivo provides a tractable model to explore this transition in a ligand-specific manner . B cells are unique hematopoietic cells because they express both BCR and TLRs. Although stimulation of naïve ex vivo B cells through both receptors elicits activation and proliferation, only LPS stimulation results in plasma cell differentiation [4, 17, 18]. We used this ex vivo model to determine how and when these ligand-specific transcriptional landscapes may diverge. We observe clear differences within 2 h post stimulation. In addition to ligand-selective differences in both protein-coding and non-coding RNAs, several other transcriptional regulatory steps differed between the activation states providing three key findings . Although BCR-induced genes show new recruitment of RNA Pol II that appeared to be paused at promoters, LPS/TLR4-induced genes exhibit enhanced transition of RNA Pol II from initiation to elongation . While the H3K4me3 (activating) mark is increased in both activation states (more so during TLR4 engagement), the appearance of the H3K27me3 (repressive) mark is reduced on BCR-responsive genes but remains relatively unchanged in LPS-responsive genes . Predicted transcription factor binding sites in the promoter proximal regions of genes also differ significantly in a ligand-selective manner. Together, our results show that B cell egress from the resting state involves a large pool of shared/common RNAs, and a small set of signal-selective RNAs that exhibit remarkable transcriptional landscape changes soon after ligand engagement.
Response dependent differential transcription during early activation
To determine how naïve B cells proceed along activation pathways in response to different signals, high-resolution RNA-seq analysis was performed after 30 and 120 min of stimulation with anti-IgM (engages BCR) or LPS (engages TLR4). Many expected genes, including the early activation marker CD69 and the PRGs c-fos, c-jun, and c-myc, were increased compared to the resting state by exposure to either anti-IgM/BCR or LPS/TLR4 as early as 30 min (Additional file 1: Figure S1). However, 120 min was required to identify genes/transcripts of the BCR and TLR4/NF-κB signaling pathways with robust P values (Additional file 1: Figure S1). Therefore, all subsequent experiments were performed at 120 min. The majority of the differentially expressed genes were shared between BCR and LPS activation henceforth referred to as shared response genes (Fig. 1a), which expectedly include c-myc, egr2, and irf4, along with more immune specific genes, ccl3, and cxcr4 (Fig. 1a). Smaller subsets of response selective transcripts were induced in each of the activated states with corresponding relevance to their specific stimuli (Fig. 1b). For example, BCR preferentially increased il-7r, which is associated with B cell development and differentiation , while LPS increased the inflammatory cytokines cxcl10 and ccl5 . Surprisingly, very few  genes were decreased in the LPS response. The significance of this observation remains unclear but this paucity of data precluded any meaningful comparison of a ligand-specific decrease in transcription. Therefore, subsequent analyses focused solely on induced mRNAs.
Figure 1c shows functional relationships between the different groups of genes with shared and response selective transcription that changed significantly relative to the resting state. c-Myc, which plays an important role in B cell proliferation in response to TLR4 and BCR ligation, was dominantly positioned in the shared response genes and also prevalent with BCR-selective genes (known Myc targets are highlighted in yellow at bottom left). As expected, preferentially increased LPS specific transcripts contained TLR signaling pathway genes and known NF-κB targets (Additional file 2). We conclude that the initial stages of B cell activation involve a large set of shared genes, despite being stimulated by distinct signaling pathways and a small but significant set of ligand-selective genes. These ligand-selective genes induced by BCR showed a predominant Myc signature, while those induced by TLR4 exhibited a prominent NF-κB signature.
We also analyzed 1315 genes whose transcription was not significantly altered by BCR or TLR4 stimulation. These genes, such as Polr2a, Max, ActB, and Dicer1, collectively annotated to biological processes, including maintenance of homeostasis, cell cycle, and apoptosis management (Additional file 1: Figure S2). In addition, transcripts associated with antigen processing, TGF-β signaling, TNF signaling, and MHC1 antigen presentation, were also unchanged (Additional file 1: Figure S2), suggesting that prior to activation, these cells are actively executing significant immune functions.
Response dependent RNA Pol II occupancy during early activation
To examine the transcriptional status of these genes, we investigated RNA Pol II occupancy at the promoters of both shared and preferentially induced genes. While RNA Pol II occupancy increased in genes, exhibiting both induced and unchanged transcription, a significant genome-wide increase in RNA Pol II recruitment to TSS(s) during BCR activation was also observed (Fig. 2a and Additional file 1: Figure S5). Statistical analysis (Additional file 1: Figure S4) indicated this difference was significant and not due to high variability. RNA Pol II occupancy did not increase during LPS activation at induced and unchanged genes, though the small increase in downstream RNA Pol II at genes preferentially induced by LPS suggested promoter-associated (paused) RNA Pol II was transiting to regions downstream of TSS upon LPS stimulation (Fig. 2a). These events occurred whether considering occupancy around the TSS(s) of all possible transcripts (Fig. 2) or at the TSS with the highest RNA Pol II occupancy (Additional file 1: Figure S3). To better determine the difference in downstream RNA Pol II occupancy between the two signals and whether these differences reflect transition from initiation to elongation, we calculated RNA Pol II traveling ratios (TR) [21, 22]. Both downstream analysis (Fig. 2b) and TR (Fig. 2c) clearly show that RNA Pol II occupancy downstream of TSS is significantly different between the two signaling pathways. However, because this assay measures the whole RNA Pol II population, it is not possible to conclusively determine whether RNA Pol II containing transcriptional complexes transition at the same or different rates. Nevertheless, these data collectively suggest that while genes induced by BCR recruited additional RNA Pol II, the transition from initiation to elongation is less, suggesting a pausing mechanism. In contrast, although the RNA Pol II association at the promoter was lower with LPS-induced genes, LPS activation results in a more robust transition from initiation to elongation mode.
Response dependent changes in chromatin
Chromatin modification is a well-defined transcriptional regulatory mechanism [14, 23]. Increased tri-methylation of promoter proximal histone 3 lysine 4 (H3K4me3) is associated with enhanced transcription complex recruitment and retention . A measurable level of H3K4me3 was observed in the resting state, which was further increased upon signaling through either BCR or TLR4. However, compared to BCR activation, LPS signaling via TLR4 increased H3K4me3 at genes with unchanged transcription (Fig. 3a, b) and in a genome-wide manner (Additional file 1: Figure S5). Increase in H3K4me3 was further enhanced at genes preferentially induced by LPS. While many of the genes regulating deposition of H3K4me3 are known , none were differentially transcribed in either BCR or TLR4 activation. We conclude that while increased H3K4me3 is associated with B cell activation in general, it is more pronounced at LPS/TLR4 responsive genes than at BCR-responsive genes.
Polycomb repressor complexes (PRCs) catalyze repressive chromatin marks [25–28]. In mammals, PRC2 activity depends on the SET domain-containing protein Enhancer of zeste homolog 2 (EZH2) that catalyzes transcriptionally repressive histone H3 methylation at lysine 27 mark (H3K27me3) . H3K27me3 plays an important role in B cell proliferation, and EZH2 expression is low in resting B cells but upregulated in activated B cells [30–32]. Hence, we analyzed H3K27me3 in response to BCR and TLR4 engagement. In genes with increased or unchanged transcription, the level of H3K27me3 at promoters was decreased during BCR activation (Fig. 3c). This decrease was most prominent around the TSS, and while there was a small decrease during LPS activation, this was restricted to an area just upstream of the TSS; otherwise, H3K27me3 in the LPS activation state remained largely unperturbed. BCR activation showed a greater genome-wide decrease in H3K27me3 than LPS (Additional file 1: Figure S5). While the boxplot comparison of median analysis did show some variation, the differences remained statistically significant. We concluded that there was a decrease in H3K27me3 during BCR engagement relative to TLR4 engagement around TSSs. Collectively, analysis of H3K4me3 and H3K27me3 illustrates that distinct chromatin modifications separate BCR and LPS activation states, regardless of whether considering all possible transcripts (Fig. 3) or the ones with the RNA Pol II highest occupancy (Additional file 1: Figure S6) was considered.
Analysis of individual genes
Having observed ligand-selective regulatory signatures at shared and response selective genes, we examined two representative genes from each category (shared, BCR- and LPS-selective genes) based on known functions in B cell activation and immunity and on the presence or absence of a Myc signature (for the shared and BCR-selective groups) or a NF-κB signature (for the TLR4/LPS-selective group). These genes collectively illustrated the characteristics observed at the global level. For example, the shared genes (Irf4 and Myc) as well as BCR-responsive genes (Il7r and Egr3) showed an increase in RNA Pol II and a decrease in H3K27me3 (Fig. 4). While Irf4 and Myc showed increased RNA expression in response to both stimuli, Il7r and Egr3 showed preferential RNA expression in response to BCR stimulation. Conversely, the LPS-responsive genes (Il12b and Tnf) showed a preferential increase in RNA expression and an increase in RNA Pol II at downstream sites, while the H3K27me3 mark around the TSS remained largely unchanged. Lastly, although H3K4me3 was increased by both stimuli, it was clearly more pronounced at Il12b and Tnf after LPS stimulation. Hence, the global transcriptional patterns observed with each category of genes are also observed at individual target genes.
Response selective transcription factor binding motif enrichment
We next examined promoter sequences to gain further insights into the observed differences in RNA Pol II recruitment or regulatory histone marks. Promoters containing “CpG islands” correlate with low nucleosome occupancy and increased RNA Pol II occupancy [33–35]. A majority of transcripts was associated with promoters containing or near (within 200 bps) a predicted CpG island; however, no preference for CpG associated promoters was evident for either response (Additional file 1: Figure S8). Although a connection between H3K27me3 and CpG islands exists , we observed no significant ligand-dependent correlation between CpG island promoters and H3K27me3 (data not shown).
To further shed light on the signal-specific activation states, transcription factor (TF) motifs within the promoter region (−1000 to +1000 of the TSS) were examined. TF motifs at induced gene promoters were clearly different between the two activation states. Most TF motifs were highly overrepresented in both the Shared- and BCR-induced gene sets, as shown by the large number of data points along the diagonal of the scatterplot (Fig. 5a). A much smaller degree of commonality was observed between the Shared- and LPS-responsive genes where many of the diagonal data points seen in the BCR vs Shared comparison migrated closer to the axes. Taking general motifs from each data set with the 30 lowest P value scores (highest association), we observed that many high-frequency motifs were common between Shared and BCR, while neither overlapped with the 30 lowest P value motifs found in LPS preferentially induced gene promoters (Fig. 5b). Direct comparison of all significantly represented TF motifs associated with the BCR- and LPS-responsive genes showed very little overlap (Fig. 5c). The general binding properties of the TFs associated with these motifs were strikingly different, while BCR activation is associated with a preponderance of helix-loop-helix containing TFs and STAT binding sites; NF-κB, GATA, and IRF binding sites dominated the LPS-selective signature (Fig. 5d). The direct role of these TFs in recruiting RNA Pol II or epigenetic marks is unknown but that the two responses exhibit such a mutually exclusive set of TFs motifs, which is noteworthy.
Expression of non-coding RNAs
Because miRNAs are known to play a pivotal role in various B cell processes such as development, inflammation, and tolerance [37, 38], we tested whether signal-selective miRNA expression could further account for the distinct transcriptional signatures. Analysis of a separate RNA sequencing set designed for improved miRNA detection revealed miRNA expression relative to the resting state that differs greatly between the two modes of activation (Fig. 6a, b). Consistent with previous observations, we found that both BCR and LPS enhanced miR-155 and miR-19b-1 expression [39, 40]. Expression of several miRNAs was decreased during BCR and LPS stimulation, including miR-125b, which when overexpressed, is capable of driving lymphoma development (Fig. 6b) . Of the 60 miRNAs with significantly changed expression relative to the resting state, the majority show a difference in expression; 20 increased by LPS and 29 decreased by BCR; the full list of miRNAs can be found in Additional file 1: Figure S9 and Additional file 3: Table S1. Two expression group representatives annotated to biological processes relevant to B cell activation, including several miRNAs decreased during BCR activation and increased by LPS, are shown (Fig. 6c). Overall, we found that a far greater decrease in miRNA expression occurred during the BCR response while LPS response showed a greater induction of miRNA compared to the resting state. An interesting time-dependent transcription of two RNA editing associated long non-coding RNAs (lncRNAs), Neat1, and Malat1 [42, 43], also suggested that lncRNA management of RNA editing might provide another level of regulation during early B cell activation (Additional file 1: Figure S10).
Deregulation of B cell activation can result in autoimmune disorders, chronic inflammation, and lymphoma. Although B cells express both BCR and TLR4, the functional consequences to these stimuli are distinct, at least under ex vivo conditions. While BCR engagement (triggered by anti-IgM stimulation) leads to proliferative expansion of activated B cells, engagement of TLR4 (triggered by LPS stimulation) leads to proliferation and eventual production of plasma cells. However, the transcriptional signatures and molecular mechanisms that distinguish these responses are relatively unknown. Here, we employed ex vivo activation of resting murine splenic B cells to examine these molecular signatures and define the transcriptional and regulatory landscape during early activation by high-resolution RNA- and ChIP-seq. We observe that at 2 h post stimulation, most genes (~90 %) induced by the two pathways are shared, despite being triggered through distinct receptors. This observation is consistent with previous studies [44, 45]. However, the rest of the transcriptome (10 %) exhibit signal-selective transcriptional programs.
While genes induced by both responses (shared) are greatly dominated by Myc targets and Myc-associated processes, groups of genes preferentially induced by LPS, as expected, show a strong correlation with NF-κB (Fig. 1). Surprisingly, BCR-induced genes showed little relationship to processes associated with the shared and LPS specific induced genes; they only overlapped with the shared-induced genes at the level of Myc and Arnt targets. Myc expression in normal lymphoid tissue is present in both activated and resting B cells in all phases of the cell cycle . We were surprised to find that the Myc co-factor Max transcription levels remained constant during activation despite a 15- to 30-fold increase in Myc, suggesting other co-factors were required to drive such a strong Myc-dominated signature. However, another Myc regulator, Mnt , was increased by BCR activation to a small degree. These data further underscore Myc’s dominant but enigmatic role in B cell activation [48–50].
Given the distinct gene sets noted in each of the responses, we examined mechanisms that could drive these observed differences. While preexisting, paused, RNA Pol II is associated with a large number of genes in diverse cell types, new RNA Pol II recruitment indicates transcriptional activity . Although preexisting RNA Pol II was found globally in resting B cells, RNA Pol II occupancy at the TSS was broadly increased during the BCR response and increased further at BCR-responsive gene promoters. This scenario suggests two general BCR-dependent mechanisms are in play , a global recruitment of Pol II to promoters and  a signal-specific guidance to intensify this general recruitment of Pol II to BCR-responsive promoters, particularly around the TSS. Conversely, LPS activation appeared to involve maintenance of steady-state RNA Pol II occupancy at the promoter relative to the resting cells. An increase in downstream RNA Pol II, possibly reflecting elongating Pol II, was noticeable, although this promoter associated RNA Pol II peak can be due to other mechanisms . To better understand this phenomenon, we took a closer look at downstream regions and calculated traveling ratios of polymerase occupancy in promoter versus downstream/coding sequences (Fig. 2c). These analyses revealed that while the BCR signal resulted in enhanced RNA Pol II at promoters, the transition to downstream region (elongation) was less/slower. In contrast, although LPS stimulation did not result in additional RNA Pol II recruitment, transition to downstream regions was greater/faster than that observed with BCR. It remains possible that although the total recruitment of RNA Pol II under two stimulations is very similar, the difference in promoter versus downstream region associated RNA Pol II reflects the rate at which the enzyme transits from initiation to elongation mode. Our preliminary analysis indicates that there are differences in the complexity of gene structure (e.g., number of exons) between the BCR- and LPS-responsive genes. Whether the difference in RNA Pol II movement between the two stimuli reflects these differences, a difference in signal strength and/or a fundamental difference in signaling pathways remain to be determined.
Epigenetic regulation of gene expression plays crucial roles in lymphoid differentiation and homeostasis [23, 33, 35, 52, 53]. Although broad chromatin modification responses were seen with both BCR (demethylation of H3K27me3) and LPS (a greater increase in H3K4me3 methylation than with BCR) activation states, both of these effects were increased at response specific promoters, suggesting signal-dependent enhancement. While increased H3K4me3 during both responses was expected, a specific decrease in H3K27me3 around the TSS of BCR-induced genes is intriguing. Because EZH2/PRC2 is involved in germinal center B cell activation, where BCR signaling would be most likely to dominate , we speculate that PRC2 plays a role in distinguishing BCR versus LPS responses (Fig. 7). lncRNAs may guide PRC2 activity (reviewed in ), indicating a point of convergence between the increased lncRNA transcription and decreased H3K27me3 that we observed during BCR activation (Additional file 1: Figure S3). Conversely, the H3K27me3-specific demethylase, JMJD3/KDM6, is involved in LPS stimulation in macrophages and is recruited to the TSS [56–58]. It is thus possible that a balance of EZH2 and KDM6 counter regulation separates these two transcriptional landscapes.
Given that helix-loop-helix (HLH) transcription factors play an important role in B cell development and differentiation [59, 60], it is intriguing that HLH TF motifs are strongly represented in BCR-responsive promoters. While much work has been done on the role of E-box binding proteins in early B cell development (reviewed in ), the role of this class of proteins in mature B cell early activation is limited . The idea that enhancer-promoter actions mediated by HLH TFs (e.g. Myc) via E-boxes might partially account for the increased RNA Pol II occupancy at promoters induced by BCR is appealing. A recent report that targeting BET proteins in high-risk acute lymphoblastic leukemia inhibits Myc and Il7r expression, both of which exhibit increased transcription in response to BCR in our experiment, also suggests Myc plays an important role during early phases of B cell activation . The fact that E-box containing promoter sequences are underrepresented in an EZH2 recruitment assay  further suggests that the decrease in H3K27me3 and preponderance of E-box sequences in BCR-induced genes observed in our study are related.
As deregulation of B cell activation is related to malignancies such as B cell lymphomas, our studies may also provide insights into lymphogenesis. Myc and NF-κB are well-established master regulators of initiation of transcriptional programs, but when deregulated, they function as oncogenic drivers in B cell lymphomas. Deregulated and increased Myc and Bcl proteins, such as Bcl2 and Bcl6, are associated with particularly aggressive lymphoma types [64, 65]. Here, we found BCR stimulation decreased Bcl6 expression; in contrast, Bcl6 remained stable during LPS activation. These and other Bcl transcription patterns (Additional file 1: Figure S11) suggest that the proper regulation of Myc and Bcl proteins is required for early activation. Further investigation into this oncogenic driver network might yield interesting relationships.
Targeting of miRNAs is complex with an average miRNA having approximately 100 target sites in addition to non-canonical miRNA binding [66, 67]. Here, we identified differential expression of many miRNAs known to regulate processes involved in B cell activation. That a greater decrease in miRNA expression occurred in the BCR response compared to LPS suggests rapid downregulation of miRNAs is necessary to orchestrate gene expression driving the adaptive immune response. Given the wide activity spectrum predicted for many miRNAs, it is possible they could coordinate some of the separate regulatory mechanisms we observed. For example, a recent report proposes a regulatory loop linking overexpression of Myc, EZH2, and miR26a repression to lymphoma growth . Our BCR activation data showing overexpression of Myc, lowered H3K27me3, and decreased miR26a highlight the multi-factorial nature and cross-dependency of regulatory systems likely to drive complicated responses such as signal-specific B cell activation. Although miRNAs have generally been associated with oncogenic pathways, targeted deletion of miR-17 cluster shows defects in B cell differentiation [39, 69]. Because LPS but not BCR signaling in splenic B cells results in differentiation, it is tempting to speculate that signal-specific regulation in the miR-17 cluster is a way of distinguishing between the two signals. The miR-15 cluster belongs to a very selective group of miRNAs enriched in the nucleus and thus capable of further directly regulating LPS specific transcription . Despite the fact that miRNAs are critical regulators of diverse biological processes, differential regulation of miRNAs to the extent observed in our analysis is very surprising. However, it is currently unknown if these miRNAs are regulated by rapid turnover of miRNAs, regulated at the level of transcription, or both.
BCR responses are slower and presumably more precise [44, 45, 71], therefore, it is tempting to speculate that a tighter regulatory environment is required to orchestrate these lengthy responses (Fig. 7). An increase in global Pol II recruitment that appears to be regulated at the level of pausing and a greater release of miRNA repression could reflect this strict regulation. Conversely, TLR4-mediated signaling is reflective of innate responses, which are generally rapid and transient, and therefore could be manifested by an enhanced transition of RNA Pol II from initiation to elongation together with a global increase in activation marks at TSSs. Additionally, one would expect the derepression (via decreasing H3K27me3) observed during the BCR response to be slower than the sharper increase in preexisting H3K4me3-activating marks observed during LPS/TLR4 signaling. How these different observations are related to each other will be the next challenging phase to understand the regulation of B cell activation. Nevertheless, our observations begin to elucidate the signal-specific signatures involved in early activation of B cells and further suggest key molecular mechanisms (Fig. 7) that govern this important process.
We conclude that the B cell commitment to different activation states is dependent upon rapid regulatory mechanisms and occurs much earlier than previously thought. Different RNA Pol II recruitment and transition from initiation to elongation, distinct activating (H3K4me3) and repressing (H3K27me3) histone signatures, mutually exclusive transcription factor binding in promoters and highly selective miRNA profiles distinguish these responses.
Cells and induction
Naïve resting B cells from splenocytes of 8-week-old male C57BL6 mice were isolated with anti-CD43 beads (Miltenyi), confirmed as 95 % CD19+ by flow cytometry (FACS Calibur), and resuspended in cold media with either 10 ug/ml anti-mouse IgM goat IgG Fab fragments (Jackson Immunology) or 25 ug/ml Salmonella typhimurium typhus LPS (Sigma) were added. The cells were rested on ice for 30 min following a previously published method  and incubated at 37 °C/5 %CO2 for the experimental times. Animal care and use in this study are covered under the “Assurance of Compliance with PHS (USA) Policy on Humane Care and Use of Laboratory animals by Awardee Institutions” and approved by the Institutional Animal Care and Use Committee of Tufts University (Animal Welfare Assurance Number A-3775-01).
Sample preparation was performed using common techniques. In general, single end, 100 bp (initial RNA), and 50 bp (secondary RNA, ChIP and miRNA) reads were mapped against the mm9/ENSEMBL build 67 genome reference using Tophat v2.0.0  and for RNA, bowtie 1.0.0 for ChIP . RNA Pol II ChIP-seq employed antibody against total RNA Pol II (Santa Cruz N-20, sc-816x), H3K4me3 with Abcam antibody ab8580, and H3K27me3 with Abcam antibody ab6002. Mapped read numbers per million and BCR or LPS time points are 120 min unless indicated. RNA , rest 75.3, BCR30 20.9, BCR120 73.6, LPS30 44.7, LPS120 40.5; RNA , rest 77, BCR 61, LPS 68.9; RNA Pol II, rest 18.1, BCR 14.2, LPS 21.3; H3K4me3, rest 13.1, BCR 16.4, LPS 19.1; K3K27me3, rest 19.0, BCR 20.0, LPS 18.1.
Differential expression analysis
Differential expression (DE) was identified by a minimal twofold difference in log ratios of normalized reads generated with Cufflinks v1.3.1. Preferentially induced or reduced genes sets included genes that were either changed by either a single response or when affected by both responses changed only two- to fourfold by one response and were changed by the preferred response at a ratio of at least twofold more than the non-referred response. A spreadsheet of the differential expression list can be found in Additional file 3. Genes were annotated to biological processes with the online Toppfun program. Gene network analysis was carried out using ToppCluster  and visualized by Cytoscape .
After Trizol isolation of RNA, TruSeq Small RNA Sample Preparation Kits were used to produce material for generating 50 bp single end reads which were then analyzed with miRDeep2  using the miRBase reference v14 with standard settings. Mapped miRNAs were confirmed by visual inspection of miRNA structure and UCSC Genome Browser tracks  and inclusion in the Ensembl data base . Differential expression from the resting state was identified by a minimal twofold difference in miRdeep2 normalized reads. Total miRNA data set reads per million are the following: rest 20.2, BCR 120 28.6, and LPS 120 14.8. Total miRdeep2 miRNA reads (per thousand) are the following: rest 55.8, BCR 16.2, and LPS 27.6.
For histograms of TSS coverage, custom R scripts were used to produce bedgraphs from mapped bam files, which were converted to BigWig files with bedgraphToBigwig for UCSC Genome Browser presentation [80, 81]. Reads per million-normalized coverage was computed for the gene sets and regions indicated, and summary statistics were calculated at each base pair, for histograms or by summing total coverage across regions, as shown for boxplots. Traveling ratios (TR) were calculated from the mean of summed transcript RPM means for each transcription group in an area representing the promoter (p) (−0.3/0.3 kbp) and downstream body (b) (0.3/2.25 kbp) of the transcript, TR = (b/p). TRs were then normalized to the resting state (TRactivation/TRrest).
CpG and TF motif analysis
Predicted CpG island locations were from preloaded USCS Genome Browser tracks and produced by common methods . Proximity of TSSs to CpG islands was analyzed with Bedtools’ IntersectBed . Enriched TF binding motifs in the promoters, defined as −1000 to +1000 regions relative to the TSS based on RNA Pol II occupancy (Fig. 2), employed the motif enrichment algorithm implemented in the HOMER tool  supplemented with the mouse TF binding motifs contained in the CisBP database (build 0.90) , resulting in a total of 3812 mouse motifs. Enrichment calculations used promoter sequences of genes whose expression did not change as our background set.
Quantitative PCR-RNA validation
The sequences have been deposited to the GEO database (NCBI/NLM/NIH)—accession number (GSE61608) (http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE61608).
Carbone A, Gloghini A, Kwong YL, Younes A. Diffuse large B cell lymphoma: using pathologic and molecular biomarkers to define subgroups for novel therapy. Ann Hematol. 2014;93(8):1263–77. Epub 2014/05/30.
Shaffer 3rd AL, Young RM, Staudt LM. Pathogenesis of human B cell lymphomas. Annu Rev Immunol. 2012;30:565–610. Epub 2012/01/10.
Sprent J. Lifespans of naive, memory and effector lymphocytes. Curr Opin Immunol. 1993;5(3):433–8. Epub 1993/06/01.
DeFranco AL, Raveche ES, Paul WE. Separate control of B lymphocyte early activation and proliferation in response to anti-IgM antibodies. J Immunol. 1985;135(1):87–94. Epub 1985/07/01.
Harwood NE, Batista FD. Early events in B cell activation. Annu Rev Immunol. 2010;28:185–210. Epub 2010/03/03.
Nutt SL, Tarlinton DM. Germinal center B and follicular helper T cells: siblings, cousins or just good friends? Nat Immunol. 2011;12(6):472–7. Epub 2011/07/09.
Doyle SL, O’Neill LA. Toll-like receptors: from the discovery of NFkappaB to new insights into transcriptional regulations in innate immunity. Biochem Pharmacol. 2006;72(9):1102–13. Epub 2006/08/26.
Pone EJ, Zhang J, Mai T, White CA, Li G, Sakakura JK, et al. BCR-signalling synergizes with TLR-signalling for induction of AID and immunoglobulin class-switching through the non-canonical NF-kappaB pathway. Nat Commun. 2012;3:767. Epub 2012/04/05.
Rawlings DJ, Schwartz MA, Jackson SW, Meyer-Bahlburg A. Integration of B cell responses through Toll-like receptors and antigen receptors. Nat Rev Immunol. 2012;12(4):282–94. Epub 2012/03/17.
Dardick I, Sinnott NM, Hall R, Bajenko-Carr TA, Setterfield G. Nuclear morphology and morphometry of B-lymphocyte transformation. Implications for follicular center cell lymphomas. Am J Pathol. 1983;111(1):35–49. Epub 1983/04/01.
Kurosaki T, Shinohara H, Baba Y. B cell signaling and fate decision. Annu Rev Immunol. 2010;28:21–55. Epub 2009/10/16.
Niiro H, Clark EA. Regulation of B-cell fate by antigen-receptor signals. Nat Rev Immunol. 2002;2(12):945–56. Epub 2002/12/04.
Kouzine F, Wojtowicz D, Yamane A, Resch W, Kieffer-Kwon KR, Bandle R, et al. Global regulation of promoter melting in naive lymphocytes. Cell. 2013;153(5):988–99. Epub 2013/05/28.
Fowler T, Sen R, Roy AL. Regulation of primary response genes. Mol Cell. 2011;44(3):348–60. Epub 2011/11/08.
Rogatsky I, Adelman K. Preparing the first responders: building the inflammatory transcriptome from the ground up. Mol Cell. 2014;54(2):245–54. Epub 2014/04/29.
Glynne R, Ghandour G, Rayner J, Mack DH, Goodnow CC. B-lymphocyte quiescence, tolerance and activation as viewed by global gene expression profiling on microarrays. Immunol Rev. 2000;176:216–46. Epub 2000/10/24.
Kearney JF, Lawton AR. B lymphocyte differentiation induced by lipopolysaccharide. I Generation of cells synthesizing four major immunoglobulin classes. J Immunol. 1975;115(3):671–6. Epub 1975/09/01.
Solvason N, Wu WW, Kabra N, Wu X, Lees E, Howard MC. Induction of cell cycle regulatory proteins in anti-immunoglobulin-stimulated mature B lymphocytes. J Exp Med. 1996;184(2):407–17. Epub 1996/08/01.
Corcoran AE, Smart FM, Cowling RJ, Crompton T, Owen MJ, Venkitaraman AR. The interleukin-7 receptor alpha chain transmits distinct signals for proliferation and differentiation during B lymphopoiesis. EMBO J. 1996;15(8):1924–32. Epub 1996/04/15.
Bandow K, Kusuyama J, Shamoto M, Kakimoto K, Ohnishi T, Matsuguchi T. LPS-induced chemokine expression in both MyD88-dependent and -independent manners is regulated by Cot/Tpl2-ERK axis in macrophages. FEBS Lett. 2012;586(10):1540–6. Epub 2012/06/08.
Reppas NB, Wade JT, Church GM, Struhl K. The transition between transcriptional initiation and elongation in E. coli is highly variable and often rate limiting. Mol Cell. 2006;24(5):747–57. Epub 2006/12/13.
Wade JT, Struhl K. The transition from transcriptional initiation to elongation. Curr Opin Genet Dev. 2008;18(2):130–6. Epub 2008/02/20.
Smale ST, Tarakhovsky A, Natoli G. Chromatin contributions to the regulation of innate immunity. Annu Rev Immunol. 2014;32:489–511. Epub 2014/02/22.
Shilatifard A. The COMPASS family of histone H3K4 methylases: mechanisms of regulation in development and disease pathogenesis. Annu Rev Biochem. 2012;81:65–95. Epub 2012/06/06.
Baxter J, Sauer S, Peters A, John R, Williams R, Caparros ML, et al. Histone hypomethylation is an indicator of epigenetic plasticity in quiescent lymphocytes. EMBO J. 2004;23(22):4462–72. Epub 2004/10/29.
Boyer LA, Plath K, Zeitlinger J, Brambrink T, Medeiros LA, Lee TI, et al. Polycomb complexes repress developmental regulators in murine embryonic stem cells. Nature. 2006;441(7091):349–53. Epub 2006/04/21.
Landeira D, Fisher AG. Inactive yet indispensable: the tale of Jarid2. Trends Cell Biol. 2011;21(2):74–80. Epub 2010/11/16.
Pereira CF, Piccolo FM, Tsubouchi T, Sauer S, Ryan NK, Bruno L, et al. ESCs require PRC2 to direct the successful reprogramming of differentiated cells toward pluripotency. Cell Stem Cell. 2010;6(6):547–56. Epub 2010/06/24.
Simon JA, Kingston RE. Mechanisms of polycomb gene silencing: knowns and unknowns. Nat Rev Mol Cell Biol. 2009;10(10):697–708. Epub 2009/09/10.
Cao R, Wang L, Wang H, Xia L, Erdjument-Bromage H, Tempst P, et al. Role of histone H3 lysine 27 methylation in Polycomb-group silencing. Science. 2002;298(5595):1039–43. Epub 2002/09/28.
van Galen JC, Dukers DF, Giroth C, Sewalt RG, Otte AP, Meijer CJ, et al. Distinct expression patterns of polycomb oncoproteins and their binding partners during the germinal center reaction. Eur J Immunol. 2004;34(7):1870–81. Epub 2004/06/24.
Velichutina I, Shaknovich R, Geng H, Johnson NA, Gascoyne RD, Melnick AM, et al. EZH2-mediated epigenetic silencing in germinal center B cells contributes to proliferation and lymphomagenesis. Blood. 2010;116(24):5247–55. Epub 2010/08/26.
Hargreaves DC, Horng T, Medzhitov R. Control of inducible gene expression by signal-dependent transcriptional elongation. Cell. 2009;138(1):129–45. Epub 2009/07/15.
Deaton AM, Bird A. CpG islands and the regulation of transcription. Genes Dev. 2011;25(10):1010–22. Epub 2011/05/18.
Ramirez-Carrozzi VR, Braas D, Bhatt DM, Cheng CS, Hong C, Doty KR, et al. A unifying model for the selective regulation of inducible transcription by CpG islands and nucleosome remodeling. Cell. 2009;138(1):114–28. Epub 2009/07/15.
Mendenhall EM, Koche RP, Truong T, Zhou VW, Issac B, Chi AS, et al. GC-rich sequence elements recruit PRC2 in mammalian ES cells. PLoS Genet. 2010;6(12), e1001244. Epub 2010/12/21.
Contreras J, Rao DS. MicroRNAs in inflammation and immune responses. Leukemia. 2012;26(3):404–13. Epub 2011/12/21.
Danger R, Braza F, Giral M, Soulillou JP, Brouard MRNAS. Major Players in B Cells Homeostasis and Function. Front Immunol. 2014;5:98. Epub 2014/03/22.
Ventura A, Young AG, Winslow MM, Lintault L, Meissner A, Erkeland SJ, et al. Targeted deletion reveals essential and overlapping functions of the miR-17 through 92 family of miRNA clusters. Cell. 2008;132(5):875–86. Epub 2008/03/11.
Rodriguez A, Vigorito E, Clare S, Warren MV, Couttet P, Soond DR, et al. Requirement of bic/microRNA-155 for normal immune function. Science. 2007;316(5824):608–11. Epub 2007/04/28.
Enomoto Y, Kitaura J, Hatakeyama K, Watanuki J, Akasaka T, Kato N, et al. Emu/miR-125b transgenic mice develop lethal B-cell malignancies. Leukemia. 2011;25(12):1849–56. Epub 2011/07/09.
Mao YS, Sunwoo H, Zhang B, Spector DL. Direct visualization of the co-transcriptional assembly of a nuclear body by noncoding RNAs. Nat Cell Biol. 2011;13(1):95–101. Epub 2010/12/21.
Tripathi V, Ellis JD, Shen Z, Song DY, Pan Q, Watt AT, et al. The nuclear-retained noncoding RNA MALAT1 regulates alternative splicing by modulating SR splicing factor phosphorylation. Mol Cell. 2010;39(6):925–38. Epub 2010/08/28.
Fambrough D, McClure K, Kazlauskas A, Lander ES. Diverse signaling pathways activated by growth factor receptors induce broadly overlapping, rather than independent, sets of genes. Cell. 1999;97(6):727–41. Epub 1999/06/25.
Zhu X, Hart R, Chang MS, Kim JW, Lee SY, Cao YA, et al. Analysis of the major patterns of B cell gene expression changes in response to short-term stimulation with 33 single ligands. J Immunol. 2004;173(12):7141–9. Epub 2004/12/09.
Cattoretti G. MYC expression and distribution in normal mature lymphoid cells. J Pathol. 2013;229(3):430–40. Epub 2012/11/21.
Link JM, Hurlin PJ. The activities of MYC, MNT and the MAX-interactome in lymphocyte proliferation and oncogenesis. Biochimica et biophysica acta. 2014. Epub 2014/04/16.
Eilers M, Eisenman RN. Myc’s broad reach. Genes Dev. 2008;22(20):2755–66. Epub 2008/10/17.
Dang CV. MYC on the path to cancer. Cell. 2012;149(1):22–35. Epub 2012/04/03.
Sabo A, Kress TR, Pelizzola M, de Pretis S, Gorski MM, Tesi A, et al. Selective transcriptional regulation by Myc in cellular growth control and lymphomagenesis. Nature. 2014;511(7510):488–92. Epub 2014/07/22.
Ehrensberger AH, Kelly GP, Svejstrup JQ. Mechanistic interpretation of promoter-proximal peaks and RNAPII density maps. Cell. 2013;154(4):713–5. Epub 2013/08/21.
Bhatt D, Ghosh S. Regulation of the NF-kappaB-mediated transcription of inflammatory genes. Front Immunol. 2014;5:71. Epub 2014/03/13.
Smith E, Shilatifard A. The chromatin signaling pathway: diverse mechanisms of recruitment of histone-modifying enzymes and varied biological outcomes. Mol Cell. 2010;40(5):689–701. Epub 2010/12/15.
Beguelin W, Popovic R, Teater M, Jiang Y, Bunting KL, Rosen M, et al. EZH2 is required for germinal center formation and somatic EZH2 mutations promote lymphoid transformation. Cancer Cell. 2013;23(5):677–92. Epub 2013/05/18.
Simon JA, Kingston RE. Occupying chromatin: polycomb mechanisms for getting to genomic targets, stopping transcriptional traffic, and staying put. Mol Cell. 2013;49(5):808–24. Epub 2013/03/12.
Chen S, Ma J, Wu F, Xiong LJ, Ma H, Xu W, et al. The histone H3 Lys 27 demethylase JMJD3 regulates gene expression by impacting transcriptional elongation. Genes Dev. 2012;26(12):1364–75. Epub 2012/06/21.
De Santa F, Narang V, Yap ZH, Tusi BK, Burgold T, Austenaa L, et al. Jmjd3 contributes to the control of gene expression in LPS-activated macrophages. EMBO J. 2009;28(21):3341–52. Epub 2009/09/26.
De Santa F, Totaro MG, Prosperini E, Notarbartolo S, Testa G, Natoli G. The histone H3 lysine-27 demethylase Jmjd3 links inflammation to inhibition of polycomb-mediated gene silencing. Cell. 2007;130(6):1083–94. Epub 2007/09/11.
Anantharaman A, Lin IJ, Barrow J, Liang SY, Masannat J, Strouboulis J, et al. Role of helix-loop-helix proteins during differentiation of erythroid cells. Mol Cell Biol. 2011;31(7):1332–43. Epub 2011/02/02.
Massari ME, Murre C. Helix-loop-helix proteins: regulators of transcription in eucaryotic organisms. Mol Cell Biol. 2000;20(2):429–40. Epub 1999/12/28.
Miyazaki K, Miyazaki M, Murre C. The establishment of B versus T cell identity. Trends Immunol. 2014;35(5):205–10. Epub 2014/04/01.
McDonald JJ, Alinikula J, Buerstedde JM, Schatz DG. A critical context-dependent role for E boxes in the targeting of somatic hypermutation. J Immunol. 2013;191(4):1556–66. Epub 2013/07/10.
Ott CJ, Kopp N, Bird L, Paranal RM, Qi J, Bowman T, et al. BET bromodomain inhibition targets both c-Myc and IL7R in high-risk acute lymphoblastic leukemia. Blood. 2012;120(14):2843–52. Epub 2012/08/21.
Ott G, Rosenwald A, Campo E. Understanding MYC-driven aggressive B-cell lymphomas: pathogenesis and classification. Hematology Am Soc Hematol Educ Program. 2013;2013:575–83. Epub 2013/12/10.
Tomita N. BCL2 and MYC dual-hit lymphoma/leukemia. J Clin Exp Hematop. 2011;51(1):7–12. Epub 2011/06/02.
Brennecke J, Stark A, Russell RB, Cohen SM. Principles of microRNA-target recognition. PLoS Biol. 2005;3(3), e85. Epub 2005/02/22.
Loeb GB, Khan AA, Canner D, Hiatt JB, Shendure J, Darnell RB, et al. Transcriptome-wide miR-155 binding map reveals widespread noncanonical microRNA targeting. Mol Cell. 2012;48(5):760–70. Epub 2012/11/13.
Zhao X, Lwin T, Zhang X, Huang A, Wang J, Marquez VE, et al. Disruption of the MYC-miRNA-EZH2 loop to suppress aggressive B-cell lymphoma survival and clonogenicity. Leukemia. 2013;27(12):2341–50. Epub 2013/03/30.
Sandhu SK, Volinia S, Costinean S, Galasso M, Neinast R, Santhanam R, et al. miR-155 targets histone deacetylase 4 (HDAC4) and impairs transcriptional activity of B-cell lymphoma 6 (BCL6) in the Emu-miR-155 transgenic mouse model. Proc Natl Acad Sci U S A. 2012;109(49):20047–52.
Taft RJ, Simons C, Nahkuri S, Oey H, Korbie DJ, Mercer TR, et al. Nuclear-localized tiny RNAs are associated with transcription initiation and splice sites in metazoans. Nat Struct Mol Biol. 2010;17(8):1030–4. Epub 2010/07/14.
Janeway Jr CA, Medzhitov R. Innate immune recognition. Annu Rev Immunol. 2002;20:197–216. Epub 2002/02/28.
Damdinsuren B, Zhang Y, Khalil A, Wood 3rd WH, Becker KG, Shlomchik MJ, et al. Single round of antigen receptor signaling programs naive B cells to receive T cell help. Immunity. 2010;32(3):355–66. Epub 2010/03/17.
Trapnell C, Pachter L, Salzberg SL. TopHat: discovering splice junctions with RNA-Seq. Bioinformatics. 2009;25(9):1105–11. Epub 2009/03/18.
Langmead B, Trapnell C, Pop M, Salzberg SL. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009;10(3):R25. Epub 2009/03/06.
Kaimal V, Bardes EE, Tabar SC, Jegga AG, Aronow BJ. ToppCluster: a multiple gene list feature analyzer for comparative enrichment clustering and network-based dissection of biological systems. Nucleic Acids Res. 2010;38(Web Server issue):W96–102. Epub 2010/05/21.
Cline MS, Smoot M, Cerami E, Kuchinsky A, Landys N, Workman C, et al. Integration of biological networks and gene expression data using Cytoscape. Nat Protoc. 2007;2(10):2366–82. Epub 2007/10/20.
Friedlander MR, Mackowiak SD, Li N, Chen W, Rajewsky N. miRDeep2 accurately identifies known and hundreds of novel microRNA genes in seven animal clades. Nucleic Acids Res. 2012;40(1):37–52. Epub 2011/09/14.
Raney BJ, Dreszer TR, Barber GP, Clawson H, Fujita PA, Wang T, et al. Track data hubs enable visualization of user-defined genome-wide annotations on the UCSC Genome Browser. Bioinformatics. 2014;30(7):1003–5. Epub 2013/11/15.
Flicek P, Ahmed I, Amode MR, Barrell D, Beal K, Brent S, et al. Ensembl 2013. Nucleic Acids Res. 2013;41(Database issue):D48–55. Epub 2012/12/04.
Quinlan AR, Hall IM. BEDTools: a flexible suite of utilities for comparing genomic features. Bioinformatics. 2010;26(6):841–2. Epub 2010/01/30.
Zhang Y, Liu T, Meyer CA, Eeckhoute J, Johnson DS, Bernstein BE, et al. Model-based analysis of ChIP-Seq (MACS). Genome Biol. 2008;9(9):R137. Epub 2008/09/19.
Gardiner-Garden M, Frommer M. CpG islands in vertebrate genomes. J Mol Biol. 1987;196(2):261–82. Epub 1987/07/20.
Heinz S, Benner C, Spann N, Bertolino E, Lin YC, Laslo P, et al. Simple combinations of lineage-determining transcription factors prime cis-regulatory elements required for macrophage and B cell identities. Mol Cell. 2010;38(4):576–89. Epub 2010/06/02.
Weirauch MT, Yang A, Albu M, Cote AG, Montenegro-Montero A, Drewe P, et al. Determination and inference of eukaryotic transcription factor sequence specificity. Cell. 2014;158(6):1431–43. Epub 2014/09/13.
Fowler T, Suh H, Buratowski S, Roy AL. Regulation of primary response genes in B cells. J Biol Chem. 2013;288(21):14906–16. Epub 2013/03/29.
We are deeply indebted to Drs. Ali Shilatifard, Deqing Hu, and Xin Gao of the Stowers Institute for performing an independent RNA-seq and ChIP-seq of RNA Pol II as well as their guidance and helpful discussions throughout the course of this work. We thank Xiaoting Chen (CCHMC) for help with the motif enrichment analysis and members of Tufts’ Computational Biology Initiative for analysis assistance. We also thank Dr. Dinah Singer (NCI) and Dr. Harinder Singh (CCHMC) for critically reading the manuscript and for their thoughtful suggestions. This work was supported in part by the Intramural Research Program of the NIH, the National Institute on Aging to RS, and by grants to STS (R01 GM086372) and ALR (AHA 12GRNT12180023).
The authors declare that they have no competing interests.
Concept was initiated by ALR. TF developed and performed experiments and analyzed data. ASG, SD, KGB, MTW, WHW, MTW, BA, STS, and RS analyzed the data. AG performed experiments. STS and RS also provided advice and material support. TF and ALR wrote the manuscript with input from all authors. All authors read and approved the manuscript.
Figure S1. 120 min is required for recognizably organized response. Figure S2. Characterization of unchanged transcription. Figure S3. Stimulus specific of high polymerase occupancy. Figure S4. RNA Pol II data for each cellular state. Figure S5. Global RNA Pol II recruitment and histone changes. Figure S6. Stimulus-specific histone methylation at high Pol occupancy TSS(s). Figure S7. Closer examination of individual BCR induced gene tracks to better show examples of H3K27me3 reduction upon BCR activation. Figure S8. Lack of CpG-dependent effect on Pol II occupancy. Figure S9. Expression of activation-responsive miRNAs. Figure S10. Observation of lncRNA Malat1 behavior with validation and correlating activation regulation mechanisms. Figure S11. Response-dependent difference in Bcl network.
Genes involved in NF-κB signaling. A list of genes induced by either BCR or LPS, or both (shared), found in the KEGG Pathway NF-κB signaling pathway.
About this article
Cite this article
Fowler, T., Garruss, A.S., Ghosh, A. et al. Divergence of transcriptional landscape occurs early in B cell activation. Epigenetics & Chromatin 8, 20 (2015). https://doi.org/10.1186/s13072-015-0012-x