p53 status determines the epigenetic response to demethylating agents azacitidine and decitabine
Emma Langdale Hands, Arndt Wallmann, Gabrielle Oxley, Sophie Storrar, Rochelle D’Souza, Mathew Van de Pette
TL;DR
This study shows that the p53 protein controls how cells respond to demethylating drugs azacitidine and decitabine, with p53 status acting as a potential biomarker for treatment choice.
Contribution
The study identifies p53 as a key determinant of drug-specific epigenetic responses and introduces p53 status as a novel biomarker for demethylating agent selection.
Findings
Azacitidine activates p53, leading to R-loop accumulation in CpG islands of p53 target genes.
Decitabine induces DNA damage in the absence of p53, while p53 removal effects on the epigenome are reversible.
p53 status determines divergent epigenetic responses to azacitidine and decitabine.
Abstract
5’-Azacitidine (Aza) and 5-Aza-2’-deoxycytidine (Dac) are widely used demethylating drugs that directly integrate into nucleic acids. They are frequently used interchangeably, surprisingly as their selectivity is unique from the other, with no predictors of response or clinical biomarkers to indicate drug preference. Using these drugs to induce demethylation, we combine DRIPc-Seq, Immunostaining, RNA-Seq and Mass spectrometry to uncover unique cellular responses. Activation of p53, exclusively by Aza, sustains accumulation of R-loops in CpG islands of p53 target genes. This effect is abolished by the removal of p53, compounded by destabilisation of heterochromatin marks. Dac treatment induces global chromatin modification, sustaining DNA damage, which is heightened in the absence of p53. Rescue experiments reverse the changes observed in the epigenome, demonstrating a direct role for…
Genes, proteins, chemicals, diseases, species, mutations and cell lines named across the full text — each resolved to its canonical identifier and authoritative record.
Click any figure to enlarge with its caption.
Figure 10
Figure 11
Figure 1
Figure 2
Figure 3
Figure 4
Figure 5
Figure 6
Figure 7
Figure 8
Figure 9
Figure 12
Figure 13- —http://dx.doi.org/10.13039/501100000265UKRI | Medical Research Council (MRC)
- —UK Health Security Agency
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Taxonomy
TopicsEpigenetics and DNA Methylation · Histone Deacetylase Inhibitors Research · Cancer-related Molecular Pathways
Introduction
The scope of our ability to target epigenetic modifications through pharmaceutical agents has expanded significantly in recent years, including Polycomb Repressor Complex II inhibitors, histone deacetylase inhibitors and drugs that specifically modify RNA methylation status (Burgess et al, 2021; Yankova et al, 2021). Many of these show great promise; however, it is only DNA methyltransferase inhibitors that have been progressed into the clinic to date. 5’-Azacitidine (Aza) and 5-Aza-2’-deoxycytidine (Dac) have been used for the treatment of acute myeloid leukaemia (AML) and myelodysplastic syndrome (MDS) for a number of years (Jabbour et al, 2008; Saiz-Rodríguez et al, 2021; Khan et al, 2012) and are currently being tested for their efficacy in solid tumour treatment, often in combination with other therapies (Luke et al, 2023; Li et al, 2015). Their mode of action differs, by virtue of the ability of Aza to integrate into RNA and DNA, inhibiting m5C methylation in both (Chen et al, 2025), whereas Dac exclusively integrates into DNA (Christman, 2002; Li et al, 1970). Within the context of DNA methylation, both drugs principally act through the formation of DNMT1-DNA adducts, consequently resulting in a depletion of the available DNMT1 pool with which to maintain DNA methylation levels (Jüttermann et al, 1994; Maslov et al, 2012; Stresemann and Lyko, 2008), with off-target modifications to the broader epigenome believed to arise from the knock-on effects of this shared mechanism (Cedar and Bergman, 2009; Seelan et al, 2018; Komashko and Farnham, 2010). While their ability to inhibit tumorigenesis is well-known, the cellular pathways through which this is achieved remain opaque, beyond their ability to reduce the hypermethylation that is generally associated with the cancer state (El Fakih et al, 2018). This is troubling as these drugs are frequently used interchangeably in the clinic, and there is no known predictive marker that would indicate which drug to prescribe to ensure maximal efficacy (Ma and Ge, 2021; Kuendgen et al, 2018).
Studies have indicated that in addition to a global removal of DNA methylation, Aza and Dac also induce broader changes to the epigenome following treatment (Manzoni et al, 2016; Komashko and Farnham, 2010; Malkaram et al, 2022; Tobiasson et al, 2017). Most recently, treatment of a cancer cell line with Dac induced an abundance of R-loops (Panatta et al, 2022). R-loops are transient and unstable tripartite complexes that form primarily when a nascent RNA hybridises with one of the DNA strands, ejecting the other. This structure carries with it risks of DNA damage through double-strand breaks, as transcriptional stalling caused by the R-loop can lead to a collision of the replication and transcription complexes (Crossley et al, 2019; Aguilera and García-Muse, 2012), thereby making any cancer treatment that induces an excess of R-loops potentially counterproductive. However, the relative abundance and distribution of R-loops in otherwise healthy cells, observed through genome-wide sequencing techniques (Sanz and Chédin, 2019; Ginno et al, 2012), demonstrates that with the risk must come the reward of such a structure. Indeed, multiple studies in recent years have revealed physiological roles of R-loops across different genes and cell types (Petermann et al, 2022; Niehrs and Luke, 2020; Yu et al, 2003; Al-Hadid and Yang, 2016; Daniels and Lieber, 1995).
The physiological function of the bulk of R-loops are thought to partly converge on transcriptional regulation, whether that be through modifying the total dosage of a gene product, or by altering the splicing pattern of the transcribed RNAs; the manner in which an R-loop achieves this functional consequence is believed to be contingent on its location within any given gene locus (Ginno et al, 2012, 2013; Skourti-Stathaki et al, 2014). In general, so-called promoter R-loops are considered to be transcriptionally activating (Niehrs and Luke, 2020), while terminator R-loops are shown to induce repressive chromatin marks and are associated with DNA damage-inducing events (Skourti-Stathaki et al, 2014; Hatchi et al, 2015). In the literature, R-loops are frequently characterised as being either physiological, and required for normal cell function, or pathological, forming aberrantly and with the potential to induce DNA damage. However, there is currently little consensus as to what constitutes, or how to distinguish, physiological or pathological R-loops. All R-loops will need to be removed however and to achieve this, cellular machinery primarily mediated by RNase H (San Martin Alonso and Noordermeer, 2021; Cerritelli and Crouch, 2009), resolves them, resulting in a highly dynamic mark.
The variation on the functional consequence of an R-loop depending on where it forms, and the overriding modification of transcription, independent of changes to genetic sequence, implicate R-loops as a component of the overall epigenetic architecture of a cell. Further to this, emerging data has highlighted that R-loop formation and removal is responsive to, and potentially responsible for, changes to distinct epigenetic marks, including DNA methylation and histone modifications (Sanz et al, 2016; Nadel et al, 2015; Skourti-Stathaki et al, 2014). In particular, R-loops are especially enriched in CpG islands that are unmethylated, with recent data suggesting that changes in methylation will consequently modify the global R-loop profile (Hartono et al, 2015; Ginno et al, 2012). With this in mind, we asked whether depletion of the methylome, through either Aza or Dac, would consequently involve the genome-wide accumulation of R-loops in newly unmethylated CpG islands, and whether this results in elevated DNA damage and cell stress. Recent advances in methodology now allows for strand-specific, single-base resolution mapping of R-loops genome-wide, through the next-generation sequencing (NGS) of the offending RNA strand (DNA-RNA immunoprecipitation followed by cDNA conversion and sequencing, DRIPc-Seq (Sanz and Chédin, 2019)). Here, we combine this approach with epigenetic analysis to uncover divergent mechanisms of cellular action for Aza and Dac.
Results
Aza and Dac treatment leads to epigenetic re-mapping at sub-toxic concentrations
Aza and Dac are well known to be substantially toxic to mammalian cells, and their anti-cancer properties are thought to occur through the demethylation and activation of tumour suppressor genes (Frontiers Production Office, 2023). While demethylation was a requirement for our model, we sought to identify a clinically relevant dosing regimen for each drug in HEK293 cells where acute toxicity was minimalized (Fig. 1A,B). Both drugs were found at 1.5 µM to satisfy this requirement, with no significant changes in viability or early apoptosis and minimal cytotoxicity after 48 h of treatment. Whole-genome bisulfite sequencing (WGBS-Seq) demonstrated that at this dose, substantial genome-wide demethylation occurred following exposure to either drug, with broadly overlapping regions of demethylation. Somewhat surprisingly, given that Aza is believed to preferentially integrate in RNA over DNA (Li et al, 1970; Christman, 2002), this demethylation of the DNA was found to be greater and more widespread in Aza samples than Dac, with 30.6% of significantly demethylated tiles unique to Aza-treated cells, compared to only 10.6% in Dac-treated cells (Fig. 1C). Analysis of the genomic location of all hypomethylated tiles revealed an enrichment in coding regions of the genome, that was predominantly shared by both drugs (Fig. 1D; Fig. EV1A–C).Figure 1. Hypomethylating conditions are sufficient to induce changes in transcription and global R-loop abundance.(A) Experimental schematic of drug treatment in HEK293 cells. Cells were treated for 48 h with either Aza, Dac or solvent control with a refresh after 24 h. (B) ApoTox Toxicity assay in HEK293 cells treated with a dose range of Aza or Dac to measure cell viability, early apoptosis and cytotoxicity levels. n = 3 biological repeats. Error bars are 士SD. An Ordinary one-way ANOVA with Šídák’s multiple comparisons test was performed between all doses vs untreated. Top panels left to right: ****P = < 0.0001, ****P = < 0.0001, ****P = < 0.0001, *P = 0.0235, ****P = < 0.0001, ****P = < 0.0001, ****P = < 0.0001, ****P = < 0.0001, ****P = < 0.0001, ****P = < 0.0001, ****P = < 0.0001, ****P = < 0.0001. Bottom Panels left to right: **P = 0.0012, **P = 0.0018, ****P = < 0.0001, ****P = < 0.0001, ****P = < 0.0001, ****P = < 0.0001, and ****P = < 0.0001. (C) Violin plots of CpG methylation of all available positions from WGBS of Aza- or Dac-treated HEK293 cells. The p values were calculated by Wilcoxon test. ***P < 2.2e-16. n = 2 biological repeats, triplicate repeats per biological repeat. Venn diagram of hypomethylated DMRs based on 1000 bp tiling of the genome. (D) Feature distribution of all tiles across the genome and shared and treatment-specific DMRs as shown in (C). (E) Volcano-plot of RNA-seq analysis of Aza- and Dac-treated HEK293 cells. Significantly differentially expressed genes (DEGs) (Aza: 111; Dac: 13) are highlighted in the respective colour. n = 3 biological repeats. (F) Venn diagram of upregulated DEGs in Aza- and Dac-treated cells. Significantly altered genes are listed in Dataset EV1. (G) Proportion of hypomethylated promoters in all genes and the upregulated DEGs under Aza treatment for Aza and Dac treatment, proportion of non-hypomethylated promoters is in grey. (H) Representative confocal images using S9.6 antibody to stain Aza- and Dac-treated HEK293 cells for R-loops (red) and DAPI as a counterstain for DNA (grey). Scale bars 10 μm. (I) Image quantification of S9.6 staining plotted as median S9.6 nuclear fluorescence and average S9.6 foci number per nuclei. n = 3 biological repeats. Statistical testing was calculated by ordinary one-way ANOVA with Dunnett’s multiple comparisons test. (R-loops) ****P = < 0.0001, (Foci) **P = 0.0015. Error bars are 士SD. Source data are available online for this figure.
To determine if this broad demethylation induced changes in transcription, RNA-sequencing was performed. Ninety-eight genes were uniquely and significantly upregulated following Aza treatment, in addition to 13 that were also increased following Dac treatment. These 13 genes represented the total of all significantly modified genes in Dac treatments (Fig. 1E,F; Dataset EV1). The fact that expression was predominantly elevated genome-wide following the treatment of Aza and Dac corresponds with the known function of methylation, and for the 98 genes that were found to be upregulated following Aza treatment, their promoters were more demethylated when compared to all promoters across the genome (62.2 vs 24.1%). Analysis of these same 98 promoter regions in the Dac-treated samples determined that a smaller proportion of these promoters became hypomethylated following that treatment (29.7 vs 15.7%) (Fig. 1G). This indicates that beyond the broad demethylation observed by both drugs, the change in methylation status of selected promoters may have a direct impact on the transcriptional profiles of the treated cells.
In light of the known enrichment of R-loops in hypomethylated CpG islands (Ginno et al, 2012; Hartono et al, 2015), and following our confirmation that Aza and Dac caused global demethylation in the absence of substantial acute toxicity, we asked whether treatment with Aza and Dac induced an accumulation of global R-loops that corresponded to the changes in DNA methylation or transcription. Immunostaining for R-loops (S9.6 clone, red, Fig. 1H) demonstrated a differential response to the two drugs, with the total fluorescence and the number of nuclear foci elevated following Dac treatment, in accordance with published work (Panatta et al, 2022), while in contrast a reduction in total fluorescence was observed in Aza-treated cells (Fig. 1I). This increase in R-loop staining in Dac was accompanied by elevated markers of DNA damage, including micronuclei formation (Fig. EV1D,E) and mitotic spindle pole errors (Fig. EV1F,G), which was not observed in Aza-treated samples. R-loop formation has been observed to coincide with a decrease in the repressive histone marker H3K27me3 by inhibiting the binding of the Polycomb Repressive Complex 2 (PRC2), the canonical writer of H3K27 methylation (Chen et al, 2015). Immunostaining for H3K27me3 in treated cells showed a significant decrease in nuclear staining intensity in Aza-treated cells, whereas there was an increase in Dac-treated cells (Fig. EV1H). In contrast, global H3K9me3 signal was significantly increased in Aza-treated cells, while no significant changes were detected in Dac-treated cells (Fig. EV1I).
Global R-loops are re-mapped following Aza and Dac treatment
Due to the substantially different staining patterns for R-loops in either Aza or Dac treatments, we chose to perform DNA:RNA immunoprecipitation followed by cDNA conversion and sequencing (DRIPc-Seq) (Sanz et al, 2016). This approach enables single-base and strand-specific resolution of R-loops, allowing a base resolution genome-wide picture of the changes that had occurred to be generated (Fig. 2A). Intra-experimental technical variation has been reported for S9.6-based R-loop detection (Fig. EV2A) (Sanz et al, 2016; Smolka et al, 2021), and so to account for this, we employed strict controls as outlined in the methods section. The size and stranded distribution of R-loops were only modestly changed by either treatment (Fig. EV2B,C). However, the genomic location of Aza or Dac R-loops were substantially disturbed, with the majority of these changes found within coding regions (Fig. 2B). These changes though were unique from one another, and while an accumulation of R-loops was found to occur in promoters or at the transcriptional start site (TSS) following Aza treatment, a notable decline in intronic R-loops was also observed. In contrast, promoter R-loops became depleted in Dac samples, while there was a substantial accumulation of R-loops near the transcriptional termination site (TTS) (Fig. 2B,C). In fact, co-accumulation or co-depletion of R-loop peaks in Aza and Dac treatment was almost completely absent, with a much stronger likelihood of an inverse response occurring (Fig. EV2D).Figure 2. Azacitidine and Decitabine treatment uniquely alter the R-loop landscape.(A–E) DRIPc-seq analysis. n = 2 biological repeats, triplicate samples were analyzed per biological replicate. (A) Circos plot showing the unstranded R-loop signals across chromosomes in 1000-kb bins for Control (Ctrl), Azacitidine (Aza), Decitabine (Dac) and RNase H-treated HEK293 cells. Schematic cartoon summarising the DRIPc-seq approach through the use of the S9.6 antibody to recognise R-loops, which can then be immunoprecipitated for sequencing. (B) Feature distribution for significantly increased (gained) and decreased (lost) R-loop signals in cells treated with Aza or Dac, compared to solvent control. (C) Peak count frequency for Aza and Dac gained R-loop peaks across the gene body. (D) Top-scoring pathway enrichment terms of genes with gained R-loop peaks from Aza-treated cells. P values were generated through the Benjamini–Hochberg procedure. (E) DRIPc-seq representative R-loop peak tracks of p53-related genes after Ctrl treatment and Aza treatment (Watson/forward-strand R-loop (+), Crick/reverse-strand R-loop (−)). The gained peaks upon Aza treatment are highlighted with grey boxes. (F,** G**) Western Blot of p53. (F) Quantification of signal intensity of Western Blot bands from phospho-p53 staining, normalised to the signal intensity of the loading control β-Actin and plotted as fold change over Ctrl. For the Ctrl-treated group, fold change was calculated against the average of all Ctrl samples vs 1. n = 8 biological repeats, p values were calculated by one-way ANOVA with Dunnett’s multiple comparisons test, ****P = 0.0001, data were presented as mean 土 SD. (G) Representative Western blot of total p53 and β-Actin (bottom two image panels), and p53 with serine 15 phosphorylation (Phospho-p53 Ser15) and β-Actin (top two image panels) in HEK293 cells treated with 1.5 μM Aza, Dac or solvent control (Ctrl) or UV (positive control for p53 activation). The β-actin images were taken from the same membrane as their corresponding p53 images; each membrane was imaged in both the 488-fluorescence channel (β-Actin) and 568 fluorescence channel (phospho-p53/Total p53). Source data are available online for this figure.
As has been already stated, R-loops are known to preferentially inhabit unmethylated CpG islands (Ginno et al, 2012); (Hartono et al, 2015). In either Aza or Dac treatment, however, the relationship between newly demethylated regions and an accumulation of R-loops was not immediately obvious at a genome-wide level. In fact, for Aza treatment, a depletion of an R-loop was more likely to occur in a newly demethylated region, rather than an enrichment, while Dac-induced demethylated regions were as likely to form R-loops as they were to lose them (Fig. EV2E). To overcome this complexity, we restricted our analysis to promoters for Aza treatment, as our earlier analysis had indicated a link between demethylation of promoter regions, and an increase in transcription of the associated gene, that was unique to Aza (Fig. 1G). Pathway enrichment for genes with a newly formed promoter R-loop, regardless of whether their expression was significantly modified, highlighted genes associated with chromatin architecture and the cell cycle, potentially indicating a global epigenetic response to the initial exposure. To our surprise, however, the top enrichment was for p53 target genes (Fig. 2D,E). p53 is a protein which is frequently named the ‘guardian of the genome’ because of its various roles in DNA damage repair, cell cycle arrest and apoptosis, in response to a variety of cell stressors, thus making it a critical tumour suppressor and by extension the most frequently mutated gene in human cancers. DNA damage as a result of increased R-loop formation has previously been linked to p53 activation (Lindström et al, 2022), and p53-deficient cells have been shown to be dependent on the DHX9 DexH-box helicase, an RNA helicase that functions in the removal of R-loops (Lee and Pelletier, 2016). Despite this, very little is known about the role of p53 directly in the formation of R-loops.
No p53-target genes were found to be differentially expressed in the RNA-Seq analysis (Figs. 1E and EV2F), while p53 activation by these drugs has not yet been clearly demonstrated. Western blotting for Phospho-p53 Ser15, a mark that represents activated p53, was performed, and determined to be activated following Aza treatment (Fig. 2F,G). In contrast, such activation was not observed following Dac treatment, and no associated changes in p53-target gene promoter R-loops were detected (Figs. 2F and EV2G). Promoter R-loops are generally considered to function through the promotion of transcription; however, our analysis had not determined transcriptional activation for the p53-target genes in question. By definition, a cis R-loop requires transcription to occur; therefore, we reasoned that it was possible our timing was missing an earlier increase in the expression of these p53-target genes. To determine if this was the cause, we performed time-course analysis of the expression of three known p53-target genes, which were enriched for R-loops in the Aza-treated samples (DDIT4, CDKN1A, and* GADD45A*) and one p53-target gene, which was not enriched for R-loops (MAGED2), as well as expression of the p53 gene (TP53) (Fig. EV2H). While modest increases in expression were detected for certain p53-target genes following Aza treatment (DDIT4 and* GADD45A*) at timepoints which would have precluded detection through our earlier analysis, a clear concerted increase in other genes (MAGED2 and CDKN1A) was not observed. Following Dac treatment, there was only a significant increase in MAGED2 expression after 48 h. Expression of TP53 was shown to significantly increase following Aza treatment, but not Dac treatment.
Knockout of p53 affects R-loop formation upon treatment
To determine if p53 had a direct role in regulating R-loop formation following Aza or Dac treatment, R-loop immunofluorescence staining was performed in HEK293 p53 KO cells. Aza and Dac-treated p53 KO cells displayed significantly increased total staining and foci number compared to untreated p53 KO cells (Fig. 3A,B). The total nuclear staining or foci formation was also significantly higher for both treatments in the KO cells when compared to the WT cells, suggesting that the initial loss of p53 results in its own increase in R-loop formation. Observations of markers of genomic instability were also significantly increased for both treatments in the p53 KO cells compared to WT cells, with a significant increase in micronuclei and abnormal mitotic events (Fig. EV3A–D). Immunofluorescence of H3K27me3 and H3K9me3 in p53 KO cells showed a significant increase in Dac-treated cells, compared to control, while Aza-treated cells demonstrated an enrichment for H3K27me3 and a depletion of H3K9me3 signal. However, when compared to WT data, treatments and the untreated control showed a significant loss of H3K27me3 and H3K9me3 staining, linking the absence of p53 to a dramatic depletion of total H3K27me3 and H3K9me3 levels (Fig. EV3E,F).Figure 3. Removal of p53 remodels the epigenetic response to demethylating agents.(A,** B**): S9.6 staining in HEK293 p53 KO cells (A) Analysis of R-loop S9.6 staining; median nuclear fluorescence and average S9.6 foci number per nucleus, in HEK293 p53 KO cells treated with 1.5 μM Azacitidine (Aza), Decitabine (Dac) or solvent control for 48 h (grey box denotes new data) vs the same treatment in HEK293 cells (WT) (data previously shown, see Fig. 1I). n = 3 biological repeats, 13 fields analyzed. P values for the new p53 KO cell data were calculated by ordinary one-way ANOVA with Dunnett’s multiple comparisons test. P values for comparisons between the same treatment across the WT vs p53 KO cells was calculated by ordinary one-way ANOVA with Šídák’s multiple comparisons test. (Left panel p53 KO) *P = 0.0411, ***P = 0.0002 (Left panel WT vs p53 KO) *P = 0.018. (Right panel p53 KO) *P = 0.0136 (Right panel WT vs p53 KO) *P = 0.0488, ***P = 0.0005, **P = 0.0057. Error bars are 士SD. (B) Representative confocal images using S9.6 antibody to stain for R-loops (red) and DAPI as a counterstain for DNA (grey) in HEK293 p53 KO cells. Scale bars 10 μm. Boxes in the merged images represent the area of the Zoom images. (C–E) DRIPc-seq Analysis in HEK293 p53 KO cells. n = 2 biological repeats, triplicate samples were analyzed per biological replicate. (C) Feature distribution for significantly increased (gained) and decreased (lost) R-loop peak signals in cells treated with Aza. (D) Overlap of gained R-loop peaks in Aza-treated WT cells vs Aza-treated p53 KO cells. (E) DRIPc-seq R-loop peak tracks of p53-related genes in p53 KO control-treated cells vs p53 KO Aza-treated cells, + and − symbols next to peak tracks denote the positive (+) and negative (−) (Watson and Crick) DNA strands. Source data are available online for this figure.
Due to these substantially altered staining patterns in both R-loop and repressive chromatin markers, DRIPc-Seq was repeated in the HEK293 p53 KO cells, demonstrating a significant change in the R-loop profile of Aza-treated KO cells compared to WT. In p53 KO Aza-treated cells, the highest percentage of R-loop gain sites were found in intronic regions (33.0%), with those gained in promoter regions only making up 13.7% (Fig. 3C), in contrast to the 29% of newly gained promoter R-loops observed in the WT cells (Fig. 2B). The enrichment of R-loops in the promoters of p53-target genes, which had been so evident in WT cells (Fig. 2D,E), was lost in KO cells (Fig. 3D,E). Our DRIPc-Seq analysis of Dac-treated p53 KO cells encountered variability, which was sustained through both biological and technical replication (Fig. EV3G–J).
R-loop interactome
Following the observed epigenetic landscape re-modelling after both the exposures of Aza and Dac, and from the absence of p53, we sought to determine the role R-loops were playing in these responses. We reasoned that a broad-scale epigenetic re-mapping was occurring in response to these drugs, and that to achieve this, recruitment of the machinery that lays down and modifies epigenetic marks would have been required. Using DNA-RNA immunoprecipitation followed by mass spectrometry (DRIP-Mass), we assessed the proteins that were binding to R-loops in the differing conditions, to test this hypothesis (Figs. 4A and EV4A). The bulk of proteins were determined to be of nuclear origin; however, an increased prevalence of cytosolic proteins was observed in KO cells following Aza treatment (Figs. 4B and EV4A). Enrichment analysis of the WT response to Aza and Dac demonstrated an accumulation of epigenetic re-modelling proteins and complexes, including those components responsible for H3K9me3 and H3K27me3 writing (Figs. 4B–D and EV4A,B). Both conditions demonstrated accumulation of DNA methyltransferase proteins, while Dac was also found to accumulate histone acetylation writers (Fig. 4D). In KO cells, while Dac exposure also induced the enrichment of epigenetic readers and writers, as had been seen in WT cells (Figs. 4B–D and EV4B), an overlapping response between Aza WT and KO cells was less evident, with no clear enrichment for the epigenetic remodelers in KO cells (Figs. 4A–C and EV4A–C). In fact, for Aza-treated KO cells, this was the only condition where more total proteins interacting with R-loops had decreased from control conditions (Figs. 4A and EV4C).Figure 4R-loops act as a recruitment site for epigenetic modifiers and damage response elements.(A) Total number of significantly enriched (up) or significantly depleted (down) proteins from mass spectrometry of the R-loop interactome in HEK293 and HEK293 p53 KO cells treated with 1.5 µM Azacitidine (Aza) or Decitabine (Dac), compared to solvent-treated cells. All mass spectrometry data analysis is n = 3 biological repeats, three samples analyzed per biological replicate. (B) Subcellular compartment analysis of significantly enriched proteins. Bar plots show the top five identified subcellular compartments by their log10(observed/expected) enrichment score. (C) Bubble plots of the top ten hits from GO Biological term enrichment analysis of significantly enriched proteins. P values were generated through the Benjamini–Hochberg procedure. (D) STRING functional protein-protein interaction networks of significantly enriched proteins from Aza/Dac-treated HEK293 cells, clustered using the Markov cluster algorithm (MCL) (Van Dongen, 2008) via clusterMaker2 (Morris et al, 2011). Clusters and proteins of interest are highlighted using colours and annotated using their top biological function hit. Note that the same protein can appear in multiple clusters. Border colours denote clusters which are unique to Aza-treated cells (yellow) or Dac-treated cells (blue), an example of a shared cluster is shown in the centre (DNA methylation). Source data are available online for this figure.
p53 rescue
It has recently been shown that loss of p53 reduces the level of S-adenosylmethionine (SAM), an important methyl-donor, resulting in a global loss of H3K9me3 (Panatta et al, 2022). This places p53 as a powerful regulator of chromatin architecture and supports the previously shown loss of H3K9me3 and H3K27me3 observed in the HEK293 p53 KO cells (Fig. EV3E, F). In order to further solidify this relationship between p53 loss and re-modelling of chromatin markers, as well as explore a role for p53 in regulating R-loop formation, we utilised the ip53 H1299 cell line. This cell line possesses a p53 Tet-on inducible system allowing for restoration of p53 expression. In standard conditions, p53 is transcriptionally silenced, but can be rescued through the addition of doxycycline (Dox) (Fig. 5A,B) (Kirschner et al, 2015). As had been evident in the HEK293 cell line, for the ip53 H1299 cells in the absence of p53, S9.6 staining intensity was observed to significantly increase following Aza and Dac exposure (Figs. 5C,E and EV5A,B), with only modest changes in staining pattern induced by the Dox exposure (Fig. EV5C,D). Nucleolin and S9.6 co-staining in Dox^-^ cells identified a correlative pattern for Aza and Ctrl cells, with lower correlation following Dac exposure (Appendix Fig. S1A,C–E). Upon restoration of p53, the Aza and Dac induced accumulation of S9.6 staining intensity was not observed, with no significant differences from either Dox^+^ or Dox^-^ control cells (Fig. 5D,E). An increased correlation between Nucleolin and S9.6 co-staining was observed for Dac following restoration of p53 (Appendix Fig. S1B–E), suggestive of at least partial reversibility of this global accumulation of R-loop burden. Restricted H3K9me3 staining was detected in Dox^-^ Ctrl cells, which were further diminished following Aza and Dac exposure (Fig. 5F). As had been seen for the R-loop staining, return of p53 to control cells was found to rapidly restore both global H3K9me3 levels and a wild-type like response to Aza and Dac exposure (Figs. 5G,H and EV3F). This effect was further observed in H3K27me3 staining, whereby rescue of p53 levels restored staining patterns and an expected response to both drugs (Figs. 5I–K and EV3E). Finally, our earlier observations regarding DNA damage in the presence (Fig. EV1D–G) and absence of p53 (Fig. EV3A–D) were tested in the H1299 cells and confirmed with γH2AX staining. Accumulation of γH2AX signal was detected in Dac-treated cells, while significance was not achieved for Aza-treated cells. Restoration of p53 abrogated the signal in all conditions, while markers of DNA damage, which had mimicked the response observed in HEK293 cells, were also rescued (Appendix Fig. S2).Figure 5. Restoration of p53 expression reverses epigenetic dysregulation and the response to demethylating agents.(A) Graphic of the TET-ON system for doxycycline-induced p53 expression. In the presence of Doxycycline, expression of the target gene (TP53) is reversibly activated. A Tet-responsive element (TRE) sits upstream of the TP53 gene, which contains an inactive promoter. The transcription factor rtTA recognises the TRE but can only induce gene expression in the presence of tetracyclines such as Doxycycline (Dox). When Dox is present, it results in a conformational change allowing rtTA to recruit RNAPII, resulting in transcription of the TP53 gene. (B) Graphic of the experimental design. ip53 H1299 cells were seeded and treated with solvent or 100 ng/uL of Dox, to induce p53 expression, this treatment was refreshed every 24 h. Twenty-four hours after seeding, cells were treated with 1.5 μM Azacitidine (Aza), Decitabine (Dac) or solvent, this treatment was then refreshed 24 h later for a total of a 48 h treatment. Cells were then fixed and stained for immunofluorescence analysis. (C–E) Immunofluorescence staining of R-loops. (C, D) Representative confocal images using S9.6 antibody to stain for R-loops (red) and DAPI as a counterstain for DNA (blue) in ip53 H1299 cells. Cells treated with Dox express p53, whilst cells treated with solvent do not. Scale bars 10 μm. (E) Analysis of R-loop S9.6 staining; median S9.6 nuclear fluorescence in ip53 H1299 cells treated with or without Dox and either 1.5 μM Aza, Dac or solvent control for 48 h. The p values were calculated by ordinary one-way ANOVA with Tukey’s multiple comparisons test. ***P = 0.0009, *P = 0.0318, ****P < 0.0001, *P = 0.0230. Error bars are 士SD. n = 5 biological repeats, two fields analyzed per repeat. (F–H) Immunofluorescence staining of H3K9me3. (F, G) Representative confocal images of H3K9me3 staining (green) and DAPI as a counterstain for DNA (blue) in ip53 H1299 cells. Cells treated with Dox express p53, whilst cells treated with solvent do not. Scale bars 10 μm. (H) Analysis of H3K9me3 staining; median corrected total nuclear fluorescence in ip53 H1299 cells treated with or without Dox and either 1.5 μM Aza, Dac or solvent control for 48 h. The p values were calculated by ordinary one-way ANOVA with Tukey’s multiple comparisons test. ***P = 0.0003, ****P = < 0.0001, ***P = 0.0005. Error bars are 士SD. n = 5 biological repeats, two fields analyzed per repeat. (I–K) Immunofluorescence staining of H3K27me3. (I, J) Representative confocal images of H3K27me3 staining (green) and DAPI as a counterstain for DNA (blue) in ip53 H1299 cells. Cells treated with Dox express p53, whilst cells treated with solvent do not. Scale bars 10 μm. (K) Analysis of H3K27me3 staining; median corrected total nuclear fluorescence in ip53 H1299 cells treated with or without Dox and either 1.5 μM Aza, Dac or solvent control for 48 h. The p values were calculated using ordinary one-way ANOVA with Tukey’s multiple comparisons test. ****P = < 0.0001, ****P = < 0.0001, *P = 0.0314, *P = 0.0384. Error bars are 士SD. n = 5 biological repeats, two fields analyzed per repeat. Source data are available online for this figure.
Discussion
Following observations that both DNA and RNA methylation content can influence the propensity to form R-loops, we chose to evaluate the chemical removal of these marks, and whether these differentially impacted upon R-loop patterns across the genome. Our finding that p53 status is a critical factor in the response of cells to both Aza and Dac (Figs. 5C–E and 6) was unexpected, as p53 activation for either of these drugs had not been directly reported, while a link between R-loop modulation and p53 had not been demonstrated. Previous work has identified that activation of CBP/p300 occurs following Aza treatment (Diesch et al, 2021), and these proteins are known to activate p53 (Grossman, 2001). This mechanism would not be predicted to activate p53 following Dac exposure, as it is known to occur through the sensing of RNA modifications and changes to protein synthesis, which Dac is not known to affect. Further work would be required to confirm if this was the mechanism for Aza-induced activation of p53. The response to Aza was, however, dependent on the activation of p53 (Fig. 2F), including the highly enriched appearance of R-loops in the promoters of p53-target genes (Fig. 2D,E). Such R-loops were not accompanied by detectable changes to expression of these loci, while for the 111 differentially expressed genes that were detected following Aza treatment, only three were found to be enriched for an R-loop (Fig. EV2F,H).Figure 6. Summary of observed Aza and Dac phenotypes in WT and p53 KO cells.Top: In WT cells (left side), Azacitidine (Aza) treatment resulted in an increase in active p53. There was also an observed enrichment of R-loops across gene promoters alongside changes in histone methylation, a significant decrease in H3K27me3 and the recruitment of histone modifiers and helicases. There was no observed significant increase in micronuclei or aberrant mitosis. It has been shown that p53 is important for maintaining levels of S-adenosylmethionine (SAM), which was shown to be linked to a global reduction in H3K9me3 (Panatta et al, 2022). Repressive histone markers H3K9me3 and H3K27me3 are also known to be important in preventing genomic instability (Caron et al, 2021). In p53 KO cells (right side), Aza treatment resulted in an enrichment of R-loops at the 3’UTR and across introns, alongside a significant increase in H3K27me3. We hypothesise that these R-loops are involved in transcription-replication conflicts (TRCs) which is driving the observed significant increase in markers of genomic instability. There was also a reduction in the recruitment of helicases to R-loops. p53 KO cells were also observed as having a significant reduction in repressive histone markers. Bottom: In WT cells (left side), Decitabine (Dac) did not significantly increase levels of active p53. Dac treatment resulted in an enrichment of R-loops across gene terminators and introns, an increase in repressive histone markers and an increase in markers of genomic instability. We hypothesise that the enrichment of these R-loops results in increased TRCs, which drives genomic instability. There was observed enriched recruitment of histone modifiers, helicases such as Aquarius (AQR), which is associated with DNA repair (Sakasai et al, 2017), and members of the polycomb repressive complex 1 (PRC1), which are important for maintaining H3K27me3 (Dobrinić et al, 2021). In p53 KO cells (right side), Dac-treatment resulted in enrichment of R-loops at gene terminators and introns, a significant decrease in H3K27me3 and H3K9me3, and a significant increase in markers of genomic instability. There was also a decrease in recruitment of histone modifiers and helicases.
The cellular response in the case of Aza treatment was sufficient to protect against the damaging impact of the drugs on genomic stability, which was lost upon removal of p53 (Figs. EV1D–G and EV3A–D). Interestingly, in the absence of p53, a larger proportion of associating proteins following Aza treatment were found to be cytosolic, rather than in the other conditions, where nuclear proteins were highly enriched (Figs. 4B and EV4B,C). Novel work has recently identified the presence of cytosolic R-loops in specific conditions (Crossley et al, 2023), with these aberrant structures associated with cell death. While our sequencing approach did not allow for the interrogation of these structures, we observed that RNase T1 treatment, which resolves all single-stranded RNA structures, including R-loops, was sufficient to ablate cytosolic S9.6 signals (Fig. EV5A,B). Further work would be required to determine if either drug was inducing an accumulation of these cytosolic R-loops.
In contrast, the response to Dac was not determined through activation of p53 (Fig. 2F,G), but by the initial state of heterochromatin and then subsequently from an ability to remodel repressive histone modifications following exposure (Figs. EV3E,F, 4C,D, and 5F–K). The fact that we see only a handful of significant transcriptional changes following Dac treatment (Fig. 1E) is testament to this re-modelling, and the ability of the epigenome to adapt and respond to stimuli, further demonstrating the multi-layered epigenetic control that is employed on gene expression. Our analysis of the proteome that was associated with these regions of the genome demonstrated a clear link between R-loops and the re-modelling of global histone marks. Following drug exposure, recruitment of these repressive histone modifications was required to mitigate against the accumulative damaging impact. In the case of Dac, such recruitment was insufficient to protect the cell from DNA damage, as exemplified by heightened micronuclei formation and errors in mitosis (Fig. EV1D–G). When the initial heterochromatin state was destabilised by the removal of p53 (Fig. EV3E,F), this damage was amplified for Dac treatment, while also seen in Aza-treated samples (Fig. EV3A–D).
While we have attempted to understand the global epigenetic consequences of exposure to these demethylating agents, the direct locus-specific function of many of the changes to the global R-loop distribution remains opaque. We found a clear gain of promoter R-loops following Aza exposure, while terminator R-loops were especially enriched in Dac-treated samples (Fig. 2B,C). Promoter R-loops are generally thought to increase transcription (Niehrs and Luke, 2020) while terminator R-loops are known to often be transcriptionally repressive and frequently damage-inducing, a feature that was detected in Dac-treated cells (Fig. EV1D–G) (Skourti-Stathaki et al, 2014; Hatchi et al, 2015). Globally, we did not detect a clear relationship between the gain of R-loops and the location of newly demethylated sites, for either drug (Fig. EV2E). This, of course, begs the question on the functionality of these R-loops, especially given the previous observation that transcriptionally active unmethylated CpG islands are enriched for R-loops (Ginno et al, 2012).
Beyond DNA methylation, while our analysis cannot preclude that genome-wide changes to splicing patterns will have occurred, we did determine that the changes in transcription that were detected, were not accompanied by a newly formed R-loop at that locus (Fig. EV2F). A time-course analysis ensured that any transcriptional activation of these p53 target genes was not missed, allowing us to conclude that the promoter R-loops in these genes were not transcriptionally activating. The fact that our DRIPc-Seq experiments were performed 48 h after the initiation of Aza exposure indicated that, contrary to the normal highly transient nature of promoter R-loops (Castillo-Guzman and Chédin, 2021), these structures were persistent. It has previously been noted that the individual elements of the histone response to Aza or Dac, including H3K9me3 and H3K27me3, shows poor correlation with the changes observed in DNA methylation or transcription (Komashko and Farnham, 2010; Lee et al, 2009; Poplineau et al, 2015). This is surprising as the primary mechanism of Aza and Dac to induce DNA demethylation occurs through DNMT1-DNA adducts, and much of the consequent cellular response is known to initiate from these adducts (Carnie et al, 2024; Borgermann et al, 2019; Zhang et al, 2025). We therefore speculate that these newly formed R-loops are not functionally modulating transcription; instead, they represent a novel component of the epigenetic response and re-modelling that the cell undergoes following the aberrant removal of nucleic acid methylation. The diverging locations of newly formed R-loops within these drugs provide a basis for the contrasting epigenetic response that has been observed following their use. The R-loop provides a recruitment point for the machinery that lays down repressive histone modifications and DNA damage response elements (Fig. 4B–D), which together serve to maintain the appropriate transcriptional dosage and epigenetic architecture of critical loci in the genome.
A limitation of our study is that we were not able to completely characterise the differential sensitivity to these drugs, including interrogation of changes to RNA m5C levels, which would be expected to be inhibited by Aza, but not Dac. Sequencing of this mark has been historically technically challenging, and the most widely used sequencing approach in mammalian systems involves the treatment of samples with Aza, precluding its applicability to our experiments (Lee and Pelletier, 2016; Khoddami and Cairns, 2013). Recent works have highlighted that RNA modifications can preferentially dispose the modified RNA to form an R-loop (Yang et al, 2023; Abakir et al, 2020), and therefore, we might hypothesise that the impaired propensity to form specific R-loops following global removal of RNA m5C may instigate the downstream epigenetic responses that we observed following Aza treatment. Future studies will be required to explore this further, with the H1299 cell line providing a potentially ideal model for this interrogation. In contrast, adenoviral integration in the different HEK293 lines can potentially modify p53 signalling, and the rescue of wild-type or expression of mutant p53 in this line is technically challenging. Similar independent observations with Dac modulating the epigenome and R-loop staining in a mouse cancer cell line, however, provide further confidence that our observations are robust (Panatta et al, 2022).
Separately, we highlighted the variability that was observed in DRIPc-Seq and the requirement for appropriate controls (Fig. EV2A). This variability was most profound in Dac-treated conditions, and was further compounded in p53 KO settings, likely indicative of a biological effect (Figs. 1E and EV3G–J). However, previous work has highlighted that the S9.6 antibody can be subject to technical challenges (König et al, 2017; Smolka et al, 2021), and our findings provide further evidence that controlled use of this antibody is imperative. In the case of S9.6 immunostaining, we have previously speculated on the possibility of cytosolic R-loop structures representing the non-nuclear signals that can be observed, as they were resolved by RNAse T1. However, to ensure reproducibility, quantification was restricted to nuclear signals. Using appropriate controls, within both sequencing and immunostaining contexts, remains an absolute requirement of S9.6-based approaches, and caution must be taken when reading across experimental replicates.
To the best of our knowledge, p53 status is not currently considered clinically when determining treatment with Aza or Dac in AML or MDS (Maslah et al, 2024; Schimmer et al, 2022). Conflicting data exists as to the impact of loss of p53 and the response to Aza and Dac treatment, however there is some consensus that loss of p53 may be an indicator of resistance to Aza treatment, while Dac induced apoptosis is thought to be enhanced by the loss of p53 (Jung et al, 2016; Kadia et al, 2016; Chang et al, 2017; Middeke et al, 2021; Takahashi et al, 2016; DiNardo et al, 2020; Döhner et al, 2018; Bories et al, 2020; Kim et al, 2021; Nieto et al, 2004). Interestingly, substantial effort is currently being employed in utilising Aza in particular, in combination therapies, including cases of MDS which are not responsive to demethylating agents in isolation. A clinical trial has demonstrated the substantially improved efficacy of Aza when used in combination with Eprenetapopt (APR-246), a first-in-class small molecule that can restore wt p53 function in null settings (Sallman et al, 2021). Our data presented here supports these observations and provides mechanistic insights as to the underlying causes (Fig. 6), linking the Aza-specific activation of p53, and subsequent re-modelling of the epigenome to the downstream apoptotic consequences. Our models were not appropriate to dissect p53 function beyond the removal and restoration of wild-type protein. Through the expression of p53 mutants that retain expression but are limited in their ability to perform specific molecular functions, we may gain greater insight into this novel function of p53. In AML, six “hot spot mutants” in the DNA-binding domain of p53 represent the most frequently observed mutations (Papaemmanuil et al, 2016; Shin, 2023), and these likely represent an ideal place for further study. However, while our data demonstrated a clear relationship between the presence of p53 and the R-loop burden of cells under control and drug-exposed conditions, our analysis of the proteome surrounding R-loops in these conditions did not, in fact, detect the presence of p53, indicating that this role may be indirect.
Perhaps our most surprising finding was the reversibility of the effects seen following p53 removal, both in the architecture of the epigenome and the response to the drugs, which further aligns with the observed efficacy of Eprenetapopt. The staining distribution of H3K9me3 and H3K27me3 was substantially ablated in the two independent p53 null cell lines, while R-loop staining was enriched. Restoration of p53 was sufficient to rapidly reverse the destabilisation of the epigenetic landscape and restore the burden of R-loops following Aza and Dac treatment to a wild-type-like response. Given the significant breadth of epigenetic therapies currently being trialled for various disease states, the dynamism in the epigenome displayed here, that enables reversibility in ways that are not possible with genetic changes, should further encourage this area of research.
Methods
Reagents and tools tableReagent/resourceReference or sourceIdentifier or catalogue number Experimental models HEK293 cellsAnne Willis Lab giftp53 KO HEK293 cellsAnne Willis Lab giftip53 H1299 cellsMasashi Narita Lab gift. Kirschner et al, 2015.H1299 parent cellsMasashi Narita Lab gift. Kirschner et al, 2015 Recombinant DNA N/A Antibodies H3K9me3Abcamab8898β ActinMerckMA1-140Total p53Cell Signaling#2527Phospho-p53 Ser15Cell Signaling#9286S9.6MerckMABE1095H3K27me3Abcamab6002ɑ TubulinAbcamab7291γ TubulinAbcamab11317NucleolinNovus BiologicalsNB600-241γH2AXAbcamab11174Alexa-Fluor Rabbit 488Thermo Fisher#A32731Alexa-Fluor Rabbit 568Thermo Fisher#A78955Alexa-Flour Mouse 488Thermo Fisher#A-21121Alexa-Fluor Mouse 568Thermo Fisher#A-11004Alexa-Fluor Mouse 488Thermo Fisher#A32723 Oligonucleotides and other sequence-based reagents
RNASEH1 This studyF: CCTCCAGTTAGCAGAGACACGTR: CCAGTAAACGCCGATTCCTGCT TP53 This studyF: GAGGTTGGCTCTGACTGTACCR: TCCGTCCCAGTAGATTACCAC GAPDH This studyF: GTCTCCTCTGACTTCAACAGCGR: ACCACCCTGTTGCTGTAGCCAA DDIT4 This studyF: GGATGGGGTGTCGTTGCCCGR: GGCAGCTCTTGCCCTGCTCC GADD45A This studyF: CGTTTTGCTGCGAGAACGACR: GAACCCATTGATCCATGTAG MAGED2 This studyF: TAGAGAAGGCAGACGCATCCR: AAGCGAGTTAGACCTGCACC p21 This studyF: CACTCCAAACGCCGGCTGATCTTCR: TGTAGAGCGGGCCTTTGAGGCCCT Chemicals, enzymes and other reagents DAPIRoche10236276001RIPA bufferThermo Fisher89900DoxycyclineMerckD98915’-AzacitidineMerckA23855-Aza-2’deoxycitidineCellGuidance SystemsSM46-5TRIzolThermo Fisher15596026Fast SYBR Green Master MixThermo Fisher4385614cOmplete™ Protease Inhibitor CocktailMerck116974980014X Laemmli BufferBio-Rad1610747β-mercaptoethanolBio-Rad1610710NuPage 4–12% Bis-Tris gelsThermo FisherNP0323BOXRapigest SF SurfactantWaters186002123RNase HThermo FisherEN0201RNase AThermo FisherEN0531Rnase T1Thermo FisherEN0541Protein A/G Magnetic BeadsThermo Fisher88803Poly-D-LysineGibcoA3890401SSC bufferMerckS6639 Kits ApoTox-Glo Triplex kitPromegaG6320RNeasy Plus Mini kitQiagen24136QIAshredderQiagen79656HS RNA Qubit kitThermo FisherQ32852ReverAid RT kitThermo FisherK1691Pierce BCA Protein Assay kitThermo Fisher23225DNeasy Blood and Tissue kitQiagen69504QIAseq Fast Select - rRNA HMR KitQiagen334386Zymo-Seq WGBS Library kitZymoD5465Select-a-size DNA/RNA Clean and Concentrator Magbead KitZymoD4085HS DNA kitAgilent5067-4626 Software ImageJFijiProteome Discoverer v2.5Thermo FisherSTRING pathway network www.string-db.org clusterMaker2BMC Bioinformatics. 2011 Nov 9;12:436.JACoP Plugin (ImageJ)FijiFastQC (version 0.11.9) https://www.bioinformatics.babraham.ac.uk/projects/fastqc/ TrimGalore (version 0.6.6) https://www.bioinformatics.babraham.ac.uk/projects/trim_galore/ Bismark (version 0.22.3) https://www.bioinformatics.babraham.ac.uk/projects/bismark/ methylKit (version 1.20.0)Genome Biol 13, R87 (2012).tileMethylCountsR packagecalculateDiffMethR packagegetMethylDiffR packageGRanges (version 1.46.1)R packagegenomation (version 1.26.0)R packageChIPseeker (version 1.30.3)R packageannotateWithFeaturesR packageannotatePeakR packageHisat2 (version 2.2)R packageedgeR (version 3.36.0)R packagelimma (version 3.50.3)R packageBowtie2 (version 2.4.2)R packageMACS3 (version 3.0.0a6)R packageDESeq2 (version 1.34.0)R packageReactomePA (version 1.9.4)R packagePicard tools https://broadinstitute.github.io/picard/
Other DRIPc-protocolNat Protoc. 2019 May 3;14(6):1734–1755.
Cell culture
The Human embryonic kidney cell line (HEK293) and a HEK293 p53 KO cell line were used for the majority of experiments in this study. Both cell lines were cultured in Dulbecco’s Modified Eagle Medium (DMEM, Gibco) supplemented with 1% glutamax (Gibco), 10% foetal bovine serum (FBS, Gibco) and 1% Penicillin and Streptomycin (Gibco). Cells were grown under standard conditions at 37 °C and 5% CO_2_. Cells were passaged every 2–3 days at a 1:10 dilution.
A p53-inducible H1299 cell line, from here on called ip53 H1299, and the parent cell line, from now on called H1299 Parent, was kindly gifted by Masashi Narita’s lab (Cancer Research UK Cambridge Institute). Information on the establishment of the ip53 H1299 cell line can be found in (Kirschner et al, 2015). ip53 H1299 and H1299 parent cells were cultured in growth medium (DMEM supplemented with 10% FBS) under standard conditions at 37 °C and 5% CO_2_. To induce the expression of wt p53 in the ip53 H1299 cells, cells were seeded into growth media containing 100 ng/mL Doxycycline (Dox) in DMSO, with cells not receiving Dox being treated with the same volume of DMSO, as a solvent control. Dox treatment was refreshed every 24 h for a total of 72 h. Dox treatment was also performed in the H1299 Parent cells in the same manner.
Drug treatment
The drugs 5’-Azacitidine (Azacitidine) (Sigma, A2385) and 5-aza-2’-deoxycitidine (Decitabine) (Cell Guidance Systems, SM46-5) were dissolved in acetic acid and water (1:1) and stored as 5 mg/mL stocks at −80 °C for a maximum of one month. HEK293, HEK293 p53 KO cells or ip53 H1299 cells were plated at 70% confluence and left for 24 h before being treated with 1.5 µM azacitidine, decitabine or solvent for 48 h, with the treatment refreshed after 24 h. The exception to this being for the time course experiments where these cells were treated with either drug for the stated amount of time and the toxicity study where these cells were treated with a range of doses of either drug.
Toxicity assay
To assess the toxicity of Aza and Dac in HEK293 cells across a dose range, the ApoTox-Glo™ Triplex Assay (Promega) was used as per the kit instructions. Measurements were taken after 48 h of drug treatment. Fold changes were calculated compared to the untreated control samples for each drug. This kit measures cell viability, cell cytotoxicity (dead cells) and cells undergoing early apoptosis.
RNA extraction
RNA was extracted from cells using the RNeasy Mini Plus Kit (QIAGEN, #24136) following kit instructions with the additional homogenisation step using the QIAshredder (QIAGEN, #79656), for the majority of experiments. RNA was quantified on a NanoDrop spectrometer or Qubit (Thermo Fisher) with the HS RNA Kit. RNA was extracted for the qPCR time-course experiment using Trizol. About 1 mL of Trizol was added directly to cells and incubated at room temperature for 5 min. Cell lysate was transferred into an 1.5 ml tube and 0.2 mL of 100% chloroform added, vortexed well and incubated at room temperature for 2–3 min. Samples were then centrifuged at no more than 12,000 × g for 15 min at 2–8 °C. The aqueous phase was retained and washed with 0.5 mL of 100% isopropanol at room temperature for 10 min. Samples were then centrifuged at 12,000 × g for 10 min at 2–4 °C to pellet the RNA. The RNA pellet was then washed twice in 1 mL of 75% ethanol before being allowed to dry for 5–10 min. RNA was then resuspended in nuclease-free water and stored in −80 °C.
qPCR
RNA was converted to cDNA using the ReverAid RT Reverse Transcription Kit (Thermo Fisher, K1691) as per kit instructions with an input of 2 µg of RNA. qPCR reactions were performed with 2 µL of cDNA, oligonucleotides targeting specific genes, prepared as 10 µM working solutions, and Fast SYBR^TM^ Green Master Mix (Thermo Fisher). Reactions were prepared as per the protocol from Thermo Fisher for the Fast SYBR Green Master Mix. qPCR was performed using the Thermo Fisher QuantStudio 6 Flex under the following cycling conditions: 95 °C 20 s, 95 °C 3 s, 60 °C 30 s; for 40 cycles, followed by a melt curve to validate the primers used. Negative controls included a no-template (no cDNA) reaction and a minus RT reaction. GAPDH was used as a housekeeping gene for all qPCR experiments. Quantification of qPCR was performed by correcting the CT mean against the housekeeping gene to produce a ΔCT. From this, ΔΔCT and RQ values were calculated compared to the control treatment group. Sequences for oligonucleotides can be found in the Reagents and Tools table.
Western blot
Protein was extracted from cells by lysing cells directly with RIPA buffer (Thermo Fisher #89900) supplemented with 1X protease inhibitor cocktail (cOmplete Protease Inhibitor cocktail, Sigma) and 1X protein phosphatase inhibitor-2 (Sigma), for 10 min on ice. To quantify protein concentration, a BCA assay was performed using the Pierce^TM^ Protein Assay Kit (Thermo Fisher) as per kit instructions and using known standards. To generate a positive control for p53 activation, additional UV-treated HEK293 samples were prepared by exposing cells with a fluorescent tubes of 312 nm shortwave ultraviolet B. Cells were irradiated in subconfluence conditions for 30 min. Protein was then extracted in the same way. About 20 µg of protein lysate was prepared with 4X Laemmli sample buffer (Bio-Rad), supplemented with 10% β-mercaptoethanol. Denatured protein was run using NuPage 4–12% Bis-Tris Gels (Thermo Fisher). Samples were transferred onto a PVDF membrane. Membranes were then blocked with either BSA or 5% milk in tris buffered saline (TBS) with 1% Tween (TBS-T) for 2 h at room temperature. Specific proteins were detected by incubating membranes in blocking buffer with primary antibodies (listed in Reagents and tools table) at 4 °C for 12–20 h. Alexa-Fluor (Thermo Fisher) secondary antibodies at 1:5000 dilution were used for 1 h at room temperature. Antibody detection was carried out on a ChemiDocMP Imaging System (Bio-Rad) with signal analysis performed in ImageJ Fiji using the Analyze-Gels function. Band intensities, for each membrane, were then normalised to loading control (Ctrl) and fold change over Ctrl samples was calculated. For the Ctrl-treated group, fold change was calculated against the average of all Ctrl samples in a membrane by membrane basis.
DRIP-mass spectrometry
To analyse the R-loop interactome, DNA:RNA immunoprecipitation (DRIP) was coupled with label-free relative quantification mass spectrometry. Three biological repeats were analyzed for each condition. Samples were obtained for mass spectrometry from HEK293 and HEK293 p53 KO cells following treatment with 1.5 µM Azacitidine or Decitabine or solvent for 48 h. Samples were then prepared as stated in (Cristini et al, 2018) with the following changes. Following sonication, one third of each sample was treated with 5.5U of RNase H per µg of DNA overnight at 37 °C. About 25 µg of the sonicated DNA or RNase H-treated DNA was diluted 1:4 in resuspension buffer (RSB) + Tween 20 (RSB + T) buffer. Samples were then subjected to immunoprecipitation with 10 µg of S9.6 antibody (or no antibody, as a negative control) and 0.1 ng RNase A per µg genomic DNA, and rotated at 4 °C for 4 h. About 50 µL protein A/G Dynabeads (Thermo Fisher) was added to each solution and incubated for a further 2 h at 4 °C with rotation. Beads were washed four times with RSB + T and then twice with RSB buffer. Samples were then eluted in RapiGest (Waters) for 30 min at room temperature. Eluted samples were then processed by filter-aided sample preparation, and LC-MS/MS, which was kindly performed by the University of Cambridge MRC Toxicology Unit Proteomics facility. In brief, samples were precipitated in a 0.07% β-mercaptoethanol buffer followed by an in-solution digestion with AMBIC and RapiGest, protein reduction in Dithiothreitol (DTT), alkylation with Iodoacetamide, and IAA quenched with DTT. Samples were then ready for injection. Samples were injected in a randomised manner. Data processing was performed by the University of Cambridge MRC Toxicology Unit Proteomics facility. Raw data were imported and processed in Proteome Discoverer v2.5 (Thermo Fisher Scientific). The raw files were submitted to a database search using Proteome Discoverer with SequestHF and Inferis re-scoring algorithm against the Homo sapiens database containing human protein sequences from UniProt/Swiss-Prot. Common contaminant proteins (several types of human keratins, BSA and porcine trypsin) were added to the database. The spectra identification was performed with the following parameters: MS accuracy, 10 ppm; MS/MS accuracy of 0.01 Da for spectra acquired in Orbitrap analyzer; up to two missed cleavage sites allowed; carbamidomethylation of cysteine as a fixed modification; and oxidation of methionine as variable modifications. The percolator node was used for false discovery rate estimation, and only rank one peptide identifications of high confidence (FDR <1%) were accepted. PD normalisation was performed using all peptide abundances for the normalisation of the samples. Significantly enriched or depleted proteins in the treated samples, identified compared to control-treated samples and not present in negative controls, were analyzed for protein localisation, protein domain enrichment and gene ontology enrichment via the STRING pathway network software (v.12) (Szklarczyk et al, 2023). Enriched and depleted protein datasets were also analyzed in the clusterMaker2 software (v.2.3.4) for Markov cluster algorithm (MCL) (Van Dongen, 2008) clustering to produce functional cluster maps of each protein network. Raw and processed data, experimental details and analysis –of DRIP-mass of AZA and DAC treated wt and p53^−/−^ HEK293 cells, can be found within the source data file for Fig. 4.
Immunofluorescence
HEK293 and HEK293 p53 KO cells were plated at 70% confluence on sanitised glass coverslips, which had been coated with 0.05 mg/mL poly-D-lysine. ip53 H1299 cells were plated at 70% confluency on uncoated, sanitised glass coverslips. All cells were fixed in 100% ice-cold methanol for 10 min on ice. RNase T1-treated controls were made in the H1299 cells; following fixation, coverslips were submerged in 0.1% BSA in TBS-T supplemented with 3 mM magnesium chloride and a 1:200 dilution of RNase T1 (Thermo Fisher) for 1 h at room temperature. Blocking and staining was then performed the same way for all coverslips. Coverslips were blocked in 3% BSA in PBS, except for S9.6 staining, in which coverslips were blocked in 3% BSA in 4xSSC. See the antibody table for dilutions used for each antibody. DAPI was used as a counterstain to identify nuclei. All image processing and analysis was performed using Fiji software (version 1.54t). Custom Macro scripts were used for image quantification and analysis (see the Data availability section for the link to the deposited ImageJ Macro scripts). Micronuclei were identified by placing a grid over the image and counting the number of micronuclei in each grid cell. The total number of micronuclei in an image was then divided by the total number of nuclei in the image. H3K9me3, H3K27me3 and yH2AX images were analyzed using Macro scripts in which only the nuclear signal is measured. Briefly, this was done by using the DAPI stain image channel to identify nuclei and create a mask, which was then applied to the other image channels to measure signal intensity (integrated density) and signal area only within the DAPI area. R-loop signal intensity (S9.6) and signal area were calculated in the same way; however, the additional measurements of foci area and foci number were also recorded. H3K9me3 and H3K27me3 signal intensity was corrected to DAPI signal intensity by dividing the target signal by the DAPI signal for each nucleus, to produce corrected total nuclear fluorescence (CTNF). This was done to correct for any potential cross-over of fluorescence excitation between the target channel secondary antibody (488) and the DAPI channel (460). Colocalisation analysis was performed in Fiji ImageJ on images from H1299 cells stained with S9.6 and Nucleolin. The JACoP plugin (version 2.1.4) was used to calculate Pearson’s and Mander’s (M1 = S9.6, M2 = Nucleolin) coefficient scores and generate cytofluorogram data. The cytofluorogram data of five replicates was then combined and plotted in R. All scale bars on images were added through the Fiji software. Scales of scale bars are described for each figure in the figure legend.
WGBS sample prep
DNA was extracted from HEK293 cells, which had been treated with 1.5 µM of azacitidine or decitabine or solvent for 48 h, using the QIAGEN DNeasy Blood and Tissue Extraction Kit, as per kit instructions. Whole-genome bisulfite sequencing (WGBS) conversion and library construction was then performed using Zymo-Seq WGBS Library Kit (Zymo Research) as per kit instructions. The Select-a-size DNA/RNA Clean and Concentrator Magbead Kit (Zymo Research) was used when needed to remove primer dimer contamination from samples. Libraries were visualised before sequencing using the High Sensitivity DNA Kit and Chip on a BioAnalyser (Agilent). All libraries were combined into a single pool and sequenced at the University of Cambridge Biochemistry sequencing facility using the Next-Seq 500 (Illumina), 300-cycles 150 bp paired-end reads.
WGBS data processing and data analysis
FastQC (version 0.11.9) was used to assess the overall quality of the sequenced samples, and TrimGalore (version 0.6.6) was used to trim low-quality bases (quality score lower than 20), adaptor sequence, and end-repair bases from the 3ʹ end of reads. Bismark (version 0.22.3) was used for alignment to the human hg19 reference genome, deduplication of reads and methylation calling (Krueger and Andrews, 2011). The methylKit R package (version 1.20.0) was used to test WGBS data for differential methylation for the Aza and Dac exposures compared pairwise to the control group (Akalin et al, 2012). Methylation calls were reported for all nucleotides with a read depth of at least 5. In addition to individual CpG sites, we used the function tileMethylCounts to determine methylated/unmethylated base counts over 1000 bp tiling windows across the genome. The calculateDiffMeth function was used to identify differentially methylated CpGs or 1000 bp tiles. To filter for differentially methylated regions (DMRs), we used the getMethylDiff function with a q-value cut-off of 0.01 and a difference cutoff (absolute value of methylation percentage change) of 25. Referencing to other sequencing data (e.g., DRIPc) was achieved through the GRanges (version 1.46.1) R package by transforming the respective regions into GRanges objects on which, e.g. overlap operations could be conducted (subsetByOverlaps) (Lawrence et al, 2013). The genomation (version 1.26.0) and ChIPseeker (version 1.30.3) R packages were employed to annotate DMRs to the TxDb hg19 genome (TxDb.Hsapiens.UCSC.hg19.knownGene) (Akalin et al, 2015; Yu et al, 2015). The functions annotateWithFeatures and annotatePeak were used to annotate the abundance of regions within custom and TxDb genetic features, respectively.
RNA-Seq library preparation, data processing and data analysis
RNA was extracted, as described above, using the kit method. The QIAseq Fast Select - rRNA HMR Kit (Qiagen) was used to prepare RNA-sequencing libraries, with 1 µg RNA as input, as per kit instructions. To remove primer dimers, the Select-a-size DNA/RNA Clean and Concentrator Magbead Kit (Zymo Research) was used, with the desired peak 300 bp. Libraries were then pooled at 10 nM in 10 µL. RNA-sequencing (RNA-Seq) was performed using the Next-Seq 500 (Illumina) for 150-cycles paired-end reads, by the University of Cambridge Biochemistry Department sequencing facility. Quality of data in fastq files was assessed using FastQC, and low-quality bases (quality score lower than 20) and adaptor sequences were removed using TrimGalore. Reads were aligned to the human hg19 reference genome using Hisat2 (version 2.2), and raw read counts were generated using featureCounts. The data were normalised for composition bias using the calcNormFactors function from the edgeR R package (version 3.36.0), and normalised again and transformed to log2 counts per million reads (CPM) using the voom function from the limma R package (version 3.50.3) (Robinson et al, 2010; Ritchie et al, 2015). To determine differential expression, the contrast.fit and eBayes limma functions were used to compute estimated coefficients and standard errors, and apply empirical Bayes moderation, respectively. The p values were adjusted by using the Benjamini–Hochberg false discovery rate (FDR) approach, applying an FDR <0.05 and log2(fold change) >1 to determine significance.
DRIPc-seq sample prep
DNA:RNA immunoprecipitation with cDNA conversion sequencing (DRIPc-Seq) was performed on HEK293 and HEK293 p53 KO cells, treated with Azacitidine, Decitabine or solvent, as described above. DRIPc-Seq was performed according to the published method (Sanz and Chédin, 2019), with the exception that the genomic DNA was fragmented using sonication to achieve a fragment size of around 300 bp. An RNase H-treated control was used as the negative control. Libraries were then prepared using the TruSeq RNA v2 kit (Illumina) by the University of Cambridge Biochemistry sequencing facility.
DRIPc-seq data processing and data analysis
All libraries were sequenced at the University of Cambridge Biochemistry sequencing facility using the Next-Seq 500 (Illumina), 150-cycles single-end reads. As described above, FastQC was used to assess the overall quality of the samples, and TrimGalore was used to trim low-quality bases and adaptor sequences. Trimmed reads were mapped to the UCSC human hg19 genome using Bowtie2 (version 2.4.2), and the mapped reads were deduplicated using the MarkDuplicates command in Picard tools. Strand-specific R-loop peak calling was performed using MACEV3 (version 3.0.0a6) with the --broad --nomodel --broad-cutoff 0.2 --max-gap 350 --extsize 150 settings. The resulting peaks were merged into union R-loop peaks in a strand-specific manner, collapsing the peaks into non-redundant and distinct peaks. We then compiled a high-confidence peak list that contained only peaks that are present in at least two datasets and had not remained in any of the RNase H-treated samples. To identify R-loop signals that emerge upon treatment, the differentially up- or downregulated R-loop signals for all high confidence peaks were determined using the R package DESeq2 (version 1.34.0) with a cutoff Benjamini–Hochberg adjusted P value (adjusted P value) <0.05 and a log2(fold change) >0.5 to identify gained R-loops or log2(fold change) <−0.5 for lost R-loops. Importantly, changes were only determined with a control treatment that was measured alongside the respective treated sample to account for any underlying changes in R-loop architecture between experiments. The resulting peaks were annotated to the hg19 genome using CHIPseeker and genomation, and pathway-enrichment was performed using the ReactomePA R package (version 1.9.4) (Yu and He, 2016; Yu et al, 2015).
Quantification and statistical analysis
Sample sizes were chosen according to standard guidelines. Sample blinding was performed post-harvesting of cells, with researchers unblinded for the final statistical assessment. Various statistical tests were used in this study; specific tests are described in the corresponding figure legends or methods sections. All error bars shown are ±SEM unless otherwise stated in the figure legend. For all experiments, the following nomenclature was used: *p < 0.05, **p < 0.01, ***p < 0.001, ****p < 0.0001. Unless otherwise stated, statistical analysis was performed in R (v.4.3.2).
Supplementary information
Appendix Peer Review File Dataset EV1 Source data Fig. 1 Source data Fig. 2 Source data Fig. 3 Source data Fig. 4 Source data Fig. 5 Figure EV1 Source Data Figure EV3 Source Data Figure EV4 Source Data Figure EV5 Source Data Expanded View Figures