Hsa-miR-1-2/miR-133a-1 Cluster Expression Silencing by DNA Hypermethylation Triggers Oncogenes Contributing to Prostate Cancer Aggressiveness

Background: Prostate cancer (PCa) is currently the most commonly diagnosed type of cancer. The incidence and mortality of PCa worldwide correlate with increasing age and bad dietary habits. Previously, we investigated the mRNA/miRNA role on PCa development and progression using high fat diet (HFD) chronically fed mice models. Here our main goal was to investigate the effect of HFD on the expression of prostate cancer-related miRNAs and their relevance in PCa patients. Methods: Using microarray data from mice prostate tumors generated by TRAMP-C1 cell inoculation, we focused on the role of three candidate miRNAs (miR-133a-3p/133b/1-3p) and their target genes. Based on data from public databases, we examined the expression levels of hsa-miR-133a-3p/133b/1-3p and their correlation with clinicopathological features in PCa patients. The biological roles of hsa-miR-133a-3p/133b/1-3p and their relevant target genes were investigated by bioinformatics approaches. The promoter methylation of hsa-miR-133a-3p/133b/1-3p host genes and their correlation with mature miRNA expression was further evaluated. Results: We identied 6 up- and 18 down-regulated miRNAs in TRAMP-C1 mice prostate tumors under HFD conditions using miRNA microarrays. Target genes (1,278) of down-regulated miRNAs involved in cancer-related biological processes were identied using DIANA-TarBase and STRING databases. Three down-regulated miRNAs: hsa-miR-133a-3p, 133b and 1a-3p showed nine common target genes that negatively correlated with miRNA expression in prostate tumors from patients. Hsa-miR-133a-3p/133b/1-3p expression levels were signicantly decreased in PCa compared to normal tissues and their low expression correlated with bad clinicopathological features in patients. We also examined the promoter region of hsa-miR-133a-3p/133b/1-3p encoding genes in PCa patients and then compared methylation at these loci with mature miRNA expression. This analysis revealed that hsa-miR-1-2/miR-133a-1

obesity have the potential to cause PCa initiation, promotion and progression [2]. The proposed mechanisms for PCa induced by dietary fats are divided into growth factor signaling, lipid metabolism, in ammation and hormonal modulation among others [2]. However, the underlying molecular mechanisms responsible for the effect of high fat diet (HFD) on PCa development and progression remained unknown.
Previously, we generated several preclinical mice models to study the in uence of HFD on PCa development and progression. We reported that C-terminal binding protein 1 (CTBP1) depletion in androgen-insensitive PCa xenografts from HFD-fed mice modulated the expression of mRNAs and microRNAs (miRNAs) involved in cancer related pathways which impacts on PCa proliferation and invasion [3][4][5][6]. Additionally, recent androgen-sensitive PCa allografts and HFD mice model demonstrated that high fat intake signi cantly increased tumor growth. Tumors developed in HFD fed mice showed overexpression of oncogenes and oncomiRs compared to control diet (CD) [4,7,8].
MiRNAs are endogenous small non-coding RNAs (18-22 nucleotides) that regulate gene expression. Compelling evidence have demonstrated that miRNA expression is deregulated in several human cancer types through numerous mechanisms, including ampli cation or deletion of miRNA genes, abnormal miRNAs transcription, deregulated epigenetic changes and defects in the miRNA biogenesis machinery [9]. MiRNAs can function either as oncogenes (oncomiRs) or tumor suppressors (tsmiRs) under certain conditions. The deregulated miRNAs have been shown to affect several hallmarks of cancer, including sustaining proliferative signaling, evading growth suppressors, resisting cell death, activating invasion and metastasis, and inducing angiogenesis [9]. In PCa, several miRNAs have been proposed to regulate cell proliferation, cell cycle, apoptosis, as well as invasion and adhesion processes [10]. In addition, diet and lifestyle factors are involved in the regulation of miRNA expression in different tissues and pathologies, including cancer [11][12][13]. Finally, an increasing number of studies have identi ed miRNAs as potential biomarkers for PCa diagnosis, prognosis and therapy.
Here our main goal was to investigate the effect of HFD on the expression of cancer-related miRNAs in prostate tumors and their relevance in PCa patients. Starting from microarray data from mice prostate tumors we focused on the role of three candidate miRNAs (miR-133a-3p/133b/1-3p) and their target genes. Using PCa patient samples from public databases, we examined the expression levels of hsa-miR-133a-3p/133b/1-3p and their correlation with clinicopathological features in PCa patients. The biological roles of hsa-miR-133a-3p/133b/1-3p and their relevant target genes were investigated by bioinformatics approaches. Finally, the promoter methylation of hsa-miR-133a-3p/133b/1-3p host genes and their correlation with mature miRNA expression was evaluated.

