ARTICLE https://doi.org/10.1038/s42003-019-0284-y OPEN Changes in long-range rDNA-genomic interactions associate with altered RNA polymerase II gene programs during malignant transformation 1234567890():,; Jeannine Diesch1,11, Megan J. Bywater1,12, Elaine Sanij1,2, Donald P. Cameron1,3,4, William Schierding5, Natalie Brajanovski1, Jinbae Son1,4, Jirawas Sornkom 1,4, Nadine Hein3, Maurits Evers3, Richard B. Pearson 1,4,6,7, Grant A. McArthur1,4,8, Austen R.D. Ganley9, Justin M. O’Sullivan 5, Ross D. Hannan1,3,4,6,7,10 & Gretchen Poortinga 1,4,8 The three-dimensional organization of the genome contributes to its maintenance and reg- ulation. While chromosomal regions associate with nucleolar ribosomal RNA genes (rDNA), the biological significance of rDNA-genome interactions and whether they are dynamically regulated during disease remain unclear. rDNA chromatin exists in multiple inactive and active states and their transition is regulated by the RNA polymerase I transcription factor UBTF. Here, using a MYC-driven lymphoma model, we demonstrate that during malignant progression the rDNA chromatin converts to the open state, which is required for tumor cell survival. Moreover, this rDNA transition co-occurs with a reorganization of rDNA-genome contacts which correlate with gene expression changes at associated loci, impacting gene ontologies including B-cell differentiation, cell growth and metabolism. We propose that UBTF-mediated conversion to open rDNA chromatin during malignant transformation con- tributes to the regulation of specific gene pathways that regulate growth and differentiation through reformed long-range physical interactions with the rDNA. 1 Cancer Research Division, Peter MacCallum Cancer Centre, Melbourne, VIC 3000, Australia. 2 Department of Pathology, University of Melbourne, Parkville, VIC 3010, Australia. 3 ACRF Department of Cancer Biology and Therapeutics, John Curtin School of Medical Research, Australian National University, Canberra, ACT 2601, Australia. 4 Sir Peter MacCallum Department of Oncology, University of Melbourne, Parkville, VIC 3010, Australia. 5 Liggins Institute, The University of Auckland, Auckland 1023, New Zealand. 6 Department of Biochemistry and Molecular Biology, Monash University, Clayton 3800 VIC, Australia. 7 Department of Biochemistry and Molecular Biology, University of Melbourne, Parkville, VIC 3010, Australia. 8 Department of Medicine, St Vincent’s Hospital, University of Melbourne, Fitzroy, VIC 3065, Australia. 9 School of Biological Sciences, The University of Auckland, Auckland 1010, New Zealand. 10 School of Biomedical Sciences, University of Queensland, Brisbane, QLD 4072, Australia. 11Present address: Josep Carreras Leukaemia Research Institute, Barcelona 08021, Spain. 12Present address: QIMR Berghofer Medical Research Institute, Brisbane, QLD 4029, Australia. Correspondence and requests for materials should be addressed to R.D.H. (email: ross.hannan@anu.edu.au) or to G.P. (email: gretchen.poortinga@petermac.org) COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio 1 ARTICLE COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y A dvances in genomics and epigenomics have provided cell lymphoma, we showed that progression from premalignancy insights into three-dimensional (3D) genome organiza- to malignancy in vivo is associated with UBTF-dependent epi- tion at an unprecedented level of detail1–3, and highlight genetic remodeling that activates a significant proportion of the dynamic relationship between genomic spatial organization previously pseudo-silent rDNA repeats. Circularized chromo- and gene regulation4,5. However, relatively little is understood some conformation capture sequencing (4C-seq) demonstrated about how genome organization and gene expression are that, concomitant with activation of rDNA during malignancy, reshaped during disease development and their impact on the the population of genomic loci interacting with the rDNA process6. The largest substructure in the nucleus is the nucleolus, changes during lymphomagenesis. Genes associated with these the site of ribosome biogenesis, which forms dynamically around rDNA-interacting loci show an inverse relationship between their transcribed ribosomal RNA (rRNA) genes (rDNA) arranged in rDNA interaction level and their gene expression. A sub- arrays of tandem repeats at chromosomal regions called nucleolar population of the rDNA-NAD interactions that change during organizer regions (NORs)7. The NORs are organized on acro- malignant progression require the active chromatin state of the centric chromosomes, with the highly variable rDNA copy rDNA repeats. Notably, of this NAD-associated loci sub-popu- number averaging more than 100 per diploid genome7,8. How- lation, those with reduced expression in malignant cells are ever, at any given time less than 50% of these rDNA copies are enriched for genes encoding proteins involved in B-cell differ- actively transcribed by the dedicated RNA polymerase I (Pol I) to entiation. Conversely, those with increased expression are enri- produce the 47S rRNA precursor9. A large body of evidence ched for growth and metabolism gene ontologies. We propose supports a model of rDNA copies existing in one of three epi- that malignancy-associated alterations in the rDNA chromatin genetic chromatin states: silent (CpG methylated at the rDNA status are linked to formation of new rDNA-genomic interactions promoter), pseudo-silent (lacking hypomethylated promoters but and co-occur with gene expression changes that contribute to the transcriptionally inactive), or active (hypomethylated and tran- malignant phenotype. scriptionally competent)10. Active genes exhibit a range of tran- scription rates depending on the cellular state, and the relative ratio of silent, pseudo-silent and active genes is modulated during Results differentiation and development11–15. These states are regulated rDNA chromatin opens during Eµ-Myc malignant progression. epigenetically and via the upstream binding transcription factor Mice harbouring the Eµ-Myc transgene display abnormal B-cell (UBTF), an architectural protein required to recruit Pol I to the development, characterized by an initial premalignant phase in rDNA promoter but also critical for binding to under-methylated, young mice (4–6 weeks) of enhanced B-cell progenitor prolifera- pseudo-silent rDNA repeats to convert them into the active, tion and growth before eventually progressing to malignant lym- transcriptionally competent state12,14,16,17. phoma54–56. Both premalignant and malignant Eµ-Myc B-cells Functions beyond production of rRNA are well documented have elevated rRNA synthesis rates compared to wild-type B cells for the rDNA and nucleolus, including regulation of genomic of the same developmental stage due to elevated expression levels stability and global gene expression18–23. Genomic sequences that of the MYC oncogene37. To investigate the status of rDNA include certain genes other than rDNA are localized to nucleoli in chromatin with respect to rRNA transcription rates during nucleolar-associated domains (NADs)24–28. The NAD nucleolar malignant progression, we compared their dynamics in wild type, chromatin compartment is enriched for repressive chromatin premalignant, and malignant cells isolated from Eµ-Myc mice. We marks and under-represented for active histone modifications, separated these changes from variables associated with differ- and NAD-associated genes are generally transcriptionally entiation status and tumor compartment by performing analyses repressed24–27. This is consistent with evidence that the nucleolus on IgM-negative pre-B cell subtypes (B220highIgMlowIgDlow) is surrounded by a facultative heterochromatic shell26,29,30, and FACS sorted from bone marrow (with individual animals as highlights a potential role for the nucleolus in dynamically reg- biological replicates). MYC overexpression induced a robust ulating global gene transcription through nucleolar colocalization. increase in rRNA transcription rates in premalignant cells com- The interplay between altered nucleolar morphology and dis- pared to cells from wild-type littermates, with no further sig- ease is well recognized and accelerated rRNA transcription and nificant increase observed between premalignant and malignant ribosome biogenesis is a common feature of many cancers31–33. cells as measured using qRT-PCR (Fig. 1a) or rRNA fluorescence This is reflected by the increased size and/or number of nucleoli in situ hybridization (FISH) (Fig. 1b). We assayed rDNA chro- in tumor cells, initially observed by pathologists over 100 years matin structure by performing psoralen crosslinking followed by ago and used as a diagnostic and prognostic marker for certain southern blot analysis. Psoralen preferentially incorporates into cancers34. Moreover, the potential of dysregulated ribosome more accessible DNA regions and distinguishes between active biogenesis as a therapeutic target in cancer has been demon- (i.e., transcriptionally competent open rDNA chromatin) and strated by the development of small molecule inhibitors of Pol I inactive (i.e., silent and pseudo-silent) rRNA genes11. In contrast transcription35–39. The MYC oncogene, a potent transcriptional to the elevated rRNA transcription rates in both premalignant and driver of growth-associated gene programs40–42 via all three RNA malignant cells (Fig. 1a, b), the fraction of active rDNA repeats polymerases (Pol I, II, and III), has been implicated in sensitizing was significantly elevated in malignant cells only (malignant, MYC-addicted cancer cells to inhibition of Pol I 38.2% active vs. premalignant, 16.6% active) (Fig. 1c, d, Supple- transcription13,43–48. However, while a substantial body of data mentary Fig. 1a-d). Furthermore, high ratios of active to inactive has provided critical insight into our understanding of rDNA as a rDNA chromatin were uniformly observed in lymphoma cell lines therapeutic target36,37,39,49–53, the precise mechanisms under- established from independent Eµ-Myc tumors (lymphoma lines lying the heightened sensitivity of tumor cells to perturbations in average, 42.8% active) as compared to wild type and premalignant Pol I transcription and the degree to which subsequent disruption cells (Fig. 1d, e). Thus, changes in the ratio of active to inactive of nucleolar integrity contributes to the therapeutic window rDNA repeats specifically occur as cells transition from pre- remain unresolved. malignant to malignant and are uncoupled from MYC-driven Here we examine whether changes to rDNA chromatin rRNA transcription rates during malignant progression. structure are associated with malignant transformation and fur- Loading of the Pol I architectural factor UBTF (consisting of two ther, are accompanied by alterations in rDNA-NAD interactions. polypeptides, UBTF1 and 2) onto rDNA is required for the Using the Eµ-Myc mouse model of spontaneous MYC-driven B- formation and maintenance of active, transcriptionally competent 2 COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y ARTICLE a b 15 *** wt Prem Mal pre-rRNA expression NS 10 ** 5 0 rRNA/DAPI em t al w M c Pr d *** e Eμ-Myc B cells 100 Eμ-Myc lymphoma cell lines % Total rDNA Active Active Active 50 Inactive Inactive Inactive 2770 bp 2770 bp 0 – + + + + Psoralen em t al w M Pr s em Pr t Li l w a ne M f UBTF ChIP g * % UBTF-bound rDNA 3 100 100 % Total rDNA wt % Total DNA Prem 2 NS * Mal Unmeth Unmeth 50 50 Meth Meth 1 0 0 0 em em t al t al w w M M CORE 28S Pr Pr h 10 MeDIP wt Prem % Total DNA Mal 5 0 S1 S S S H E E 28 IG ET C R EN IT U O 5′ C Fig. 1 rDNA chromatin undergoes remodeling as Eμ-Myc cells transition to malignancy. a Relative pre-rRNA expression (rRNA transcription rate) in wild type (wt), premalignant (prem) and malignant (mal) Eμ-Myc cells determined by qRT-PCR using primers to the rapidly degraded 5’ external transcribed spacer (5’ETS) region. Pre-rRNA levels were normalized to B2M mRNA and represented relative to wt (x ± s.e.m., n = 3 mice/cell stage). b rRNA levels in cells from a as determined by rRNA FISH using a probe (red) to the degraded internal transcribed spacer (ITS1) region. Nuclei were counterstained with DAPI (blue). Scale bars, 10 μm. c Representative Southern blot of psoralen cross-linked genomic rDNA from wt, prem and mal Eμ-Myc cells; arrows indicate active and inactive rDNA repeats. d Four independent psoralen Southern blot experiments as in c were quantitated with ImageQuant (x ± s.e.m., n = 4). Also shown is quantitation (lines) of e Southern blot of psoralen cross-linked genomic rDNA from Eμ-Myc lymphoma cells lines (lane 1: no psoralen control, 4242 line; lanes 2–5: lines 226, 107, 4242, and 102); arrows indicate active and inactive rDNA repeats. f ChIP analysis of UBTF binding at the rDNA promoter (CORE) and transcribed (28S) regions. The percentage (%) of DNA immunoprecipitated with anti-UBTF or rabbit sera (RS) was calculated relative to the unprecipitated input control and the % of DNA associated with RS was subtracted from corresponding UBTF samples (x ± s.e.m., n = 3). g CpG methylation at the rDNA core promoter determined by HpaII digestion followed by qPCR of total rDNA (left panel) and UBTF-bound rDNA (right panel). The percentage of unmethylated and methylated rDNA is plotted (x ± s.d., total rDNA n = 3, UBTF ChIP-CHOP n = 2). h CpG methylation frequency across the rDNA unit determined by methylated DNA immunoprecipitation (MeDIP) followed by qPCR using amplicons spanning the rDNA repeat. All amplicon primers as labeled except: 5’ETS (ETS1), 28S (28S-2), and IGS (IGS2). Samples are analyzed by qPCR with results expressed as percent (%) of total DNA (minus RS % total) as described in f (x ± s.d., n = 2). The Student's t-test was used for all statistical analyses: *p-value < 0.05; **p-value < 0.01; ***p-value < 0.001; NS, not significant rDNA12,17,57. Therefore, we used quantitative ChIP to determine While CpG dinucleotide methylation of the rDNA promoter UBTF enrichment at rDNA in all stages of Eµ-Myc pre-B cells. core region is implicated in silencing of murine rDNA14,58, we Consistent with the psoralen results, UBTF enrichment was previously demonstrated that UBTF regulates the transition from significantly increased in malignant cells compared to premalignant pseudo-silent to open rDNA chromatin in the absence of this cells (Fig. 1f), with a modest increase in premalignant cells over wild methylation12,13. To investigate the methylation status at the type. These results suggest that UBTF binding at rDNA mediates rDNA promoter compared to the rDNA chromatin state, we the conversion of rDNA to the active chromatin state as cells performed methylation-sensitive digestion at the −143 CCGG transition from premalignancy to malignancy and this occurs in the HpaII site of genomic DNA from wild type, premalignant, and absence of further increase in the rRNA transcription rate. malignant cells followed by qPCR with primers spanning this site COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio 3 ARTICLE COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y and the adjacent −133 CpG methylated residue at the rDNA association with the rDNA51. Following treatment of lymphoma promoter12,59. We detected no significant change in methylation cells with CX-5461, we observed no significant change to the ratio at the −143 CpG (Fig. 1g, left), despite the striking shift in rDNA of active to inactive rDNA despite repression of rRNA transcrip- chromatin status during malignant progression (Fig. 1c, d). tion (Supplementary Fig. 2). Taken together, these data indicate Similarly, we did not identify any significant changes in that survival of malignant Eμ-Myc cells is dependent on the ability enrichment of CpG methylation at amplicons spanning the of UBTF to maintain an active rDNA chromatin state indepen- rDNA unit between cell types using quantitative methylated DNA dent of the rRNA transcription rate. immunoprecipitation (qMeDIP) (Fig. 1h). These findings indicate that the switch from inactive to active rDNA state during the rDNA repeats engage in long-range chromatin interactions. transition to malignancy does not correlate with changes in The observation that the rDNA chromatin status affects tumor rDNA CpG methylation. Because UBTF binding is associated cell survival independent of rRNA transcription rate per se sug- with active rDNA repeats, we examined the methylation status of gests that changes in rDNA chromatin might impart functional UBTF-bound regions by performing UBTF ChIP followed by effects. Non-rDNA-genomic regions, called NADs, localize to the HpaII digestion and qPCR at the −143 CpG site (ChIP-CHOP), nucleolus, which forms around actively transcribed rDNA and found that the majority of UBTF-bound rDNA (~80%) was repeats24–28. Changes in rDNA structure/activity can influence unmethylated in all cell stages (Fig. 1g, right). Collectively, these genome-wide chromatin structure and gene expression19,20, results demonstrate that the MYC-driven premalignant state, therefore we examined whether rDNA class switching promotes compared to wild type, is characterized by increased rRNA reorganization of NADs. We used 4C-seq to identify long-range production without major alterations in the rDNA chromatin chromatin interactions that occur with rDNA sequences during state. In contrast, the transition to malignancy is accompanied by lymphomagenesis. Multiple bait sequences located in the tran- the UBTF-associated activation of pseudo-silenced rDNA scribed 18S and 28S rDNA regions (Fig. 3a) were used to capture (termed herein rDNA class switching), which is uncoupled from interactions between the rDNA and rest of the genome in isolated changes in rRNA transcription rates. wild type, premalignant, and malignant pre-B cells. Cell-type specific rDNA interactions were included in subsequent analyses if they were observed in two biological replicates (i.e., animals). Eµ-Myc lymphoma survival requires UBTF-mediated active We combined all replicated rDNA interactions for each cell type rDNA. Given the observation that clonal Eμ-Myc lymphomas and calculated the average read count per interaction in 5-kb show significantly elevated ratios of active rDNA repeats, we non-overlapping windows as a semi-quantitative measure of hypothesised that rDNA class switching is favorable to the interaction frequency (read counts listed in Supplementary malignant state. To test this, we investigated whether tumor cell Data 1). The interactions were not restricted to a few specific loci survival is dependent on the elevated active:inactive rDNA ratio and did not show a bias towards the rDNA-containing chro- independent of rRNA transcription rate. We previously estab- mosomes (12, 15, 18, and 19 in the Eμ-Myc C57BL/6 genetic lished that Ubtf knockdown leads to a reduction in the number of background60), but were uniformly distributed across all chro- active rDNA repeats without affecting overall rRNA transcription mosomes (Supplementary Fig. 3a). We categorized the rDNA- rate due to compensatory mechanisms that maintain rRNA genome interactions into two categories: those that show no synthesis at a constant level12. Here, we confirmed this observa- significant change (FDR > 0.1) between wild type, premalignant, tion using a shRNAmiR targeting Ubtf (Ubtf1/2)37 in an Eμ-Myc and malignant cells (constitutive rDNA interactions, Supple- lymphoma cell line (Eμ-Myc-shUbtf). Knockdown of Ubtf led to a mentary Data 2); and those that are significantly changed (FDR < robust reduction of UBTF compared to control cell line (Eμ-Myc- 0.1) between any pair of cell stages (differential rDNA interac- LMP) (Fig. 2a) that was accompanied by a pronounced reduction tions). We plotted the top 10% of differential rDNA interactions in active rDNA repeats (Fig. 2b) but only a modest reduction in (regardless of fold-change direction) and while qualitative only, relative rRNA transcription rates (Fig. 2c, d). circos plots revealed that the major overall interaction re- To investigate whether the reversal of rDNA class switching by configuration occurs between premalignant and malignant cells Ubtf knockdown impacts cell proliferation and survival, we used (Fig. 3b). These results suggest that B-cell lymphomagenesis is an in vitro GFP-based cell competition assay37. Mock-transduced associated with changes in the pattern of interactions that rDNA Eμ-Myc lymphoma cells were seeded in culture with equal forms with the rest of genome. numbers of cells retrovirally transduced with either empty vector or shUbtf expressing vector, and survival of the vector-driven GFP expressing cells was monitored by FACS. Long-term (8 days) rDNA-genome interactions correlate with transcriptome competition revealed that shUbtf-expressing cells have an acute changes. Specific gene families are reported to be enriched within survival disadvantage compared to control cells (Fig. 2e). More- NAD sequences24,25,27, therefore we sought to characterize over, this survival handicap was only partially rescued by whether this occurred during malignant progression. We deter- overexpression of the anti-apoptotic factor BCL2. In comparison, mined the closest gene for each rDNA-interacting fragment based we performed knockdown of the Pol I transcription initiation on the nearest transcription start site (TSS). Analysis of gene factor RRN3 (shRrn3)37, which leads to reduced rRNA transcrip- ontologies (GOs) of protein coding genes associated with con- tion but does not regulate rDNA chromatin remodelling57. Rrn3 stitutive rDNA-interacting regions showed an over- knockdown also resulted in a proliferation disadvantage compared representation of GO terms that comprise zinc finger domain to control cells, but unlike Ubtf knockdown, the cell death in proteins, transmembrane channel proteins and G-protein- response to RRN3 loss was completely rescued by BCL2 coupled receptors (Fig. 3c, Table 1 and Supplementary Data 2). overexpression (Fig. 2e). This is consistent with our previous Inspection of all genes mapping to constitutively interacting data showing that Eμ-Myc lymphoma cells treated with CX-5461, regions revealed an enrichment of non-coding RNA (ncRNA) a highly selective inhibitor of Pol I-mediated transcription36, are genes (e.g., small nucleolar RNAs (snoRNAs), spliceosomal small over 100-fold less sensitive to apoptotic cell death following BCL2 nuclear RNAs (snRNAs)) and genes encoding olfactory and overexpression37. CX-5461 inhibits rDNA transcription by vomeronasal receptors (Table 1, Supplementary Data 2). These preventing recruitment of the transcription competent Pol I results are in keeping with previous studies demonstrating complex (including RRN3) while minimally perturbing UBTF nucleolar enrichment of specific gene classes including certain 4 COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y ARTICLE a b c LM cl2 -B f 1.5 f bt bt 42 P P U U * LM sh sh 42 LMP expression 47S rRNA 1.0 shUbtf α-UBTF 90 kD Active rDNA Inactive rDNA 0.5 2770 bp α-GAPDH 36 kD 0.0 d e LMP shUbtf 8 d post-transduction Fold change %GFP+ ** LMP 1.0 shUbtf * shRrn3 0.5 rRNA/DAPI 0.0 – + – + – + Bcl2 Fig. 2 UBTF loss in lymphoma cells reduces active rDNA and cell survival independent of rDNA transcription. a Eμ-Myc-Bcl2 cells (4242-Bcl2) were stably transfected with empty vector LMP (Eμ-Myc-LMP cells) or shUbtf (Eμ-Myc-shUbtf cells) and Ubtf knockdown verified by Western blot of UBTF protein compared to GAPDH loading control. b Representative Southern blot of psoralen cross-linked genomic rDNA from Eμ-Myc-LMP and Eμ-Myc-shUbtf cells as in Fig. 1c. c Relative pre-rRNA (47SrRNA) expression in Eμ-Myc-LMP or Eμ-Myc-shUbtf cells determined by qRT-PCR. Pre-rRNA levels were determined as in Fig. 1a, at 72 h post transduction. RNA levels were normalized to B2M mRNA and are represented relative to Eμ-Myc-LMP (x ± s.e.m., **p-value < 0.01; n = 3). d rRNA levels in Eμ-Myc-LMP or Eμ-Myc-shUbtf cells determined by rRNA FISH as in Fig. 1b. Scale bars, 10 μm. e GFP competition assay using Eμ- Myc-LMP, Eμ-Myc-shUbtf, and Eμ-Myc-shRrn3 (shRrn3) cells co-cultured with mock-transduced Eμ-Myc cells either in the absence (−) or presence (+) of Bcl2 overexpression. Data are expressed as percent fold change in GFP+ cells 8 days post transduction (x ± s.e.m., *p-value < 0.05; **p-value < 0.01; n = 3) ncRNAs and zinc finger and olfactory receptor encoding UBTF-dependent interactions involve specific gene pathways. genes24,25,27. Interestingly, while we detected 5S rRNA pseudo- To identify interactions that are reliant upon UBTF-dependent genes within the constitutive interaction set, we did not detect any rDNA class switching as opposed to those that are intrinsically a rDNA interactions with the canonical 5S rRNA array on chro- consequence of malignant transformation, we performed 4C-seq mosome 8, consistent with reported Hi-C data28,61. in the stable Ubtf knockdown lymphoma cells (Eμ-Myc-shUbtf) The nucleolar periphery is characterized by a facultative compared to control cells (Eμ-Myc-LMP). Again, the distribution heterochromatic shell and transcriptionally inactive regions30, of the Ubtf-dependent rDNA interactions (Supplementary thus we determined if there is an association between physical Data 5) was spread across all chromosomes (Supplementary interactions with the rDNA and gene expression. We performed Fig. 3b). rDNA interactions with genomic DNA in mouse RNA-seq on wild type, premalignant, and malignant cells and embryonic stem cells (mESCs) have been reported using the determined the average gene expression (reads per kilobase of recently developed split-pool recognition of interactions by tag transcript per million reads mapped [RPKM]) for all genes extension (SPRITE) method62. We interrogated the overlap of (Supplementary Data 3). We compared the density of RPKM mESCs SPRITE rDNA contacts with our 4C-seq rDNA interac- distributions for genes associated with rDNA-interacting (con- tions identified in the Eμ-Myc lymphoma cell lines. This analysis stitutive and differential) and non-interacting regions for all indicated a strong concurrence between the SPRITE-identified stages of malignant progression (Fig. 3d). Genes associated with rDNA hub genomic connections and our 4C rDNA-NAD con- rDNA-interacting regions (Fig. 3d, red lines) in the premalignant nections for both the LMP and shUbtf conditions, with bootstrap and malignant cells, but not wild-type cells, had significantly samples confirming the non-random nature of the intersection lower expression levels than non-interacting genes (Fig. 3d, blue (Supplementary Fig. 3c). We intersected the rDNA interactions lines). This observation is consistent with previous findings24–27 detected in the Ubtf knockdown cells with those from the Eμ-Myc that association with the rDNA correlates with a repressive effect pre-B cells and identified reciprocal interactions that either on gene expression. increase upon premalignant to malignant transition and decrease Global gene expression analysis identified 4908 differentially with UBTF loss (1822 interactions) or decrease upon pre- expressed genes during malignant progression (wild type vs. malignant to malignant transition and increase with UBTF loss malignant cells, FDR ≤ 0.1; −0.5 ≥ logFC ≥ 0.5) (Supplementary (1246 interactions) (Fig. 4a, Supplementary Data 6). These Data 3). We intersected these with the 23,476 rDNA interactions UBTF-dependent rDNA class switch interactions correlated with that changed during progression from premalignant to malig- the changes in rDNA chromatin state that occur during transition nant, coincident with rDNA class switching (Supplementary to malignancy and are reversed by UBTF knockdown (Fig. 4a). Data 4). The majority of differentially expressed genes (3961 out A substantial fraction (28%: 864 of 3068 total interactions) of of 4908 genes) were proximal to differentially interacting regions. the rDNA class switch interacting regions were located within Notably, we observed a robust inverse relationship between 5–20 kb of the associated gene’s TSS (Supplementary Fig. 4a). rDNA interaction and gene expression for rDNA-interacting Considering that changes in rDNA interaction correlated with genes whose expression changes significantly from premalignant altered expression levels of the nearest gene (Fig. 3e), we to malignant (746 genes) (Fig. 3e). These results showed a examined whether these rDNA-interacting fragments are correlation between interaction with the rDNA and transcrip- enriched for putative regulatory elements by intersecting them tional repression of associated genes as cells transition from with regions marked by the enhancer-associated H3K4me1 premalignancy to malignancy, suggesting that rDNA-NAD histone modification in Eμ-Myc tumor cells63. Thirteen percent interactions contribute to Pol II gene regulation during the of loci that differentially interacted with the rDNA during the development of malignancy. premalignant to malignant transition were located in H3K4me1- COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio 5 ARTICLE COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y a rDNA repeat 18S 5.8S 28S IGS 5′ETS ITS1 ITS2 3′ETS Baits b rDNA rDNA rDNA 1 1 1 x x x 25 25 25 50 50 50 150 150 150 100 100 100 0 0 0 0 0 0 75 75 75 125 125 125 0 0 0 125 125 125 100 100 100 150 150 150 175 175 175 75 75 75 50 50 50 19 2 19 2 19 2 25 25 25 0 0 0 25 25 25 0 0 0 50 50 50 50 50 50 75 75 75 25 25 25 0 0 0 10 10 10 18 18 18 0 0 0 5 5 5 12 12 12 75 75 75 0 0 0 15 15 15 50 50 50 5 5 5 25 17 25 17 25 17 0 0 0 0 0 0 25 25 25 50 50 50 17 17 17 75 75 75 3 3 3 50 75 50 75 50 75 25 100 25 100 25 100 0 125 0 125 0 125 150 150 150 75 16 75 16 16 75 0 0 0 50 50 50 25 25 25 25 25 25 50 50 50 4 4 4 0 0 0 75 75 75 100 100 100 100 100 100 75 75 75 15 15 15 125 125 125 50 50 50 150 150 150 25 0 25 0 25 0 0 0 0 25 25 25 125 125 125 50 50 50 100 100 100 5 5 5 75 75 75 14 14 14 75 75 75 100 100 100 50 50 50 125 125 125 25 25 25 150 150 150 0 0 0 0 0 0 25 25 25 100 100 100 50 50 50 75 75 75 75 75 75 13 13 13 50 50 50 6 6 6 100 100 100 25 25 25 12 12 12 0 5 0 5 0 5 0 0 0 0 0 0 10 10 10 25 25 25 75 75 75 50 50 50 50 50 50 12 12 12 75 75 75 25 25 25 7 7 7 10 10 10 0 0 0 0 0 0 12 12 12 5 5 5 100 100 100 15 15 15 0 0 0 0 0 0 75 75 75 11 11 11 50 50 50 25 25 25 25 25 25 50 50 50 0 0 0 75 75 75 125 125 125 100 100 100 100 125 100 125 100 125 8 8 8 75 0 75 0 75 0 50 50 50 25 25 25 25 25 25 50 50 50 100 75 100 75 100 75 0 0 0 10 9 10 9 10 9 wt Prem Mal c d wt 0.25 Not interacting Metal ion binding 0.20 Interacting Transmission of nerve impulse Synaptic transmission Density 0.15 Purine nucleoside binding p < 0.412 Phosphotransferase activity 0.10 Nervous system development p < 0.4486 Substrate-specific channel activity 0.05 Cytoskeleton organization Ion transmembrane transporter activity 0.00 Regulation of cell communication Cytoskeleton 0 2 4 6 8 10 Cation transport RPKM (log2) Actin binding Protein modification process 0.25 Prem Dendrite Not interacting 0 5 10 15 20 0.20 Interacting –log (p value) Density 0.15 p < 0.003 0.10 p < 0.006054 0.05 e Prem vs. mal 0.25 0.00 All DE genes Increased interaction 0 2 4 6 8 10 0.20 Decreased interaction RPKM (log2) 0.15 0.25 Mal Density Not interacting 0.20 Interacting 0.10 p < 1.9 × 10–12 Density 0.15 0.05 p < 0.002 0.10 0.00 0.05 –10 –5 0 5 10 0.00 Fold change (log2) 0 2 4 6 8 10 RPKM (log2) Fig. 3 rDNA repeats engage in long-range chromatin interactions during lymphomagenesis. a Schematic representation of a representative rDNA repeat indicating the 4C-seq bait positions sequences listed in Supplementary Table 1). b Circos plots illustrating the top 10% significantly differential interactions (between a representative 28S rDNA bait and the rest of the genome) in wild type (wt), premalignant (prem) and malignant (mal) Eμ-Myc cells. The depicted differences in rDNA interactions between cell stages is qualitative only. Interactions with each chromosome are highlighted in a different color. c GO analysis of protein coding genes in regions constitutively interacting with rDNA between wild type, premalignant, and malignant cells. Constitutively interacting regions associated with protein coding genes (9382 unique genes, from Supplementary Data 2) were subjected to gene set over-representation analysis using ConsensusPathDB (mouse) to identify enriched GO terms90–92, with p values referring to enrichment of the terms in the described dataset. The genes contributing to these GOs include those encoding zinc finger domain, transmembrane channel and receptor proteins. p-value (x-axis) denotes the significance of the number of interacting genes compared to the total number of genes in the GO classification. d Density histograms of gene expression (RPKM) distributions for genes associated with regions that do interact (red line) and do not interact (blue line) with the rDNA in wild type, premalignant and, malignant cells. Significance between groups was assessed using the nonparametric Kolmogorov-Smirnov (K-S) test. e Density histogram of fold change (log2) in differentially expressed (DE) genes (premalignant vs. malignant) for all genes (including not interacting, black line) and genes associated with significantly increased interactions (red line) and significantly decreased interactions (blue line) with the rDNA in malignant vs. premalignant cells. Significance between increased and decreased interactions was assessed using the K-S test 6 COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y ARTICLE Table 1 Constitutive rDNA-genome interactions between all cell states and associated genes (GO term enrichment, Fig. 3c) Window Mean norm Mean norm counts: Mean norm Gene Gene description (GO term) counts: wild type premalignant counts: malignant name chr8:111250001–111255000 76.8 66.5 80.2 Zfhx3 zinc finger homeobox (metal ion binding) chr7:87180001–87185000 73.7 70.5 55.5 Zfp710 zinc finger protein (metal ion binding) chr2:62840001–62845000 99.7 61.5 72.0 Kcnh7 potassium voltage-gated channel (substrate- specific channel activity) chr11:48850001–48855000 62.3 70.5 77.1 Olfr56 olfactory receptor chr15:95775001–95780000 64.4 95.3 86.4 Ano6 Ca2+ activated chloride channel (ion transmembrane transporter activity) chr8:41560001–41565000 72.7 52.6 80.2 Zdhhc2 zinc finger, DHHC-type containing (protein modification process) chr7:49120001–49125000 159.9 170.7 186.1 Vmn2r58 vomeronasal receptor chr6:136000001–136005000 83.0 92.3 72.0 Grin2b NMDA ionotropic glutamate receptor (transmission of nerve impulse) chr1:87570001–87575000 76.8 85.3 97.7 n-R5s215 nuclear encoded rRNA 5S 215 chr12:92015001–92020000 76.8 63.5 126.5 n-R5s65 nuclear encoded rRNA 5S 65 chr13:20960001–20965000 145.3 127.0 170.7 Gm25605 U1 spliceosomal snRNA (U1) chr15:94530001–94535000 174.4 229.2 232.3 SNORA17 Small nucleolar RNA SNORA17 chr17:88855001–88860000 152.6 158.8 164.5 snoU13 Small nucleolar RNA U13 enriched peaks, and this increased two-fold to 27% (835 of 3,068 ribosome function, RNA processing/metabolism and mitochon- interactions) when we intersected UBTF-dependent rDNA class drial/energy metabolism processes, e.g., Parn, Eif3d, Dhodh, and switch interactions only (Supplementary Fig. 4b). H3K4me1 and Pdp2 (Fig. 4d, Table 3 and Supplementary Data 8). Both RNA H3K4me3 are co-located at transcriptionally poised and/or active processing and mitochondrial function pathways are transcrip- genes64. Therefore, we intersected the H3K4me1 overlapping tionally regulated by oncogenic MYC and dysregulated during rDNA class switch regions with identified H3K4me3 enriched malignancy41,42, and are consistent with GOs that were recently regions63. The majority (66% of 835) of the UBTF-dependent identified as being enriched in rDNA Hi-C contacts detected in rDNA class switch regions were enriched in H3K4me1 only, human lymphoblastoid and leukemic cell lines69. These findings while the remaining 34% were enriched for both H3K4me1 and reinforced the concept that UBTF-dependent changes in rDNA H3K4me3 (Supplementary Fig. 4b). In addition to its roles in interaction are associated with altered gene expression, potentially regulating Pol I transcription and rDNA chromatin status, UBTF through enhancer sequences located within the rDNA-interacting has also been shown to directly regulate a sub-fraction of Pol II fragments and, in a minor subset, through direct UBTF gene transcription65. Thus, we performed UBTF ChIP- regulation. Collectively, these results are consistent with a model sequencing in the malignant progression setting to determine in which reformed rDNA-NAD interactions that occur with whether UBTF is enriched at UBTF-dependent rDNA class UBTF-driven changes in rDNA chromatin during the transition switch interaction regions. As previously reported65,66, the to B-cell malignancy coincide with gene expression changes that majority of UBTF enrichment occurred at the promoter and impact specific functional pathways that promote the malignant transcribed sequences of the rDNA repeat (Fig. 4b). We phenotype. intersected non-rDNA repeat UBTF-enriched peaks from pre- malignant and malignant cells (Supplementary Data 7) with the rDNA-interacting loci undergo nuclear re-localization. Our corresponding UBTF-dependent rDNA interaction regions and results demonstrated that specific loci undergo altered interac- found only ≤1% overlap (Fig. 4c). However, when we intersected tions with rDNA during malignant progression, which is asso- the nearest-located genes to UBTF peaks with those associated ciated with differential gene expression (Fig. 3e). Such changes in with UBTF-dependent rDNA-interacting regions, we observed an interaction frequency are likely the result of physical changes in ~10% overlap of genes in both premalignant and malignant cells the relative positions of rDNA and the locus in question. To test (Fig. 4c, Supplementary Data 7). While this degree of overlap this, we measured the center-to-center distances between repre- might, at least in part, reflect chance occurrence, it nonetheless sentative interacting loci and the rDNA using dual-labeled 3D suggests the potential for UBTF loss to directly impact expression DNA-FISH (3D-FISH)70 (Fig. 5, Supplementary Fig. 5). The of a small subset of rDNA class switch associated genes. rDNA FISH specificity was validated by comparing its localiza- We performed GO analysis on the genes associated with rDNA tion with that of the rDNA-associated pre-rRNA processing class switch interactions that exhibited either increased rDNA protein Fibrillarin (FBL) (Supplementary Fig. 5a). We measured interaction/decreased gene expression or conversely, decreased the distance between rDNA and Ebf1, a gene identified by 4C-seq rDNA interaction/increased gene expression (wild type to to have increased rDNA interaction in malignant vs. pre- malignant) during malignant transformation. B-cell differentia- malignant cells, which is decreased upon Ubtf knockdown in tion and lineage specification ontologies—including genes such as malignant cells. As predicted, the distance between Ebf1 and Ebf1, Pbx1, and Runx1—were significantly enriched in the regions rDNA decreased in malignant compared to premalignant cells that exhibited UBTF-dependent increased interaction and (representative images Fig. 5a, top; quantitation Fig. 5b, left) decreased expression (Fig. 4d, Table 2 and Supplementary while conversely, the distance increased in Eμ-Myc-shUbtf com- Data 8). This observation is consistent with the established role pared to control cells (Eμ-Myc-LMP) (Fig. 5a, bottom; Fig. 5b, of compromised B-lineage-specific gene function in hematologic right). The significant increase in interaction frequency malignancies via impaired differentiation67,68. Furthermore, the (decreased distance) between the Ebf1 locus and rDNA during regions that displayed UBTF-dependent decreased interaction malignant progression was further validated by 3C-qPCR and increased expression were enriched in genes associated with (Fig. 5c). We also confirmed the opposite 4C-seq interaction COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio 7 ARTICLE COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y a 0 10 20 30 40 50 60 70 80 90 100 110 b 18S 5.8S 28S IGS SP rDNA Interaction UP 5′ETS ITS1/2 3′ETS UCE/CPE (0 -50000) Range: 50,000 0 10 20 30 40 50 60 70 80 90 100 110 Mal Interaction (0 -50000) DOWN shUbtf Prem UBTF (0 -50000) wt Significantly (0 -1000) Range: 1000 Input changed interactions 0 10 20 30 40 50 60 70 80 90 100 110 prem vs. mal 5 kb Interaction DOWN Interaction UP shUbtf c d increased interaction Leukocyte activation UBTF-dependent Lymphocyte activation Positive regulation of immune response Hemopoiesis T-cell activation Prem Mal Immune system development Positive regulation of cell differentiation UBTF UBTF 4C-seq 4C-seq B-cell activation ChIP-seq ChIP-seq B-cell differentiation 0 2 4 6 8 10 –log (p value) Peaks 3068 15 1150 3068 30 1739 Ribosome Ribonucleoprotein Ribosomal protein 0.5% 1.0% Structural constituent of ribosome Ribonucleoprotein complex decreased interaction ncRNA processing UBTF-dependent UBTF UBTF Nucleolus 4C-seq 4C-seq ncRNA metabolic process ChIP-seq ChIP-seq rRNA processing rRNA metabolic process 0 5 10 15 Genes 1654 130 997 1654 189 1466 –log (p value) Mitochondrion Mitochondrial part 8% 11% Mitochondrial envelope Mitochondrion inner membrane 0 5 10 –log (p value) Fig. 4 UBTF-dependent rDNA class switch interactions involve specific gene pathways. a Schematic illustrating the filtering workflow for identification of UBTF-dependent rDNA-genome interactions that change from premalignant (prem) to malignant (mal) cells. Interactions that were not UBTF-dependent or not identified in Eμ-Myc-LMP or Eμ-Myc-shUbtf cells were removed. b Integrated Genome Viewer screenshots of UBTF ChIP-seq enrichment at the mouse rDNA repeat (one replicate shown), indicated by the schematic (top); ETS and ITS, external and internal transcribed spacers, respectively; IGS, intergenic spacer; SP, spacer promoter; UCE, upstream control element; and CPE, core promoter element, along with the 18S, 5.8S, and 28S sequences depicted. Below the schematic are ChIP-seq tracks from wild type (wt), premalignant, and malignant cells (set to display the same data range) with the input control track shown at the bottom (data range adjusted for visualisation) and scale depicted underneath. c Intersection of UBTF-dependent rDNA- NAD interactions with UBTF ChIP-seq. Venn diagrams indicating the overlap of UBTF-dependent rDNA-interacting regions (selected as shown in (a)) with UBTF binding sites (top) and nearest associated unique genes for both sets of regions (bottom) in premalignant and malignant cells. d The enriched GO terms were identified using standard DAVID analysis93,94, with p values referring to enrichment of the terms in the described datasets as in Fig. 3c. These enriched ontologies include B-cell differentiation in UBTF-dependent increased interactions or ribosome metabolism and energy metabolism ontologies in UBTF-dependent decreased interactions Table 3 Decreased rDNA class switch interacting regions (window) in malignant compared to premalignant cells Table 2 Increased rDNA class switch interacting regions corresponding to RNA metabolism or energy metabolism (window) in malignant compared to premalignant cells genes (increased expression from wild type to malignant corresponding to B-cell differentiation genes (decreased cells) expression from wild type to malignant cells) Window LogFC LogFC Gene name Window LogFC LogFC Gene interaction expression interaction expression name chr14:55520001–55525000 −1.734 0.851 Pabpn1 chr16:35870001–35875000 2.829 −1.198 Parp14 chr17:56755001–56760000 −0.758 1.771 Rpl36 chr11:33930001–33935000 1.892 −4.911 Lcp2 chr6:71885001–71890000 −0.719 2.377 Polr1a chr16:19925001–19930000 1.250 −1.850 Klhl6 chr15:77855001–77860000 −0.662 1.415 Eif3d chr11:44660001–44665000 0.909 −1.925 Ebf1 chr8:107125001–107130000 −2.623 0.980 Pdp2 chr16: 92765001–92770000 1.514 −1.504 Runx1 chr8:112115001–112120000 −0.657 1.402 Dhodh chr11:44790001–44795000 0.443 −3.657 Il7r chr1:60380001–60385000 −0.118 1.621 Cyp20a1 chr19:34395001–34400000 0.751 −2.090 Fas chr3:97415001–97420000 −3.280 1.570 Chd1l chr3:132755001–132760000 −1.026 1.770 Gstcd 8 COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y ARTICLE pattern, with 3D-FISH showing that the distance between the facilitate UBTF enrichment and rDNA class switching in malig- Eif3d locus and rDNA increased during the premalignant to nant cells. malignant transition and decreased with UBTF knockdown We recently reported a link between UBTF and Pol II gene (Fig. 5d; Supplementary Fig. 5b shows additional 3D-FISH transcription by identifying a role for UBTF at the promoters of assayed loci: Cyp51, Parn, Chd1l; Supplementary Fig. 5c shows highly Pol II-transcribed genes65. Our genome-wide analysis of additional 3C-qPCR assayed loci: Pdp2 and Pbx1). 3D-FISH also UBTF enrichment in Eµ-Myc cells showed minimal overlap confirmed that no significant change in distance was measured (≤1%) between UBTF-bound (non-rDNA) sequences and the between the U1 small spliceosomal RNA gene and the rDNA in rDNA-interacting regions that require UBTF (Fig. 4c). Moreover, either cell condition (Fig. 5e), consistent with its constitutive the rDNA class switch regions do not intersect with previously association with the rDNA. Notably, some of the validated relo- identified UBTF-bound regions such as the histone gene clus- calized genes (e.g. Ebf, Pbx1) were also associated with UBTF ters65. However, when we extend the analysis to the nearest enrichment (Supplementary Data 7). Taken together, these results associated genes, we observe that ~10% of these rDNA- demonstrate that malignant progression coincides with a subset interacting loci are potentially regulated by UBTF. While fur- of specific loci that are ontologically associated with the malig- ther investigation is required, these data suggest that UBTF does nant phenotype undergoing a concurrent program of changes in not regulate Pol II genes directly through rDNA-interacting physical localization, rDNA interaction and gene expression, and regions, but may influence the expression of a small subset of the that these changes directly correlate with UBTF-mediated rDNA UBTF-dependent rDNA-interacting genes either via proximity class switching. while bound to rDNA or enrichment at these genes exclusive of rDNA. Notably, several of the UBTF-enriched Pol II genes (Supplementary Data 7) were validated by 3D-FISH and/or 3C- Discussion qPCR to relocalize with respect to the rDNA during malignant We demonstrate that the rDNA chromatin state is decoupled progression (e.g., Ebf1 and Pbx1). Thus, we propose that a from rRNA transcription during the progression from pre- complex regulatory network including rDNA interactions and, in malignancy to malignancy through a process of pseudo-silenced some instances, direct UBTF regulation of Pol II gene tran- rDNA activation that we term rDNA class switching. This scription, contributes to gene expression reprogramming during observation is consistent with the demonstration that rDNA malignant progression. Nonetheless, silencing Ubtf as a means to chromatin status does not strictly correlate with rRNA tran- reverse rDNA class switching revealed that the majority of scription rates12,13. Indeed, complete loss of UBTF via condi- interactions rely on the altered rDNA state independent of UBTF tional deletion and subsequent eradication of the active rDNA enrichment. Ongoing studies are focused on identifying factors fraction is required before a significant impact on rRNA tran- that directly mediate these dynamic rDNA interactions during scription is observed17. Interestingly, the occurrence of rDNA malignant transformation. The NoRC subunit TIP5/BAZ2A and class switching during MYC-driven lymphomagenesis mirrors the MYC mediate interactions between the rDNA and the nuclear reverse observation as granulocytic cells undergo differentiation matrix74,75 while the chromatin regulatory factor CTCF interacts accompanied by a reduction in MYC, which regulates UBTF with rDNA and is proposed to facilitate UBTF loading onto expression12,13,44. We speculate, given that MYC is dysregulated rDNA, thereby contributing to the activation of silent rDNA in upwards of 70% of human cancers41,71, that rDNA class repeats57,66,72,76,77. The observed enrichment of the enhancer- switching is a more general mechanism underlying malignant associated histone H3K4me1 modification at rDNA class switch- transformation. dependent regions associated with altered gene expression sug- Concomitant with rDNA class switching we observed a marked gests a gene regulatory function for these interactions. Ultimately, repatterning of rDNA-NAD interactions. While chromosome a comprehensive interrogation of these genomic loci via their conformation capture techniques have previously identified genetic disruption and the identification of their interacting genomic regions associated with the rDNA28,62,72, the potential protein components will be required to elucidate their functional dynamic regulation of rDNA-genome contacts and their func- significance. tional significance for disease are not well understood. Here we MYC overexpression drives the Eμ-Myc premalignant pheno- demonstrate that as preneoplastic cells become malignant, type, which we have characterized by increased rRNA tran- changes in the rDNA-interactome correlate with changes in scription of existing active rDNA repeats. We propose that proximal gene expression, an observation supported by previous subsequent epigenetic events during the transition to malignancy reports that established a link between NAD-associated genes and result in increased UBTF association with rDNA and rDNA class transcriptional repression24–27. Importantly, our findings identify switching. The question of whether the captured genomic inter- a subset of altered rDNA-interacting loci and their associated actions occur with the active rDNA repeats only or also the silent gene expression changes that are dependent on UBTF-mediated rDNA remains. Technical challenges preclude our ability to rDNA class switching. Elucidation of this dynamic rDNA inter- directly address this question. However, it is predicted that the action network in the Eµ-Myc cancer model revealed that rDNA malignancy-associated pronounced increase in active rDNA class switching is associated with functional consequences for coincides with altered physical properties of the rDNA chroma- programs of Pol II transcriptional activity that are relevant to tin. Recent studies implicate a role for the biophysical properties malignant cell fitness. of phase separation between sub-nuclear compartments leading rDNA class switching during the malignant transition is to the nucleolar sequestration of certain proteins and driving accompanied by a significant increase in UBTF loading onto the functional heterochromatin formation78–80. Another intriguing rDNA. We and others have shown that UBTF binding correlates consideration is that the epigenetic and/or transcriptional activity with psoralen accessibility12,17,57,73; however, UBTF alone is of the rDNA might impact the transcriptional state of the unlikely sufficient to drive rDNA class switching. Specifically, interacting genomic regions through titration of epigenetic Ubtf knockdown does not significantly alter repressive chromatin modifying or Pol II transcription factors. modifications and enforced Ubtf expression induces only a Our results provide compelling evidence that the formation of modest increase in the proportion of active rDNA12. These specific rDNA interactions contributes to transcriptional repro- observations suggest that prior to association of UBTF, additional gramming of functionally relevant gene pathways. For example, epigenetic modifications occur at the rDNA that are required to genes associated with increased rDNA interaction and decreased COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio 9 ARTICLE COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y a b Ebf1 c Ebf1 Prem Mal 6 3C-qPCR 4 * * ** 400 Distance (μM) 4 % of prem 300 Prem 2 200 2 Mal 100 0 0 0 Ebf1 rDNA em P al f bt LM M Pr U sh LMP shUbtf d Eif3d e U1 6 6 *** 4 4 NS Distance (μM) Distance (μM) **** 4 NS 4 2 2 2 2 Ebf1 rDNA 0 0 0 0 em em P P al al f f bt bt LM LM M M Pr Pr U U sh sh Fig. 5 rDNA class switch-dependent interacting loci undergo nuclear re-localization. a Representative 2D-FISH images taken from the 3D-FISH Z stack images of nuclei depicting a decrease in distance between the early B cell factor 1 (Ebf1) gene (cyan, Alexa 488 signal; cyan arrows) and the rDNA (magenta, Alexa 555 signal; magenta arrows) between premalignant (prem) and malignant (mal) cells. The distance decrease is reversed from LMP (Eμ- Myc-LMP) to shUbtf (Eμ-Myc-shUbtf) cells. Nuclei were counterstained with DAPI (blue). Scale bars, 2 μm. b Boxplots of distances measured from 3D- FISH analysis between the rDNA and the Ebf1 locus in premalignant (red) compared to malignant (gray) cells and LMP (blue) compared to shUbtf cells (cyan). Distances between the closest two foci in each nucleus were measured using the Volocity software. Whiskers correspond to the 10th and 90th percentiles of the data range, with data points outside this range not depicted. Significance was assessed using the Wilcoxon–Mann–Whitney test (*p- value < 0.05, **p-value < 0.01, NS = not significant, n > average 170 distances). c 3C-qPCR validation of the increased interaction frequency between the Ebf1 locus and the rDNA in malignant compared to premalignant cells. The interaction frequency was normalised to a single copy locus (Ercc3 gene) and a multi copy locus (28S rDNA) and is represented as interaction frequency relative to premalignant cells (% of premalignant) (x ± s.e.m., *p-value < 0.05; n = 3). Boxplots of distances measured from 3D-FISH analysis between the rDNA and the d eukaryotic translation initiation factor 3, subunit D (Eif3d) locus between premalignant and malignant cells (increased distance) and reversed with Ubtf knockdown; and the e constitutively interacting U1 spliceosomal snRNA (official MGI gene symbol, Gm25605) locus as described for b (***p-value < 0.001, ****p-value < 0.0001, n > average 170 distances) expression during malignancy are enriched for B-cell lineage and benefit. If so, an increased fraction of active rDNA repeats may differentiation genes, i.e., pre-B cell leukemia homeobox 1 (Pbx1), function as a biomarker for the sensitivity of cancer cells to Early B-cell Factor 1 (Ebf1), and Runx1. The dysregulation of inhibition of rRNA transcription. lineage-associated genes leads to impaired differentiation in hematologic malignancies, such as acute leukemias and some Methods lymphomas67,68,81. While perturbation of differentiation or line- Animal experiments and B-cell purification strategy. All animal experiments age specification is often due to genetic alterations including were performed with approval from the Peter MacCallum Cancer Centre Animal sequence mutation or chromosomal translocation, our findings Experimentation Ethics Committee. Eμ-Myc mice were maintained as heterozygote indicate that dysregulation of B-cell lineage genes can be influ- on a C57BL/6 background and monitored twice a week for lymphoma development. Lymphomas can arise from expansion of either B220+/IgM− or B220+/IgM+ enced by relocalization of genomic regions such as NADs. The clones. Thus, to reduce the potential for variability in our study, we focused our demonstration that a developmentally orchestrated nuclear re- analysis on the IgM-negative B-cell subtypes sorted via flow cytometry from the localization of the Ebf1 locus contributes to Ebf1 transcriptional bone marrow niche from mice with normal B-cells, and from both premalignant activation during B-cell differentiation82 corroborates such a and fully malignant B cells obtained from Eµ-Myc mice, i.e., the B220high- IgMlowIgDlow progenitor B (pro-B)/precursor B (pre-B) cells from premalignant Eµ- mechanism. Our data support a model where structural changes Myc mice and their wild-type littermate controls, and the B220+IgM−IgD− in rDNA chromatin promote gene expression changes that transformed pre-B cells from mice that have developed spontaneous disease influence clonal selection of a malignant cell population. (malignant). Female mice were used for all wild type and premalignant cell col- Together, our data are consistent with the notion that the lections while both female and male mice were harvested for malignant cells due to rDNA repeats play a key regulatory role in cellular function the spontaneous development of the lymphomas. Therefore, this strategy allowed us to perform our experiments largely outside of variables due to differentiation and during disease development, independent of their role in deter- tumor niche. Wild type and premalignant mice were sacrificed at age 4–6 weeks and mining ribosome content and growth rates. We propose that this malignant mice as soon as they established spontaneous lymphoma (age novel role is mediated through the formation and maintenance of 10–20 weeks). dynamic programs of spatial genome organization that regulate distinct patterns of gene expression. Thus, developing tools to Lymphoma cell line studies. For all in vitro experiments, the tissue cultured target rDNA chromatin regulatory factors (such as UBTF) might adapted Eμ-Myc clone #4242 was used and maintained in Anne Kelso DMEM provide a therapeutic window in certain cancers. Previous studies supplemented with 10% fetal bovine serum (FBS), 100 mM L-asparagine (Sigma), on Pol I transcription inhibitors as cancer therapies37,39,50,52 0.5% ß-mercaptoethanol, and penicillin/streptomycin/glutamine (GIBCO). A stable Ubtf knockdown (the validated shUbtf hairpin37 targets both UBTF1 and suggest that there are multiple mechanisms via which rDNA UBTF2 protein isoforms) cell line was generated in a 4242 Eμ-Myc line37 over- transcription and chromatin may be targeted for therapeutic expressing BCL2 to maintain sufficient cell survival with reduced UBTF (also 10 COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y ARTICLE generated without BCL2 for GFP competition assay). Cells were retrovirally temperature for 30 min Cross-links were reversed by the addition of proteinase K transduced with shUbtf-MSCV-LTRmiR30-PIG or empty LMP retroviral vectors (300μg) and incubation at 65 °C (overnight) before treatment with RNase A (also shRrn3-MSCV-LTRmiR30-PIG for GFP competition assay; the shRrn3 (300μg, incubation at 37 °C for 1h). Ligation products were purified by phenol- hairpin was also previously validated37) (Open Biosystems) and transduced chloroform extraction and ethanol precipitation. For 4C-seq the DNA fragments populations selected by sorting for high GFP expression using flow cytometry. For were amplified by PCR using inverse primers in the rDNA loci (Supplementary the GFP-based cell competition assays37, the indicated Eμ-Myc-LMP, Eμ-Myc- Table 1) in the presence of two PCR blockers (Supplementary Table 2), which have shUbtf, or Eμ-Myc-shRrn3 cells were mixed with equal numbers of mock- a C3 spacer at the 5’ and 3’end, located on either side of the specific restriction transduced (no vector) cells, cultured over 8 days and the percentage of GFP+ cells fragments in the rDNA84. The PCR blockers inhibit elongation and thus amplifi- determined by FACS. Eμ-Myc lymphoma cells were treated with 50 nm CX-5461 cation of rDNA sequences85. By using PCR blockers, amplification of adjacent and NaH2PO4 vehicle for 1 h prior to harvesting for psoralen cross-linking and fragments, caused by incomplete digestion or re-ligation of bait and adjacent qRT-PCR analysis. regions, will be avoided (see schematic, Supplementary Fig. 6). All primer sequences are listed in Supplementary Tables 1 (bait primers) and 2 (blockers). PCR products were mixed in equimolar ratios and sequenced on an Illumina 3D DNA-FISH. Cells were fixed with 4% paraformaldehyde and cytospun on poly- HiSeq2500 platform at Peter MacCallum Cancer Centre (150bp, paired-end). l-lysine-coated slides. Prior to hybridization, cells were permeabilized using 0.5% Triton X-100 (20 min at room temperature), incubated in 20% glycerol for at least 1 h and subsequently freeze–thawed for three cycles in liquid nitrogen. Cells were 4C-seq analysis. The sequencing reads were demultiplexed by splitting the reads then treated with 0.1 N HCl (in 0.7% Triton X-100), incubated in RNase ONE according to the 5’ bait sequence followed by computational removal of undigested (Promega) for 30 min, dehydrated through a 70–80–100% ethanol series and (or digested then re-ligated, both of which are minimized by PCR blockers, Sup- incubated in 50% formamide, 2 x SSC for at least 2 h prior to hybridization. Probes plementary Fig. 6) or self-ligated reads. Bait sequences were trimmed using cuta- were hybridized to pretreated slides at 85 °C for 2 min on a heat block and incu- dapt and the trimmed reads were aligned to the mm9 genome assembly bated at 37 °C in a humidified chamber overnight. Post-hybridization washes were (ENSEMBL Mus musculus NCBIM37.67) using bowtie286 with default parameters. performed at 42 °C in 2 x SSC/50% formamide followed by washes in 0.1 x SSC at Reads for which the bait sequences (up to the DpnII site) were not an exact match 60 °C. Slides were incubated with Streptavidin-Alexa fluor 555 (Invitrogen) for were removed and thus interactions originating from 18S and 28S pseudogenes biotinylated probes or anti-digoxigenin Alexa fluor 488 (Invitrogen) for were not considered for downstream analysis. Aligned reads were divided into non- digoxigenin-labeled probes. After washing, nuclei were counterstained using DAPI overlapping 5000-bp windows and the average read number calculated as semi- and slides mounted with Prolong Gold antifade reagent (Invitrogen). Z-sections quantitative measure of interaction frequency. To determine the expected inter- were captured using a Nikon 90i eclipse microscope with equivalent exposure times action frequency relative to TSS a simulated dataset was generated by in silico for all images. An average of more than 170 nuclei were analyzed for each probe DpnII digestion, and analysed using the same parameters as for the real dataset mix using the Volocity 3D image software (PerkinElmer) to measure center-to- (Supplementary Fig. 7). The average distance to TSS was then calculated and center distances. The individual populations (cyan signals for the interacting compared to the distance measured in wild type, premalignant, and malignant regions and magenta signals for the rDNA) were detected separately and only spots cells. While we observed a similar overall distribution of interacting regions within the nucleus (as determined by the DAPI signal) were considered. Fur- compared with our dataset, our findings deviated in that we saw a greater thermore, only cells in which both populations were present in the same nucleus enrichment of regions around approximately 2kb upstream of the closest TSS, but were used for further analysis. By defining the average cell size, individual cells less so immediately at the TSS (Supplementary Fig. 7), corroborating the non- could be discriminated. All distances were measured from the centroids of each random nature of our interactions. population and the minimum distance for each spot selected for statistical analyses, To focus on rDNA interactions independent of bait location within the rDNA thus accounting for allelic differences. repeat, we combined the reads for all bait sequences and removed interactions that Significance was assessed using the Wilcoxon–Mann–Whitney test. The BACs only appeared in less than four of the 11 bait sequences, and then the average read used in the experiments were RP23–225M6 (rDNA), RP23–312-I5 (Ebf1), number was calculated as semi-quantitative measure of interaction frequency. RP23–443-E23 (Chd1l), RP23–16-E24 (Parn), RP23–407-C12 (Cyp51), RP23–449- Subsequently, the counts were normalized to the effective library size, which takes F4 (Eif3d), RP23–190-N22 (Pbx1), RP23–471-F21 (Pdp2), and RP23–448-J8 (U1). the median of the ratio of the count for each gene to the geometric mean of the For the generation of probes, BAC DNA was labeled by nick translation with counts as the scaling factor for each library (see Supplementary Data 1, raw and biotin-16-dUTP or digoxiginin-11-dUTP (Roche) using a nick translation kit normalized to effective library size provided). After removing a small number of (Roche). interactions that displayed very high interaction frequency in only one of the rDNA bait sequences, all normalized interactions were mapped to each chromosome to Psoralen cross-linking assay. Psoralen cross-linking12 was performed in isolated examine their distribution (Supplementary Figure 3a and b). nuclei that were irradiated in the presence of 10 μg/ml 4,5,8’-trimethylpsoralen Differential interactions were determined using the normalized read count data (Sigma–Aldrich). Genomic DNA was isolated and 10 μg was digested with SalI across replicates and applying a FDR cutoff (FDR ≤ 0.1). The R package DEseq was (Promega). Fragmented DNA was separated on a 0.9% agarose gel, transferred used for differential interaction analysis. DEseq is considered more stringent than onto nylon membranes, and immobilized by UV irradiation at 1875 x 100 μJ/cm2 EdgeR, calculating a scaling factor based on mean ratio, which has been shown to using a UV cross-linker (Stratalinker 2400; Agilent Technologies). The membrane perform better in regards to false-positives87. Non-differential interactions (no was hybridized with 32P-labelled rDNA (generated using a nick translation labeling significant change between wild type, premalignant, and malignant; FDR > 0.1), i.e., kit, Roche) and visualized on a PhosphoImager (GE Healthcare). Quantification constitutively associated with the rDNA at all stages, were also identified was performed using ImageQuant (TLv2005.04; GE Healthcare). Original psoralen (Supplementary Data 2). We then focused on differential interactions between blot scans that correspond to quantitated and/or edited blot images shown in premalignant and malignant cells (coincident with rDNA class switching, Figs 1, 2 and Supplementary Fig. 2 are shown in Supplementary Fig. 1a, b and Supplementary Data 4). Only interactions also identified in Eμ-Myc-LMP/Eμ-Myc- Supplementary Fig. 8. A representative depiction of the psoralen blot quantitation shUbtf cells (Supplementary Data 5) were used for further analysis, allowing us to is shown in Supplementary Fig. 1c, d. identify and examine UBTF-dependent interactions. We categorized interactions identified in the Eµ-Myc pre-B cells into those that increase upon transition from premalignant to malignant and decrease in the absence of UBTF, and interactions 4C-library preparation. Chromatin libraries for 4C-seq were prepared as pre- that decrease upon transition from premalignant to malignant and increase in the viously described83,84 and detailed here. Isolated wild type, premalignant and absence of UBTF (Fig. 4a and Supplementary Data 6). malignant pre-B or harvested Eμ-Myc lymphoma cells (LMP and shUbtf), all Interacting regions were visualized using circos plots88, which were generated performed in two biological replicates, were crosslinked with 1% formaldehyde for using the most significantly differential interactions, where interactions identified 10 min, quenched with 125 mM glycine before being pelleted by centrifugation in both replicates with one representative bait in the 28S rDNA (read counts > 10) (320g, 8 min, 4 °C). Cells were suspended in ice cold lysis buffer (10 mM Tris pH were selected and those with the highest interaction frequency, i.e., the top 10% of 8.0, 10 mM NaCl, 0.2% NP-40, protease inhibitor cocktail-EDTA-free (Roche)) and interactions (232 in wild type, 493 in premalignant and 1546 in malignant), were lysed (10 min). Cell nuclei were collected by centrifugation (600g, 5 min, 4 °C), used to generate the circos plots. Comparisons of 4C-seq data and RNA-seq data supernatant removed, and nuclei resuspended in SDS (0.3% final conc., 1 h, 37 °C, for associated regions were performed by displaying the RPKM distribution of with shaking). Triton X-100 (1.8% final conc.) was added to sequester excess SDS genes interacting and genes not interacting with rDNA as a kernel density estimate. and the chromatin solution incubated (1 h, 37 °C, with shaking). Restriction buffer To compare differentially interacting regions with differentially expressed genes, all (NEB) was added (1x final conc.) before the chromatin was digested with DpnII differentially expressed genes between malignant vs. wild type and malignant vs. (100 U, New England Biolabs, 37 °C, o/n, with shaking). DpnII was inactivated premalignant (FDR ≤ 0.1, −0.5 ≥ logFC ≥ 0.5) were grouped into genes associated (SDS [1.3% final conc.], 20 min, 65 °C) and the chromatin solution diluted (10- with increased interacting regions and genes associated with decreased interacting fold) in water. Triton X-100 (1% final conc.) was added and the digested chromatin regions and the kernel density of fold change (log2) in expression estimated. incubated (1 h, 37 °C, with shaking) to sequester the unreacted SDS. As a control The SPRITE method interrogates crosslinked complexes of DNA segments in for possible inter-molecular ligation events, an external DNA fragment from physical proximity by marking all reads in a cluster (nuclear components that had Escherichia coli (E.coli K12 MG1655, sequence: 29555–29765) was added, at the been crosslinked together) with the same barcode62. SPRITE reads were mapped same concentration as a single copy locus, to each reaction. Fragments were ligated (mapq > 10) to NCBIM37.67 using Bowtie2 (v2.3.1) with the default parameters using T4 ligase (20 U; Invitrogen) for 5h at 16 °C before being moved to room (GEO: GSE114242). To find the concordance to our 4C data, we selected all COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio 11 ARTICLE COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y clusters that included rDNA regions 18S and 28S and merged overlapping bins/ experiments (wild type, premalignant, and malignant cells) were prepared and regions. We then compared the 4C interactions (5-kb bins, merged baits, and sequenced using the Illumina HiSeq2500 platform at Peter MacCallum Cancer replicates) to the regions containing SPRITE interactions. Using the R package Centre. The 50-bp paired-end sequencing reads were mapped to the mouse gen- ChIPseeker (v1.8.6)89, we searched 100-kb up- and downstream of the SPRITE ome database mm9 containing one copy of the rDNA (GenBank: BK000964.1) peak (getTagMatrix with windows=SPRITE +/− 100 kb) and plotted the average using bowtie286 and duplicate reads were removed using the picard tool Mark- profile (plotAvgProf, 200-kb region) for the Eµ-Myc lymphoma cell lines (LMP, Duplicates. Reads were visualized using the Integrative Genomics Viewer (IGV)98. shUbtf). For the bootstrap control analyses, SPRITE regions were replaced with the Peaks were called by the R package SPP99 and then used to perform an Irrepro- same number (and size) of random regions from across the genome. The 4C ducible Discovery Rate (IDR) analysis100, which measures the consistency between regions were then mapped to these random regions with the same parameters as ChIP-seq replicates. Common peaks that passed the IDR reproducibility threshold with the SPRITE-identified interactions (window size +/− 100 kb around the of <0.05 were annotated to the mouse genome (mm9) and used for further random bootstrapped regions). intersection/comparison analyses. All analyses were performed using the default GO analysis was performed using ConsensusPathDB90–92 for the large set (9382 options. unique genes) of constitutive (between wild type, premalignant and malignant) interacting regions and DAVID93,94 for the UBTF-dependent differential (between premalignant and malignant) interacting regions. All statistical analyses were RNA-seq. Sequencing libraries were prepared from two biological replicate performed using the R package (R team-Foundation for Statistical Computing experiments (wild type, premalignant, and malignant cells) using the TruSeq RNA 2010). sample preparation kit (Illumina) and sequenced on an Illumina HiSeq2500 platform at Peter MacCallum Cancer Centre (50 bp, SE). The generated 50-bp reads were aligned to the genome using TopHat v2.0.8b with default parameters 3C-qPCR. 3C libraries were prepared as described for 4C-seq but ending with the and the reads counted using HTSeq101. The differential expression between wild extraction of DNA fragments. Quantification of ligated products was performed by type vs. malignant and premalignant vs. malignant was then calculated utilizing the Taqman qPCR95. The primers used for 3C sample amplification are listed in DEseq package102 in R (version 3.0.2) (R team-Foundation for Statistical Com- Supplementary Tables 3 and 4. Primer efficiency and ligation efficiency were puting 2010) and RPKM values calculated (see Supplementary Data 3). Only dif- determined using randomly ligated BAC DNA containing genomic regions of ferentially expressed genes with FDR ≤ 0.1 and −0.5 ≥ logFC ≥ 0.5 were considered interest (BACs listed in 3D DNA-FISH methods), recommended as a control for for further analysis. more reliable quantification in the context of complex mouse/human genomes96. Further, we employed standard protocols for the generation of our control BAC Immunofluorescence. Immunofluorescence followed by rDNA FISH was per- library95. Cross-linking frequencies were normalized to a single copy loci (Ercc3 formed as previously described37,65. Briefly, cells were fixed in a 2% paraf- gene) and multi copy loci (primer in 28S rDNA) (relative interaction frequencies ormaldehyde suspension, spun onto slides and stained with a fibrillarin (FLB) were determined by qPCR compared to standard curves as previously determined). antibody (Abcam, ab5821) followed by the 594-conjugated goat anti-rabbit IgG As internal control, the crosslinking frequency of fragments separated by 2 kb at secondary (Molecular Probe, A21442). Cells were fixed again with 4% paraf- the Ercc3 gene was analysed97. Results were expressed as interaction frequency in ormaldehyde and subjected to rDNA FISH (detailed in Methods, 3D DNA-FISH) malignant cells relative to premalignant cells. and then visualized on an Olympus BX-51 microscope, with images captured using the SPOT Advance image acquisition software (Diagnostic Instruments). FACS analysis and sorting. The following antibodies were used for flow cyto- metric analysis of cells isolated from bone marrow: APC-conjugated anti-B220/ Immunoblotting. Western blot analysis was performed with equal total protein CD-45R (eBioscience), fluorescein isothiocynate (FITC)-conjugated anti-mouse from whole cell lysate samples resolved by SDS-PAGE and transferred onto PVDF IgM, and PE-conjugated anti-mouse IgD (BD Pharmingen). Stained cells were membranes (Millipore). The in-house polyclonal rabbit antisera anti-UBTF1/2 resuspended in buffer containing 2 μmol/L FluoroGold (hydroxystilbamidine, antibody (also used for ChIP-seq) was detected following incubation with anti- Molecular Probes, Invitrogen) before sorting on a FACSAria II flow cytometer (BD rabbit secondary antibody conjugated with horseradish peroxide (HRP) with Biocsiences, San Jose, USA). enhanced chemiluminescence (ECL) reagent (GE Healthcare). RNA-FISH. Cells were fixed with 2% paraformaldehyde and cytospun on glass slides. Slides were washed twice in 10% deionized formamide/2 x SSC, then Data availability All data (4C-seq, RNA-seq and UBTF ChIP-seq) have been deposited in the NCBI hybridized for 3 h at 37 °C in wash buffer containing 10% dextran sulphate, Gene Expression Omnibus (GEO, http://www.ncbi.nlm.nih.gov/geo/) and are 50 μg/ml BSA, 500 ng/ml tRNA, and 10 μg/ml Cy3-conjugated RNA probe (Sigma) available under accession number GSE70226. complementary to part of the rapidly degraded ITS region (thus capturing rDNA transcription rate) of the ribosomal DNA. The slides were then washed twice in wash buffer, then three times in 0.1xSSC at 32 °C, before rinsing in PBS and Received: 4 March 2018 Accepted: 28 December 2018 mounting in VectaShield containing DAPI. rDNA methylation analyses. HpaII digestion and/or ChIP-CHOP assays12 were performed by digesting genomic DNA or DNA isolated following UBTF ChIP with HpaII before qRT-PCR. The relative amount of DNA resistant to HpaII digestion was calculated after normalization to mock-digested DNA. The MeDIP assay was References performed as described12 by incubating 4 μg of sonicated single stranded (i.e., heat- 1. Gibcus, J. H. & Dekker, J. The hierarchy of the 3D genome. Mol. Cell 49, denatured) genomic DNA with 4 μg anti-5mC antibody for 2 h at 4 °C in 0.14 M 773–782 (2013). NaCl, 16.7 mM Tris, pH 8.0, and 0.05% Triton X-100. DNA-antibody complexes 2. Pombo, A. & Dillon, N. Three-dimensional genome architecture: players and were isolated with protein A-Sepharose beads (Millipore)/salmon sperm DNA for mechanisms. Nat. Rev. 16, 245–257 (2015). 2 h at 4 °C and following three washes, the DNA was eluted in 50 mM Tris, pH 8.0, 3. Fraser, J., Williamson, I., Bickmore, W. A. & Dostie, J. An overview of genome 10 mM EDTA, and 0.5% SDS, then purified (including proteinase K digestion) and organization and how we got there: from FISH to Hi-C. Microbiol. Mol. Biol. analyzed by qRT-PCR. The percentage of bound DNA was calculated (subtracting Rev. 79, 347–372 (2015). rabbit serum control) after normalization to 20 ng of input DNA. rDNA amplicon 4. Gorkin, D. U., Leung, D. & Ren, B. The 3D genome in transcriptional primer sequences are listed in Supplementary Table 5. regulation and pluripotency. Cell Stem Cell 14, 762–775 (2014). 5. Dekker, J. & Mirny, L. The 3D genome as moderator of chromosomal RNA extraction and quantitative real-time PCR. RNA was isolated from sorted communication. Cell 164, 1110–1121 (2016). pre-B cells using the NucleoSpin RNA extraction kit (MACHEREY-NAGEL). 6. Krijger, P. H. L. & de Laat, W. Regulation of disease-associated gene cDNA was prepared and qRT-PCR was performed and normalized to B2M mRNA expression in the 3D genome. Nat. Rev. 17, 771–782 (2016). levels44. qPCR analysis of rDNA transcription rates by assaying amplicons to 7. Stults, D. M., Killen, M. W., Pierce, H. H. & Pierce, A. J. Genomic architecture rapidly processed regions of the pre-rRNA was previously validated in Eμ-Myc cells and inheritance of human ribosomal RNA gene clusters. Genome Res. 18, by comparison to [32P] orthophosphate pulse labeling followed by northern ana- 13–18 (2008). lysis, i.e., direct measurement via metabolic labeling of the newly synthesized 8. Gibbons, J. G., Branco, A. T., Godinho, S. A., Yu, S. & Lemos, B. Concerted rRNA37. Primer sequences are listed in Supplementary Table 5. copy number variation balances ribosomal DNA dosage in human and mouse genomes. Proc. Natl Acad. Sci. USA 112, 2485–2490 (2015). 9. McStay, B. & Grummt, I. The epigenetics of rRNA genes: from molecular to ChIP and ChIP-seq. Standard UBTF ChIP was carried out as previously chromosome biology. Annu. Rev. Cell Dev. Biol. 24, 131–157 (2008). described12,13. The polyclonal rabbit antiserum to UBTF1/2 was generated in- 10. Sanij, E. & Hannan, R. D. The role of UBF in regulating the structure and dynamics house and ChIP amplicon primer sequences are listed in Supplementary Table 5. of transcriptionally active rDNA chromatin. Epigenetics 4, 374–382 (2009). ChIP-sequencing of UBTF-immunoprecipitated and input genomic DNA was performed as previously reported65. In brief, libraries from two biological replicate 12 COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y ARTICLE 11. Conconi, A., Widmer, R. M., Koller, T. & Sogo, J. M. Two different chromatin 42. Dang, C. V. MYC, metabolism, cell growth, and tumorigenesis. Cold Spring structures coexist in ribosomal RNA genes throughout the cell cycle. Cell 57, Harb. Perspect. Med. 3, a014217 (2013). 753–761 (1989). 43. Schlosser, I. et al. A role for c-Myc in the regulation of ribosomal RNA 12. Sanij, E. et al. UBF levels determine the number of active ribosomal RNA processing. Nucleic Acids Res. 31, 6148–6156 (2003). genes in mammals. J. Cell Biol. 183, 1259–1274 (2008). 44. Poortinga, G. et al. MAD1 and c-MYC regulate UBF and rDNA transcription 13. Poortinga, G. et al. c-MYC coordinately regulates ribosomal gene chromatin during granulocyte differentiation. EMBO J. 23, 3325–3335 (2004). remodeling and Pol I availability during granulocyte differentiation. Nucleic 45. Grandori, C. et al. c-Myc binds to human ribosomal DNA and stimulates Acids Res. 39, 3267–3281 (2011). transcription of rRNA genes by RNA polymerase I. Nat. Cell Biol. 7, 311–318 14. Grummt, I. & Längst, G. Epigenetic control of RNA polymerase I transcription (2005). in mammalian cells. Biochim. Biophys. Acta 1829, 393–404 (2013). 46. Grewal, S. S., Li, L., Orian, A., Eisenman, R. N. & Edgar, B. A. Myc-dependent 15. Chandrasekhara, C., Mohannath, G., Blevins, T., Pontvianne, F. & Pikaard, C. regulation of ribosomal RNA synthesis during Drosophila development. Nat. S. Chromosome-specific NOR inactivation explains selective rRNA gene Cell Biol. 7, 295–302 (2005). silencing and dosage control in Arabidopsis. Genes Dev. 30, 177–190 (2016). 47. Shiue, C. N., Berkson, R. G. & Wright, A. P. c-Myc induces changes in higher 16. Mais, C., Wright, J. E., Prieto, J. L., Raggett, S. L. & McStay, B. UBF-binding order rDNA structure on stimulation of quiescent cells. Oncogene 28, site arrays form pseudo-NORs and sequester the RNA polymerase I 1833–1842 (2009). transcription machinery. Genes Dev. 19, 50–64 (2005). 48. Poortinga, G., Quinn, L. M. & Hannan, R. D. Targeting RNA polymerase I to 17. Hamdane, N. et al. Conditional inactivation of upstream binding factor reveals treat MYC-driven cancer. Oncogene 34, 403–412 (2015). its epigenetic functions and the existence of a somatic nucleolar precursor 49. Colis, L. et al. DNA intercalator BMH-21 inhibits RNA polymerase I body. PLoS Genet. 10, e1004505 (2014). independent of DNA damage response. Oncotarget 5, 4361–4369 (2014). 18. Boisvert, F.-M., van Koningsbruggen, S., Navascués, J. & Lamond, A. I. The 50. Devlin, J. R. et al. combination therapy targeting ribosome biogenesis and multifunctional nucleolus. Nat. Rev. 8, 574–585 (2007). mRNA translation synergistically extends survival in MYC-driven lymphoma. 19. Paredes, S. & Maggert, K. A. Ribosomal DNA contributes to global chromatin Cancer Discov. 6, 59–70 (2016). regulation. Proc. Natl Acad. Sci. USA 106, 17829–17834 (2009). 51. Quin, J. et al. Inhibition of RNA polymerase I transcription initiation by CX- 20. Paredes, S., Branco, A. T., Hartl, D. L., Maggert, K. A. & Lemos, B. Ribosomal 5461 activates non-canonical ATM/ATR signaling. Oncotarget 7, DNA deletions modulate genome-wide gene expression: ‘rDNA–sensitive’ 49800–49818 (2016). genes and natural variation. PLoS Genet. 7, e1001376 (2011). 52. Hein, N. et al. Inhibition of Pol I transcription treats murine and human AML 21. Zhang, Q., Shalaby, N. A. & Buszczak, M. Changes in rRNA transcription by targeting the leukemia-initiating cell population. Blood 129, 2882–2895 influence proliferation and cell fate within a stem cell lineage. Science 343, (2017). 298–301 (2014). 53. Xu, H. et al. CX-5461 is a DNA G-quadruplex stabilizer with selective lethality 22. Ganley, A. R. D. & Kobayashi, T. Ribosomal DNA and cellular senescence: in BRCA1/2 deficient tumours. Nat. Commun. 8, 1–18 (2017). new evidence supporting the connection between rDNA and aging. Fems. 54. Adams, J. M. et al. The c-myc oncogene driven by immunoglobulin enhancers Yeast. Res. 14, 49–59 (2014). induces lymphoid malignancy in transgenic mice. Nature 318, 533–538 (1985). 23. Diesch, J., Hannan, R. D. & Sanij, E. Perturbations at the ribosomal genes loci are 55. Langdon, W. Y., Harris, A. W., Cory, S. & Adams, J. M. The c-myc oncogene at the centre of cellular dysfunction and human disease. Cell Biosci. 4, 1–11 (2014). perturbs B lymphocyte development in Eµ-myc transgenic mice. Cell 47, 24. Németh, A. et al. Initial genomics of the human nucleolus. PLoS Genet. 6, 11–18 (1986). e1000889 (2010). 56. Harris, A. W. et al. The Eµ-myc transgenic mouse. A model for high-incidence 25. van Koningsbruggen, S. et al. High-resolution whole-genome sequencing spontaneous lymphoma and leukemia of early B cells. J. Exp. Med. 167, reveals that specific chromatin domains from most human chromosomes 353–371 (1988). associate with nucleoli. Mol. Biol. Cell 21, 3735–3748 (2010). 57. Herdman, C. et al. A unique enhancer boundary complex on the mouse 26. Kind, J. et al. Single-cell dynamics of genome-nuclear lamina interactions. Cell ribosomal RNA genes persists after loss of Rrn3 or UBF and the inactivation 153, 178–192 (2013). of RNA polymerase I transcription. PLoS Genet. 13, e1006899 (2017). 27. Pontvianne, F. et al. Identification of nucleolus-associated chromatin domains 58. Santoro, R. & Grummt, I. Molecular mechanisms mediating methylation- reveals a role for the nucleolus in 3D organization of the A. thaliana genome. dependent silencing of ribosomal gene transcription. Mol. Cell 8, 719–725 (2001). Cell Rep. 16, 1574–1587 (2016). 59. Santoro, R., Li, J. & Grummt, I. The nucleolar remodeling complex NoRC 28. Yu, S. & Lemos, B. A portrait of ribosomal DNA contacts with Hi-C reveals 5S mediates heterochromatin formation and silencing of ribosomal gene and 45S rDNA anchoring points in the folded human genome. Genome Biol. transcription. Nat. Genet. 32, 393–396 (2002). Evol. 8, 3545–3558 (2016). 60. Dev, V. G., Tantravahi, R., Miller, D. A. & Miller, O. J. Nucleolus organizers in 29. Huang, S., Deerinck, T. J., Ellisman, M. H. & Spector, D. L. The dynamic mus musculus subspecies and in the RAG mouse cell line. Genetics 86, organization of the perinucleolar compartment in the cell nucleus. J. Cell Biol. 389–398 (1977). 137, 965–974 (1997). 61. Fedoriw, A. M., Starmer, J., Yee, D. & Magnuson, T. Nucleolar association and 30. Németh, A. & Längst, G. Genome organization in and around the nucleolus. transcriptional inhibition through 5S rDNA in mammals. PLoS Genet. 8, Trends Genet. 27, 149–156 (2011). e1002468 (2012). 31. Derenzini, M. et al. Nucleolar size indicates the rapidity of cell proliferation in 62. Quinodoz, S. A. et al. Higher-order inter-chromosomal hubs shape 3D cancer tissues. J. Pathol. 191, 181–186 (2000). genome organization in the nucleus. Cell 174, 744–757 (2018). 32. Hein, N., Hannan, K. M., George, A. J., Sanij, E. & Hannan, R. D. The 63. Sabò, A. et al. Selective transcriptional regulation by Myc in cellular growth nucleolus: an emerging target for cancer therapy. Trends Mol. Med. 19, control and lymphomagenesis. Nature 511, 488–492 (2014). 643–654 (2013). 64. Rye, M., Sætrom, P., Håndstad, T. & Drabløs, F. Clustered ChIP-Seq-defined 33. Montanaro, L., Treré, D. & Derenzini, M. The emerging role of RNA transcription factor binding sites and histone modifications map distinct polymerase I transcription machinery in human malignancy: a clinical classes of regulatory elements. BMC Biol. 9, 80 (2011). perspective. Onco. Targets Ther. 6, 909–916 (2013). 65. Sanij, E. et al. A novel role for the Pol I transcription factor UBTF in 34. Derenzini, M., Montanaro, L. & Treré, D. What the nucleolus says to a tumour maintaining genome stability through the regulation of highly transcribed Pol pathologist. Histopathology 54, 753–762 (2009). II genes. Genome Res. 25, 201–212 (2015). 35. Drygin, D., Rice, W. G. & Grummt, I. The RNA polymerase I transcription 66. Zentner, G. E., Saiakhova, A., Manaenkov, P., Adams, M. D. & Scacheri, P. C. machinery: an emerging target for the treatment of cancer. Annu. Rev. Integrative genomic analysis of human ribosomal DNA. Nucleic Acids Res. 39, Pharmacol. Toxicol. 50, 131–156 (2010). 4949–4960 (2011). 36. Drygin, D. et al. Targeting RNA polymerase I with an oral small molecule CX- 67. Mullighan, C. G. et al. Genome-wide analysis of genetic alterations in acute 5461 inhibits ribosomal RNA synthesis and solid tumor growth. Cancer Res. lymphoblastic leukaemia. Nature 446, 758–764 (2007). 71, 1418–1430 (2011). 68. Iacobucci, I. & Mullighan, C. G. genetic basis of acute lymphoblastic leukemia. 37. Bywater, M. J. et al. Inhibition of RNA polymerase I as a therapeutic strategy J. Clin. Oncol. 35, 975–983 (2017). to promote cancer-specific activation of p53. Cancer Cell 22, 51–65 (2012). 69. Yu, S. & Lemos, B. The long-range interaction map of ribosomal DNA arrays. 38. Bywater, M. J., Pearson, R. B., McArthur, G. A. & Hannan, R. D. PLoS Genet. 14, e1007258 (2018). Dysregulation of the basal RNA polymerase transcription apparatus in cancer. 70. Cremer, M. et al. Multicolor 3D fluorescence in situ hybridization for imaging Nat. Rev. 13, 299–314 (2013). interphase chromosomes. Methods Mol. Biol. 463, 205–239 (2008). 39. Peltonen, K. et al. A targeting modality for destruction of RNA polymerase I 71. Gabay, M., Li, Y. & Felsher, D. W. MYC activation is a hallmark of cancer that possesses anticancer activity. Cancer Cell 25, 77–90 (2014). initiation and maintenance. Cold Spring Harb. Perspect. Med. 4, a014241 (2014). 40. Ji, H. et al. Cell-type independent MYC target genes reveal a primordial 72. Tchurikov, N. A. et al. Hot spots of DNA double-strand breaks and genomic signature involved in biomass accumulation. PLoS One. 6, e26057 (2011). contacts of human rDNA units are involved in epigenetic regulation. J. Mol. 41. Dang, C. V. MYC on the Path to Cancer. Cell 149, 22–35 (2012). Cell Biol. 7, 366–382 (2015). COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio 13 ARTICLE COMMUNICATIONS BIOLOGY | https://doi.org/10.1038/s42003-019-0284-y 73. Hamperl, S. et al. Chromatin states at ribosomal DNA loci. Biochim. Biophys. 98. Thorvaldsdóttir, H., Robinson, J. T. & Mesirov, J. P. Integrative Genomics Acta 1829, 405–417 (2013). Viewer (IGV): high-performance genomics data visualization and exploration. 74. Zillner, K. et al. Large-scale organization of ribosomal DNA chromatin is Brief. Bioinform. 14, 178–192 (2013). regulated by Tip5. Nucleic Acids Res. 41, 5251–5262 (2013). 99. Kharchenko, P. V., Tolstorukov, M. Y. & Park, P. J. Design and analysis of 75. Shiue, C.-N., Nematollahi-Mahani, A. & Wright, A. P. H. Myc-induced ChIP-seq experiments for DNA-binding proteins. Nat. Biotechnol. 26, anchorage of the rDNA IGS region to nucleolar matrix modulates growth- 1351–1359 (2008). stimulated changes in higher-order rDNA architecture. Nucleic Acids Res. 42, 100. Li, Q., Brown, J. B., Huang, H. & Bickel, P. J. Measuring reproducibility of 5505–5517 (2014). high-throughput experiments. Ann. Appl. Stat. 5, 1752–1779 (2011). 76. van de Nobelen, S. et al. CTCF regulates the local epigenetic state of ribosomal 101. Anders, S., Pyl, P. T. & Huber, W. HTSeq--a Python framework to work with DNA repeats. Epigenetics Chromatin 3, 19 (2010). high-throughput sequencing data. Bioinformatics 31, 166–169 (2015). 77. Huang, K. et al. Ribosomal RNA gene transcription mediated by the master 102. Anders, S. & Huber, W. Differential expression analysis for sequence count genome regulator protein CCCTC-binding factor (CTCF) is negatively data. Genome Biol. 11, R106 (2010). regulated by the condensin complex. J. Biol. Chem. 288, 26067–26077 (2013). 78. Mitrea, D. M. et al. Nucleophosmin integrates within the nucleolus via multi- modal interactions with proteins displaying R-rich linear motifs and rRNA. Acknowledgements eLife 5, 1–33 (2016). We thank Richard Tothill for helpful discussions. We thank Tim Semple and Gisela Mir 79. Feric, M. et al. Coexisting liquid phases underlie nucleolar subcompartments. Arnau and the Peter Mac Molecular Genomics Core Facility for their contributions to Cell 165, 1686–1697 (2016). this study. We also thank Kerry Ardley and the Peter Mac Animal Core Facility for their 80. Strom, A. R. et al. Phase separation drives heterochromatin domain formation. assistance with the mice used in this study. We thank Jason Ellul and Maria Doyle from Nature 547, 241–245 (2017). the Peter Mac Research Computing Facility for their assistance with the bioinformatics 81. Garraway, L. A. & Sellers, W. R. Lineage dependency and lineage-survival associated with this work. This work was supported by project and program grants (R.D. oncogenes in human cancer. Nat. Rev. 6, 593–602 (2006). H., G.A.M., and R.B.P.) from the National Health and Medical Research Council 82. Lin, Y. C. et al. Global changes in the nuclear positioning of genes and intra- (NHMRC) of Australia and also Cancer Council of Victoria and Leukaemia Foundation and interdomain genomic interactions that orchestrate B cell fate. Nat. of Australia grants in aid (G.P., G.A.M., and R.D.H.). The researchers were funded by the Immunol. 13, 1196–1204 (2012). following: NHMRC fellowships (R.D.H., R.B.P., and G.A.M.) and a Sir Edward Weary 83. Rodley, C. D. M., Pai, D. A., Mills, T. A., Engelke, D. R. & O'Sullivan, J. M. Dunlop fellowship (G.A.M.). tRNA gene identity affects nuclear positioning. PLoS ONE. 6, e29267 (2011). Author contributions 84. Doynova, M. D. et al. Interactions between mitochondrial and nuclear DNA J.D., J.M.O., R.D.H. and G.P. were responsible for the overall concept and design of in mammalian cells are non-random. Mitochondrion 30, 187–196 (2016). experiments. J.D., M.J.B., E.S., D.P.C., W.S., N.B., Jb.S., J.S., J.M.O., R.D.H. and G.P. were 85. Vestheim, H. & Jarman, S. N. Blocking primers to enhance PCR amplification responsible for the collection and assembly of data. J.D., M.J.B., E.S., D.P.C., W.S., N.H., M. of rare sequences in mixed samples—a case study on prey DNA in Antarctic E., R.B.P., G.A.M., A.R.D.G., J.M.O., R.D.H. and G.P. were involved in data analysis and krill stomachs. Front. Zool. 5, 12 (2008). interpretation. J.M.O., R.D.H., and G.P. supervised the research. J.D., M.J.B., E.S., R.B.P., G. 86. Langmead, B. & Salzberg, S. L. Fast gapped-read alignment with Bowtie 2. Nat. A.M., A.R.D.G., J.M.O., R.D.H. and G.P. were involved in writing the manuscript. Methods 9, 357–359 (2012). 87. Dillies, M.-A. et al. A comprehensive evaluation of normalization methods for Illumina high-throughput RNA sequencing data analysis. Brief. Bioinform. 14, Additional information 671–683 (2013). Supplementary information accompanies this paper at https://doi.org/10.1038/s42003- 88. Krzywinski, M. et al. Circos: an information aesthetic for comparative 019-0284-y. genomics. Genome Res. 19, 1639–1645 (2009). 89. Yu, G., Wang, L.-G. & He, Q.-Y. ChIPseeker: an R/Bioconductor package for Competing interests: The authors declare no competing interests. ChIP peak annotation, comparison and visualization. Bioinformatics 31, 2382–2383 (2015). Reprints and permission information is available online at http://npg.nature.com/ 90. Kamburov, A., Wierling, C., Lehrach, H. & Herwig, R. ConsensusPathDB--a reprintsandpermissions/ database for integrating human functional interaction networks. Nucleic Acids Res. 37, D623–D628 (2009). Publisher’s note: Springer Nature remains neutral with regard to jurisdictional claims in 91. Kamburov, A. et al. ConsensusPathDB: toward a more complete picture of cell published maps and institutional affiliations. biology. Nucleic Acids Res. 39, D712–D717 (2011). 92. Herwig, R., Hardt, C., Lienhard, M. & Kamburov, A. Analyzing and interpreting genome data at the network level with ConsensusPathDB. Nat. Open Access This article is licensed under a Creative Commons Protoc. 11, 1889–1907 (2016). Attribution 4.0 International License, which permits use, sharing, 93. Huang, D. W., Sherman, B. T. & Lempicki, R. A. Bioinformatics enrichment adaptation, distribution and reproduction in any medium or format, as long as you give tools: paths toward the comprehensive functional analysis of large gene lists. appropriate credit to the original author(s) and the source, provide a link to the Creative Nucleic Acids Res. 37, 1–13 (2009). Commons license, and indicate if changes were made. The images or other third party 94. Huang, D. W., Sherman, B. T. & Lempicki, R. A. Systematic and integrative material in this article are included in the article’s Creative Commons license, unless analysis of large gene lists using DAVID bioinformatics resources. Nat. Protoc. indicated otherwise in a credit line to the material. If material is not included in the 4, 44–57 (2009). article’s Creative Commons license and your intended use is not permitted by statutory 95. Hagège, H. et al. Quantitative analysis of chromosome conformation capture regulation or exceeds the permitted use, you will need to obtain permission directly from assays (3C-qPCR). Nat. Protoc. 2, 1722–1733 (2007). the copyright holder. To view a copy of this license, visit http://creativecommons.org/ 96. Dekker, J. The three ‘C’ s of chromosome conformation capture: controls, licenses/by/4.0/. controls, controls. Nat. Methods 3, 17–21 (2006). 97. Splinter, E. et al. CTCF mediates long-range chromatin looping and local histone modification in the beta-globin locus. Genes Dev. 20, 2349–2354 (2006). © Crown 2019 14 COMMUNICATIONS BIOLOGY | (2019)2:39 | https://doi.org/10.1038/s42003-019-0284-y | www.nature.com/commsbio