Genome-wide 5hmC profiling of the genitourinary system in health and cancer
To precisely profile genome-wide 5hmC patterns in the genitourinary system, we performed hMeDIP-seq of 10 tumors (3 kidney cancer, 3 prostate cancer and 4 urothelial cancer) and paired normal tissues (Fig. 1A). As expected, unsupervised clustering uncovered the existence of tissue-specific 5hmC patterns in healthy and cancerous GU tissues (Additional file 1: Fig. S1A, B). Interestingly, we found that normal renal tissues had the highest 5hmC level, while normal urothelial tissues showed the lowest 5hmC modification at the gene body and 1-kb flanking region (Fig. 1B, Additional file 1: Fig. S1C). The global 5hmC pattern of healthy genitourinary samples was further validated by dot blot using normal tissues and three normal GU cell lines (Fig. 1C and Additional file 1: Fig. S1D). Consistent with previous reports [2, 3, 7], the meta-gene analysis suggested a global loss of 5hmC in all of the 3 types of GU cancers compared with the matched normal tissues, which was also confirmed by the reduced number of valid 5hmC peaks (Fig. 1B–D). By combining the normal and tumor peaks, we found that the kidney cancer had decreased hydroxymethylation at the 3’-UTR, TES and introns, while the urothelial cancer had a widespread dehydroxymethylation (Additional file 1: Fig. S1E).
Tumors are disease of high heterogeneity. Although cancers of the kidney, prostate and urothelium had distinct DNA 5hmC landscape, we still wondered whether there is a universal change of biological process underlying the tumorigenesis of different GU cancers. Therefore, we listed the 5hmC peaks for the three normal and cancerous GU tissues and excluded those consecutively active for each type of tissue (Additional file 2: Fig. S2A, B). The remaining normal and tumor-specific peaks were intersected (Additional file 2: Fig. S2A–D). To this end, we identified 3787 and 1363 5hmC peaks that were shared by cancerous and healthy GU tissues, respectively (Fig. 1E). While 5hmC peaks shared by normal GU tissues influenced genes involving tissue development and homeostasis, common 5hmC peaks in GU cancers were enriched for several oncogenic GO terms including the regulation of insulin-like growth factor and Notch signaling[37] (Fig. 1F). Specifically, a consecutive hydroxymethylation at the promoter and gene body of HOXB8 and HOXB9 was seen in normal GU tissues (Additional file 2: Fig. S2E, F). Similarly, the promoter of the oncogene WNT7B was uniformly hydroxymethylated in GU cancers but not normal tissues, which could also be validated by hMeDIP-qPCR (Fig. 1G, H). Collectively, these results offered a genome-wide 5hmC landscape of the genitourinary system in both health and cancer, which revealed common 5hmC alterations contributing to the progression of GU tumors.
Aberrant gain and loss of 5hmC coexist in genitourinary cancers and are associated with different oncogenic pathways
To systemically investigate the disruption of DNA hydroxymethylation during tumorigenesis, we identified differentially hydroxymethylated regions (DhMRs) in cancerous genitourinary tissues (Additional file 3: Fig. S3A). A total of 24,718, 13,090, and 15,233 peaks were subsequently obtained as DhMRs for the kidney, prostate, and urothelial cancers, of which the majority were hypo-DhMRs (Fig. 2A, Additional file 3: Fig. S3B, and Additional file 10: Table S2). In addition, the preference of hypo-hydroxymethylation over hyper-hydroxymethylation was observed at almost all autosomes in the three GU cancers we investigated (Additional file 3: Fig. S3C). Notably, both hypo- and hyper-DhMRs we collected were highly tissue-specific (Fig. 2B, C).
The annotation of hyper-DhMRs suggested the presence of substantial differently enriched gene sets in the three GU cancers (Fig. 2G–I). In kidney cancer, the hyper-DhMRs were significantly associated with hypoxia (Fig. 2G). Consistently, hyper-hydroxymethylation of HIF1A was observed (Fig. 2D). In prostate cancer, there was a strong enrichment of androgen signaling encompassing hyper-hydroxymethylated ARHGEF16 (Fig. 2E–H). Additionally, the urothelial carcinoma harbored elevated 5hmC at the promoter and gene body of FGFR3, which was reported to play a central role in its tumorigenesis [38] (Fig. 2F–I). However, hypo-DhMRs were linked to genes involving biosynthesis and metabolism irrespective of tumor type (Fig. 2G–I).
Functionally, we found that these DhMRs affected genes were positively changed in mRNA expression during the tumorigenesis of the corresponding GU cancers by analyzing the RNA-seq data from The Cancer Genome Atlas (TCGA) (Additional file 4: Fig. S4A–C). Collectively, our results suggested that the tumorigenesis in the genitourinary system is accompanied by both gain and loss of DNA 5hmC, which is tissue-specific but also indicative of a universal tumor biology.
Genes associated with 5hmC alterations predict clinical outcome in genitourinary cancers
State of 5hmC is closely related to gene expression. By integrating the RNA-seq data of GU cancers from TCGA with our hMeDIP-seq data, we found a significant positive correlation between 5hmC alterations and changes in mRNA expression during the tumorigenesis of renal cancer, prostate cancer and urothelial carcinomas (Fig. 3A and Additional file 5: Fig. S5A, Additional file 11: Table S3). Since many gene signatures have been reported to be prognostic biomarkers for GU cancers [39], we sought to evaluate the efficacy of genes associated with hypo- and hyper-hydroxymethylation in predicting patient survival. To generate such gene signatures, we collected these hypo- and hyper-DhMRs in each type of cancer and annotated them to genes, from which top 10 candidates were selected according to the expression fold change in tumors (Fig. 3B–D, Additional file 5: Fig. S5B–D). Surprisingly, patients with high hyper-hydroxymethylated signature had significantly worse overall survival or progression-free survival (PFS) in the TCGA KIRC, BLCA and PRAD cohorts, which illustrates that genes associated with 5hmC alterations are of potential clinical significance (Fig. 3E–G, Additional file 5: Fig. S5E–G).
Increased stemness and DNA 5hmC in 3D-cultured prostate cancer cells
Cancer progression to an aggressive phenotype often co-opts gain of stem cell-like states and epigenome reprogramming [18, 40, 41]. To gain more insights into such epigenetic reprograming processes, we applied a previously described mechanical method to obtain aggressive cancer stem-like cells by culturing prostate cancer cells in fibrin matrices of ∼90 Pa stiffness [42, 43], which resulted in the formation and growth of multicellular tumor spheroids (Fig. 4A). Using real-time reverse transcription-PCR (RT-qPCR), we discovered that several canonical stem cell makers, including NANOG, CD133, and OCT4, were up-regulated in the soft-fibrin-gel-selected cancer cells (Fig. 4B, C), which confirmed their cancer stem cell identity. Additionally, a significant increase in tumor volume and mass for the 3D-cultured group compared with two-dimensional (2D) controls was observed in the mouse xenograft model, suggesting that 3D-cultured prostate cancer cells acquired an aggressive cancer stem cells (CSCs) phenotype (Fig. 4D–F).
Having generated CSCs from prostate cancer cells, we sought to better resolve cellular 5hmC dynamics introduced by 3D-culture. Interestingly, using dot blot we identified a global increase of 5hmC in tumor spheres induced by 3D-culturing (Fig. 4G), which was also evident in subsequent hMeDIP sequencing (Fig. 4H, I). A total of 1,214 and 1,671 genes were identified as being hyper-hydroxymethylated in the 3D-cultured 22RV1 and PC3 cells compared with 2D controls, respectively (Additional file 12: Table S4). By annotating them to relevant biological processes, we noticed significant enrichments in signaling pathways involving mTOR, WNT, Notch, C-myc and hypoxia, which might contribute to the maintenance of stemness and promote tumor progression (Fig. 4J-K). Intriguingly, these enriched GO terms coincide well with what describes hyper-DhMRs in prostate tumors (Fig. 2H), suggesting that specific 5hmC gains might facilitate the initiation and progression of prostate cancer and are associated with cancer cell stemness.
Stemness-associated 5hmC gain defines a highly malignant cancer cell population in prostate tumors
We next sought to further investigate the clinical implication of the stemness-associated 5hmC gain. By coupling with RNA-seq, we found that the 3D-culturing gave rise to synchronous enhancements in the gene expression and DNA 5hmC of KITLG, SOX21, CDH1 and WNT9A, which are known stem cell maintainers [44,45,46] (Additional file 6: Fig. S6A). Moreover, enrichment scores of gene sets comprising simultaneously up-regulated and hyper-hydroxymethylated genes in 3D-cultured PC3 (n = 166) and 22RV1 cells (n = 171) positively correlated with the Gleason score in the TCGA PRAD cohort (Additional file 6: Fig. S6A-B), which resembles to the adult stem cell signature, naive hESC signature and primed hESC signature described in previous studies [47], implying the critical role of the regional 5hmC gain in the tumorigenesis and progression of prostate cancer (Additional file 6: Fig. S6C).
Among genes exhibiting hyper-hydroxymethylation and elevated expression in both PRAD tissues and 3D-cultured cells, 7 genes (NVL, LUC7L3, TTC13, RHOU, CABLES2, ZFPM1 and ZSWIM4) were associated with poor PFS in the TCGA PRAD cohort (Fig. 5A). Using hMeDIP-qPCR and RT-qPCR, we validated that higher 5hmC and expression levels of these genes were readily detected in 3D-cultured 22RV1 cells (Additional file 6: Fig. S6D, E). A malignant signature comprising the 7 genes was hence developed for the quantification of the stemness-associated 5hmC gain of individual tumors. As expected, while the signature was positively associated with the previously established hyper-hydroxymethylation signature (Fig. 5B), it was also indicative of the PFS and the Gleason score in the TCGA PRAD cohort (Fig. 5C, D and Additional file 6: Fig. S6F). Given that tumors are heterogeneous entities, we next explored the presence of tumor cells with stemness-associated 5hmC gain in prostate cancer. The Scissor [33] algorithm was applied to a published single-cell RNA-seq dataset of prostate cancer comprising basal/intermediate, luminal and CellCycle subtypes of tumor cells (Fig. 5E) [34]. A total of 1248 Scissor-positive (Scissor+) tumor cells were subsequently identified, while the rest 1059 tumor cells were assigned as Scissor-negative (Scissor−) (Fig. 5F). The Scissor+ rate was higher in the CellCycle subtype than the basal/intermediate subtype (Fig. 5G), which is consistent with previous reports that CellCycle and basal/intermediate cells were associated with worse and better clinical outcomes, respectively [34]. Additionally, Scissor+ tumor cells showed both a higher malignant score and a higher adult stem cell signature than Scissor− cells (Fig. 5H). Gene set enrichment analysis also suggested that stem cell-related pathways, such as hypoxia, MTORC1 signaling, Notch signaling and somatic stem cell maintenance, were enriched in Scissor+ tumor cells (Fig. 5I). Collectively, the identification of Scissor+ tumor cells emphasized the clinical relevance of the malignant score and confirmed the presence of tumor cell subpopulations with stemness-associated 5hmC gain in prostate tumors.
Inhibition of the growth of tumor stem cell-like cells by vitamin C
We recently showed that the oxidation-resistant vitamin C derivative, APM, can restore the 5hmC pattern by acting as a cofactor of TETs in kidney cancer and bladder cancer [3]. However, it remains unclear whether vitamin C can also inhibit prostate cancer cell growth especially cancer stem cell-like cells by reprogramming 5hmC. As anticipated, APM treatment restored the 5hmC level and inhibited the proliferation of PC3 and 22RV1 cells (Fig. 6A–C, Additional file 7: Fig. S7A, B). Interestingly, in the 3D culture system, treatment of prostate cancer cells with APM impaired not only the initiation of 3D tumor spheres but also the growth of established ones in a dose-dependent manner (Fig. 6D–F), suggesting that APM may be able to inhibit the proliferation of prostate cancer stem-like cells.
To investigate whether any epigenetic events were involved in this process, we performed hMeDIP-seq, which in turn revealed that APM-treated PC3 and 22RV1 cells acquired 5hmC landscapes more similar to normal prostate tissues (Fig. 6G). For example, the 5hmC level and expression of CARD14 were remarkably increased upon APM treatment to match what was observed in the normal prostate (Fig. 6H–J). Of note, genes restoring 5hmC level upon APM treatment were enriched in apoptosis-related cellular processes (Fig. 6H and Additional file 7: Fig. S7C), implying a pivotal role of apoptosis in APM-induced inhibition. The result of Gene Ontology (GO) analysis on RNA-seq data of APM-treated cells and 3D gel cultured cells also showed activated apoptotic process in APM-treated 22RV1 and PC3 cell lines, but negative regulation in 3D-cultured cell lines (Additional file 7: Fig. S7D), suggesting the potential regulatory role of APM in apoptosis-related processes. Furthermore, we performed flow cytometry analysis, which confirmed that APM treatment increased the proportion of apoptotic prostate cancer cells in both 2D and 3D culture systems (Additional file 7: Fig. S7E, F). Collectively, these results suggested that APM could be novel epigenetic therapeutics targeting cancer stem cell-like cells and inducing apoptosis for prostate cancers.