Cell culture
Murine TRAMP-C1 cell line was obtained from the American Type Cell Culture Collection (Manassas, VA) (ATCC: CRL-2730). Cells were cultured in DMEM medium (GIBCO, Thermo Scienti c, Massachusetts, USA) supplemented with 10% of fetal bovine serum, antibiotics and 0.25 IU/µl of human recombinant insulin in a 5% CO2 humidi ed incubator at 37°C.

PCa allograft and HFD murine model
We previously reported that HFD signi cantly increased TRAMP-C1 tumor growth, oncogene and oncomiRs expression [4,7]. Here, we used the same model to assess miRNA expression in these tumors.
Brie y, four-weeks-old C57BL/6J male mice (N = 16) were housed under pathogen-free conditions following the IBYME's animal care guidelines. Mice were randomized into two dietary groups and fed ad libitum during 27 weeks with CD (3,120 kcal/kg, 5% fat) or HFD (4,520 kcal/kg, 37% fat). After 15 weeks of diet, 3 × 10 6 TRAMP-C1 cells were subcutaneously injected. Animals were sacri ced in the 27th week and tumor samples were collected. The metabolic state of the animals and tumor volume were previously reported [4].

RNA isolation and microarrays analysis
For microarray analysis, total RNA from CD or HFD allografts was isolated using TriReagent (Molecular Research Center) and hybridized with GeneChip® miRNA 4.0 Array (Affymetrix) (N = 3 per group). For miRNA expression analysis, we employed the Limma and pd.mirna.4.0 packages in the R/Bioconductor environment. For differential expression analysis we used the Rank Product Method for two class unpaired data and a fold discovery rate (FDR) < 0.05 [14].

Functional enrichment analysis
To investigate the functional role of differentially expressed miRNAs between CD and HFD tumors we performed Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis using DIANA-miRPath v3 tool (http://snf-515788.vm.okeanos.grnet.gr/) and STRING database (https://string-db.org/). Brie y, we obtained a list of experimentally validated target genes derived from the differentially expressed mice miRNAs and their human orthologous from DIANA-TarBase v7 that were entered into the STRING database for KEGG analysis.

Principal component analysis (PCA)
Publicly available log2(norm_count + 1) gene expression values for TCGA Prostate Cancer (PRAD) and the GTEx project patient samples were downloaded from UCSC's Xena Browser (http://xena.ucsc.edu/) [15]. PCA was performed to determine samples distribution based on the expression of the target genes identi ed from the downregulated miRNAs. We included in the analysis normal prostates (n = 100) (GTEx), normal adjacent tissue (NAT) (n = 52) and prostate tumors (n = 497) (TCGA-PRAD) samples. For PCA plots, the R function "prcomp" from stats package (version 4.0.2) was used.

Correlation matrix
Expression levels of down-regulated miRNAs in HFD tumors and their respective target genes were obtained from TCGA-PRAD patient cohort data available in UCSC Xena. Expression data of miRNA mature strand [miRNA-Seq (IlluminaHiSeq_miRNASeq)] and genes [RNAseq (IlluminaHiSeq)] were downloaded as log2 (RPM + 1) values. PCa samples from 497 patients were included in the present analysis. We generated a correlation matrix for each miRNA and their target genes applying the Spearman correlation coe cient using the Hmisc R package. For all down-regulated miRNAs, target genes with a negative correlation coe cient rho <-0.2 plus p-Value < 0.05 were selected for further analysis.
Single-sample gene-set enrichment analysis (ssGSEA) A ssGSEA was performed to analyze the coordinated regulation of a de ned gene signature and the activation of speci c biological processes in PCa samples. Log2 (RPM + 1) gene expression values of 497 PCa samples were obtained from TCGA-PRAD cohort and loaded into GenePattern web-tool (https://www.genepattern.org). A gene set enrichment pro le was obtained for each sample using a ninegene signature obtained from the correlation matrix described above. Then, Spearman correlation analysis was performed between the expression of three down-miRNAs and the enrichment score obtained from the ssGSEA (GraphPad Prism 8.0.1). Afterwards, to analyze the contribution of the 9 genes in speci c biological processes a compute overlaps was performed against C1, C2, C3, C4, C5, C6, C7, C8 and H collections in Molecular Signatures Database (MSigDB) (http://www.gseamsigdb.org/gsea/msigdb/). The top 10 genesets with FDR q-value less than 0.05 were selected to graphed the gene/geneset overlap matrix.

Bioinformatics analysis databases
MiRNA and gene expression levels in addition to clinical-pathological data (e.g. Gleason score, pathological state N or T, etc.) were obtained from two independent PCa patient datasets. The RNA-seq and miRNA-seq data from TCGA-PRAD cohort was downloaded using UCSC Xena resource. PCa samples (n = 497) and NAT (n = 52) with expression data of miRNA mature strand [miRNA-Seq For statistical analysis normality of data was assessed using Shapiro-Wilk test and homogeneity of variances was analyzed by boxplot. We used a paired sample t-test or a Sign Median test to compare the expression between PCa and NAT samples (N = 52 per group) and a Mann Whitney test to compare the expression between PCa and normal prostate samples (normal N = 28, PCa N = 111). Also, one-way ANOVA followed by Tukey or Kruskal-Wallis test followed by Dunn were used to analyze the statistical differences between PCa samples according to Gleason grade (N = 497 TCGA-PRAD, N = 111 Taylor).
Progression free interval curves were plotted using the Kaplan-Meier method and compared by log-rank test using Xena resource from TCGA-PRAD cohort.

HFD modulate the expression of multiple miRNAs in murine prostate tumors
Previously, we reported that HFD signi cantly increased TRAMP-C1 tumor growth, oncogenes and oncomiRs expression in C57BL/6J male mice [4,7]. In this work, we investigated the impact of HFD on miRNA expression pro le from these tumors using a high-throughput platform. GeneChiP miRNA 4.0 Affymetrix was hybridized with total RNA isolated from CD or HFD TRAMP-C1 tumors. After data normalization, we selected differentially expressed miRNAs with a FDR<0.05 and identi ed 18 down-and 6 up-regulated miRNAs in HFD tumors compared to CD group ( Figure 1, Table 1). Particularly, mmu-miR-133a-3p, 133b and 1a-3p showed the most robust expression changes between the two groups. As miR-133a-3p and 1a-3p are expressed from the same miR-1/133a genomic clusters (located on chromosomes 2 and 18 in mouse and 18 and 20 in human and mmu-miR-133a-3p and 133b are part of the same miRNA gene family they were selected for further analyses.
Since each miRNA can control several genes, we selected only the deregulated miRNAs that had orthologous in human to obtain a more robust result. We employed DIANA-TarBase v7 of the miRPath tool, which is a reference database of experimentally supported miRNA targets, and followed the Figure   2A work ow. Eight of the down-modulated mouse miRNAs (8 out of 18) and 10 of its orthologous human miRNAs (10 out of 18) showed validated target genes ( Figure 2AI). Additionally, two (2 out of 6) of the upregulated mouse miRNAs and 2 of its human orthologous presented validated target genes. Next, we excluded common target genes among down-and up-regulated miRNAs from the same species and selected common target genes between species (human and mice). This analysis revealed that downregulated miRNAs showed 1,278 target genes common to human and mice, while up-regulated miRNAs did not present common target genes (Table S1). We continue working with down-modulated miRNAs and their target genes. For this group of genes, KEGG signaling pathways were identi ed with a FDR< 0.05 using STRING database ( Figure 2B, Table S2). Target genes of the down-modulated miRNAs showed enriched processes mainly associated with cancer-related and signal transduction pathways ( Figure 2B).
To determine the relevance of down-regulated miRNAs-target genes in human samples, we performed a PCA of the 1,278 target genes using the normalized expression data from normal prostate (GTEX), NAT and prostate tumors samples (TCGA-PRAD). The 2-dimensional scatterplot of the rst two principal components revealed marked differences in overall gene expression between normal prostate and PCA samples ( Figure 2C).
In order to found relevant target genes for the down-modulated miRNAs, within the 1,278 genes, we performed a correlation matrix, using the expression data of the down-modulated miRNAs and their target genes of prostate tumors from patients available in TCGA-PRAD (UCSC Xena). From this analysis, we selected the target genes that showed a negative Spearman correlation coe cient rho <-0.2 and a p-Value <0.05 ( Figure 2AII) with its regulatory miRNA. Next, to shorten the number of relevant miRNAs and target genes, we choose common target genes that negatively correlate with the expression of the previously selected miRNAs hsa-miR-133a-3p, 133b and 1-3p. We found nine common target genes (BCL7A, CENPF, CPSF6, ELAVL1, NRP1, PAICS, RALA, RBM15B and WHSC1) with a signi cant negative correlation with the mentioned miRNAs ( Figure 3A, Table S3). Also, hsa-miR-133a-3p, 133b and 1-3p showed a strong positive correlation between them with a signi cant p-value ( Figure 3A). Additionally, we performed a ssGSEA to explore the degree to which the 9 common target genes were coordinately up-or down-regulated within the prostate tumor samples in two different PCa datasets: TCGA-PRAD and Taylor (GSE21032). As shown in Figure 3B-C we found that the 9 target genes showed a positive ssGSEA-enrichment in most of PCa samples in both datasets. Moreover, a strong negative correlation between the ssGSEA-enrichment and the expression of miRs 133a-3p, miR-133b and miR-1-3p were found in the TCGA-PRAD dataset. A similar result was found for 133a-3p, miR-133b in Taylor dataset. Finally, we analyzed the contribution of the 9 target genes to speci c gene sets annotated in MSigDB collections. As shown in Figure 3D, three of the nine target genes (WHSC1, CENPF and PAICS) were found to be signi cantly enriched in cancerrelated gene sets.
Hsa-miR-133a-3p, miR-133b and miR-1-3p expression are decreased in PCa In order to further understand the pattern expression of the three selected miRNAs in prostate tumors from patients, we performed a bioinformatic analysis using the TCGA-PRAD, and Taylor (GSE21032) datasets. First, we evaluated the expression of the selected miRNAs in prostate primary tumors in comparison to NAT and normal prostates. As shown in Figure 4A, we found hsa-miR-133a-3p, 133b and 1-3p signi cantly downregulated in PCa compared to NAT (TCGA-PRAD cohort). Additionally, hsa-miR-133a-3p, 133b and 1-3p expression were signi cantly decreased in PCa compared to normal prostate (Taylor cohort) ( Figure 4B). These results proposed that hsa-miR-133a-3p, miR-133b and miR-1-3p might act as tumor suppressors in PCa.
CENPF and WHSC1 are increased in primary prostate tumors from patients and positively correlated with Gleason grade We found that hsa-miR-133a-3p, 133b and 1-3p showed nine common target genes with negative correlation. Then, we analyzed the expression pattern of the nine target genes in prostate tumors and normal samples from TCGA-PRAD and Taylor (GSE21032) datasets. As shown in Figure 7A, we found that CENPF, WHSC1, CPSF6, RBM15B, PAICS and ELAVL1 expression was signi cantly increased in primary prostate tumors compared to NAT (TCGA-PRAD dataset). Moreover, CENPF, WHSC1 and PAICS expression were signi cantly increased in prostate tumors compared to normal prostates (Taylor dataset, Figure 7B). With these results we analyzed the correlation of CENPF, WHSC1 and PAICS with PCa stage. CENPF and WHSC1 were found to be signi cantly up-regulated according as Gleason score increases in TCGA-PRAD and Taylor datasets ( Figure 8A-B). Moreover, a high expression of these genes correlates with a worst progression free interval ( Figure 8A-B). These results collectively suggest that common target genes of hsa-miR-133a-3p, 133b and 1-3p act as oncogenes favoring PCa aggressiveness.

Discussion
Here we demonstrated that HFD markedly reduces the expression of potential tsmiRs in TRAMP-C1 tumors developed as allograft in C57BL/6J mice. Target genes modulated by these tsmiRs regulate processes mainly associated to cancer related pathways. Among these, mmu-miR-133a-3p, miR-133b and miR-1-3p were dramatically diminished in tumors by HFD. Additionally, human orthologous miRNAs were signi cantly downregulated in human prostate tumors compared to NAT and normal prostates. Besides, a low expression of hsa-miR-133a-3p, 133b and 1-3p correlated with more aggressive pathological features of prostate tumors, such as Gleason score and progression free-interval in PCa patients.
Hsa-miR-133a-3p is probably the most studied miRNA of this family and has been extensively reported as down-regulated in several types of cancer and predicted a poor prognosis [21][22][23][24][25][26][27]. In PCa, a low expression of hsa-miR-133a-3p has been associated with the recurrence and distant metastasis of PCa [28][29][30][31][32]. Likewise, a recent study from Tang et al demonstrated that hsa-miR-133a-3p expression is reduced in PCa tissues compared with the NAT and benign prostate lesion tissues, particularly in bone metastatic PCa tissues. Also, low expression of miR-133a-3p is signi cantly correlated with advanced clinicopathological characteristics and shorter bone metastasis-free survival in PCa patients [33].
As referred to hsa-miR-133b and miR-1, different studies demonstrated that both miRNAs are among the most down-regulated miRNAs in PCa compared non-cancerous prostate tissues and signi cantly decreased in recurrent PCa specimens in comparison to non-recurrent PCa samples [34][35][36]. MiR-1 is further down-regulated in cancer progression and alone can predict disease recurrence [35]. Also, miR-133b and miR-1 have su cient power to distinguish recurrent PCa specimens from non-recurrent [34,35,37]. Lastly, both miRNAs have been reported to target the oncogenic function of purine nucleoside phosphorylase (PNP) in PCa [31]. Thus, both miR-133b and miR-1 may exert similar tumor suppressor activities and coordinately regulate the expression of oncogenes controlling PCa initiation and progression.
Based on our and other groups' ndings, we propose hsa-miR-133a-3p, 133b and 1-3p as a promising miRNA to be studied as potential biomarkers for PCa diagnosis and prognosis.
In this work, we also aim to nd relevant common target genes for hsa-miR-133a-3p, miR-133b and miR-1-3p. We performed a correlation analysis using data from PCa patients available in public algorithms. This analysis allowed us to nd nine relevant common target genes for the three miRNAs with a signi cant negative correlation. Among the nine target genes, CENPF, WHSC1 and PAICS were signi cantly increased in primary prostate tumors compared to NAT and normal prostates, and positively correlated with more aggressive pathological features of PCa as Gleason grade and a worst progressionfree interval. Also, the three target genes were signi cantly enriched in several cancer-related biological processes, suggesting that may act as oncogenes in PCa.
Centromere protein F (CENPF) is a key component of the kinetochore complex and plays a crucial role in chromosome segregation and cell cycle progression. Different studies have investigated the role of CENPF in PCa. Two independent reports demonstrated that CENPF overexpression in PCa tissues is linked to higher Gleason grade, advanced pathological tumor stage, the presence of metastasis, worst overall survival and prostate-speci c antigen (PSA) failure [39,40]. CENPF is also a critical regulator of cancer metabolism potentially through its role in the mitochondria. Shahid et al. reported that CENPF silencing increased inactive forms of pyruvate kinase M2, a rate limiting enzyme needed for an irreversible reaction in glycolysis, and reduced global bio-energetic capacity, Acetyl-CoA production, histone acetylation, and lipid metabolism [41]. Further untargeted metabolomics analysis revealed that silencing of CENPF alters the global metabolic pro le of PCa cells and inhibits cell proliferation [42].
Meanwhile, Wolf-Hirschhorn syndrome candidate 1 (WHSC1/NSD2) is a histone methyltransferase enzyme that targets dimethyl and trimethyl H3K36. WHSC1 has been shown to be up-regulated in prostate tumors and promote malignant growth and progression [43,44]. WHSC1 silencing in DU145 and PC-3 tumor cells decreased cell proliferation, colony formation and strikingly diminished cell migration and invasion [45]. Additionally, pharmacologically inhibition of WHSC1 with MCTP-39 inhibited DU145 prostate tumor growth in vivo [46]. WHSC1 levels positively correlate with the presence of an immunosuppressive microenvironment. Genetic and pharmacological inhibition of WHSC1 restores antigen presentation via an epigenetic regulation of gene expression at the chromatin and DNA methylation levels [44].
Therefore, based on our and previous studies, CENPF and WHSC1 have a critical role in PCa pathogenesis and progression.
In conclusion, our results demonstrated that HFD dramatically reduces the expression of tsmiRs in androgen-sensitive prostate tumors. Additionally, 3 miRNAs: hsa-miR-133a-3p, -133b and − 1-3p are epigenetically silenced by promoter hypermethylation and functions as a tsmiRs in PCa. Moreover, their low expression signi cantly correlated with advanced clinicopathological characteristics of PCa patients, including Gleason score and progression free interval. Besides, the expression of hsa-miR-133a-3p, -133b and − 1-3p negatively correlates with CENPF and WHSC1, two PCa driver oncogenes. Although evaluations of the three miRNAs and their target genes expression in larger populations are still needed, our results indicate that hsa-miR-133a-3p, -133b, -1-3p CENPF and WHSC1 are functional drivers of PCa and may be target for advanced PCa treatment.

Conclusions
HFD modulates the expression of a substantial number of miRNAs in PCa. Attenuation of hsa-miR-133a-3p, -133b and − 1-3p expression by promoter methylation in prostate tumors may enhance PCa development, in part by targeting CENPF and WHSC1. Hsa-miR-133a-3p, -133b, -1-3p, CENPF and WHSC1 are functional drivers of PCa and might be potential targets for additional therapeutic opportunities for lethal PCa.      133b (B) and miR-1-3p (C) in primary PCa tissues from TCGA-PRAD cohort. Patients were categorized into two groups based on the median expression levels of miRNA. Data were analyzed using log-rank Test.