- Research
- Open access
- Published:
Transcriptomic profile of RNA pseudouridine modification as a biomarker for cellular senescence associated with survival outcomes in colorectal cancer
BMC Biology volume 23, Article number: 61 (2025)
Abstract
Background
Colorectal cancer (CRC) is considered as an age-related disease, and cellular senescence (CS) plays a crucial role in cancer development and progression. Previous studies have shown the role of epigenetic changes in aging and cancer development, but the role of RNA pseudouridine (Ψ) modification in aging and cancer remains to be explored.
Results
Using bulk RNA sequencing, CRC cells with low Ψ writers expression levels have higher CS levels. We developed the Psi Score for assessing the transcriptomic profile of RNA Ψ modification regulation and found that the Psi Score correlates with CS. Furthermore, Psi-related senescence may be mediated by mTOR, TGF-β, TNF-α, and inflammatory response signaling pathways. Meanwhile, Psi Score could predict the anti-cancer treatment outcomes of anti-aging interventions and could be used to predict the response to immunotherapy.
Conclusions
Overall, these findings reveal that RNA Ψ modification connected aging and cancer and provided novel insights into biomarker-guided cancer regimens.
Background
Cellular senescence (CS), the state of stable and terminal cell cycle arrest, has been associated with processes such as tumor suppression and organismal aging [1,2,3]. CS can be classified into developmentally programmed senescence, replicative senescence (RS), and stress-induced premature senescence (SIPS) triggered by stimuli including DNA damage, oncogene activation, and oxidative stress [4]. Senescence has been considered as a potent barrier to prevent tumor progression [5, 6]. Therefore, there is increasing interest in finding new markers of senescence that may have prognostic potential for aging and cancer. Notably, epigenetic alterations affect gene expression and various essential cellular processes contributing to senescence-related diseases and cancer [7]. Understanding the mechanisms of these changes is important for the development of new tumor treatment strategies and prevention. However, the relationship between RNA modifications and senescence in cancer remains poorly understood. For instance, whether RNA Ψ modifications regulate the senescence in cancer has not been explored in the field of epitranscriptomics.
Pseudouridine (Ψ) modification is one of the most abundant modifications in RNA epitranscriptomics [8,9,10,11,12]. The incorporation of Ψ into synthetic mRNA plays a significant role in translation, enhancing the efficiency of protein synthesis and improving mRNA stability [13,14,15]. Notably, the ability of Ψ and its derivatives N1-methyl-Ψ in mRNA to evade immune recognition has contributed significantly to the COVID-19 mRNA vaccines [16,17,18,19]. RNA modifications play an important role in biological processes such as CS and aging and influence the progression of aging-related diseases, including cancer [8, 15, 20,21,22]. However, due to the low abundance of Ψ in mRNA and the limitations of Ψ modification techniques, further studies are needed to investigate the role of Ψ modification in senescence in cancer.
Colorectal cancer (CRC) is the third most common type of cancer and the second leading cause of cancer-related death globally [23, 24]. Both genomic instability and epigenetic alterations have been identified as hallmarks of aging and contribute significantly to CRC pathogenesis [7, 25, 26]. Thus, identifying more genetic and epigenetic loci associated with CRC risk is essential. Studies have shown that the role of pseudouridine synthases (PUS) in RNA Ψ modification is multifaceted and is not only involved in RNA stabilization and processing, but may also affect gene expression and tumor development [27, 28]. It has been shown that three of these enzymes are associated with cancer and influence the phenotype of cancer cells [28]. Nevertheless, the biological functions of PUS in CRC remain undefined.
Here, we constructed the Psi Score to quantify the regulation of RNA Ψ modification in CRC transcriptome and found a strong association between Psi Score and CS. We then showed the RNA Ψ modification is involved in the transcriptional alterations underlying CS in the tumor development by constructing the Psi Score model. Our findings highlight the potential relationship between RNA Ψ modification and CS, offering new insights into the identification of biomarkers and treatment related to senescence-based anti-cancer.
Results
Genetic and transcriptional alterations of 13 RNA Ψ modification writers in CRC
Since CRC carcinogenesis is frequently associated with genome instability and susceptibility to mutations [29, 30], we examined 13 RNA Ψ modification writers for somatic mutations in CRC to determine the prevalence of genetic alterations (Additional file 1: Table S1). Among 455 COAD samples, 66 (14.51%) exhibited mutations in the RNA modification writers (Fig. 1A). In CRC samples, PUS7 exhibited the highest mutation frequency (4%), followed by PUS3, PUS10 and PUS7L, while RPUSD3 and PUSL1 were not present. Our cohort showed multiple pairwise mutation correlations (Fig. 1B). As anticipated, there was significant co-occurrence among RNA Ψ modification writers with no mutual exclusion. Analysis of hallmark pathway gene signatures [31] highlighted that the mutation group was significantly enriched not only in carcinogenic pathways but also in inflammatory pathways, such as KRAS, TNF-α, TGF-β, IL6-JAK-STAT3, and inflammatory response signaling pathways (Fig. 1C). As a pleiotropic cytokine produced by macrophages, IL-6 is involved in immune and inflammatory responses via the JAK3/STAT3/SOCS3 pathway. Furthermore, CRC is closely related to dysregulation of the IL-6-mediated JAK/STAT3/SOCS3 signaling pathway [32, 33]. It appears that mutations in RNA Ψ modification writers may influence the survival prognosis of CRC patients. Mutations in RNA Ψ modification writers may accompany biological alterations through changes in the above pathways.
Genetic and transcriptional alterations of RNA Ψ modification writers in CRC. A The mutation frequency of 13 RNA Ψ modification writers in 455 CRC patients from the TCGA cohort. Each column represents individual patients. The upper bar graph shows TMB; the number on the right indicates the mutation frequency in each “writer”. The right bar graph shows the proportion of each variant type. B The pairwise association of Ψ writer mutations in CRC patients is shown, with the degree of association indicated by the color of the box. Blue color indicates co-occurence and red color indicates a mutually exclusive pattern. P < 0.05 was considered statistically significant. C Differences in pathway activities scored by GSVA between mutation and non-mutation CRC samples. Shown are t values from a linear model; blue: mutation; red: non-mutation. D Box plots show the expression distribution of 13 writers of RNA Ψ modification between paired normal (blue) and CRC (red) tissues. (* p < 0.05, ** p < 0.01, *** p < 0.001)
Furthermore, in order to elucidate the potential impact of genetic mutations on the transcriptional regulation of genes and their role in cancer progression, we examined the mRNA alterations of Ψ writers between paired CRC and normal tissue samples and found that all writers except PUS10 were significantly increased in CRC (Fig. 1D). Consistent with previous reports that upregulation of PUS expression in CRC may cause cancer development, and downregulation of DKC1 in lung adenocarcinoma (LUAD) cells induces senescence has been reported previously [34]. Together, we hypothesized that the expression of Ψ modification writers is relevant to CRC development and cellular senescence.
RNA Ψ modification pattern is associated with cancer hallmarks
For further analysis of the expression patterns of the Ψ writers involved in tumorigenesis, 1293 CRC patients were selected (Additional file 1: Table S2). Based on pairwise correlations between the expression of 13 Ψ writers in CRC, we observed a more frequent presence of positive correlations than negative correlations (Fig. 2A). Notably, the expression of DKC1 was positively correlated with other writers, consistent with significant mutational co-occurrence between DKC1 and others (Fig. 1B). Meanwhile, DKC1 showed significant variability between paired normal and tumor samples, suggesting that it may play a unique and important role in tumor development (Fig. 1D).
Patterns of RNA Ψ modification and biological characteristics of each pattern and the construction of Psi Score. A The heatmap displays the positive (red) and negative (blue) correlations among RNA Ψ modification writers in CRC, as determined by Spearman correlation analysis. B Kaplan–Meier curves depict overall survival differences between two distinct RNA Ψ modification patterns, Cluster 1 (pink) and Cluster 2 (orange), based on the GSE39582 dataset. The grouping of CRC samples is shown at the bottom of the chart. p < 0.05 in the two-sided log-rank test was considered statistically significant. C A heatmap visualizing the GSVA enrichment analysis shows the activation states of biological pathways in distinct RNA Ψ modification patterns. Red, activated pathways; blue, inhibited pathways. D Unsupervised clustering of 13 RNA Ψ modification writers. The names of CRC cohorts were used as sample annotations. Red, high expression of DEGs; blue, low expression. E Overall methodology. Workflow for integrative analysis of the Ψ modification levels in colorectal cancers using GEO transcriptome data
To classify the CRC patients from the integrated GEO-CRC cohorts with qualitatively different RNA modification patterns (Additional file 1: Table S3), we used consensus clustering based on the expression profiles of 13 RNA Ψ modification writers. Cluster 1 consisted of 490 CRC patients, while Cluster 2 comprised 803 patients, with the latter demonstrating a significant survival advantage (Fig. 2B; log-rank test, p = 0.0173). GSVA enrichment analysis emphasized the biological significance of these distinct patterns, with Cluster 1 enriched in p53, TNF-α, and interferon-gamma response pathways, and Cluster 2 exhibiting enrichment in DNA repair, G2M checkpoint, MYC targets, and E2F targets pathways (Fig. 2C).
We next further investigated the functional significance of these two Ψ modification patterns by enrichment analysis of cluster-specific differentially expressed genes (DEGs). These DEGs were enriched in biological processes such as cell cycle, phagosome, and DNA replication (Additional file 2: Fig. S1). They were also enriched in signaling pathways, especially the TGF-beta and BMP signaling pathways (Additional file 2: Fig. S2). We performed unsupervised clustering analyses based on the cluster-specific DEGs to categorize patients into geneClusterA and geneClusterB. The clustering of DEGs corresponded with the Ψ modification patterns, resulting in the classification of 837 out of 1,293 CRC patients into geneClusterA, associated with Cluster2, and 456 into geneClusterB, associated with Cluster1 (Fig. 2D).
Evaluation of RNA Ψ modification in transcriptomic regulation with Psi Score
In order to quantify the RNA modification patterns of individual CRC patients, we developed the Psi Score to quantify Ψ modification levels by analyzing differential gene expression between two Ψ clusters and calculating the expression activity of DEGs (Fig. 2E). As anticipated, Cluster1 exhibited significantly higher Psi Score than Cluster2 (Additional file 2: Fig. S3; Wilcoxon test, p < 2.2 × 10−16), and geneClusterB showed a significantly higher Psi Score than geneClusterA (Additional file 2: Fig. S3; Wilcoxon test, p < 2.2 × 10−16).
In order to verify the differential regulation, we analyzed the overlap between these three different classifiers. A significant overlap was observed between the Psi Score-high group and geneClusterB, with 456 samples (70.59%) found in both categories (Additional file 2: Fig. S4). Conversely, geneClusterA encompassed all 647 samples (100%) of the low Psi Score group. Additionally, 379 out of 646 (58.69%) CRC samples in the Psi Score-high group were associated with Cluster1, and a significant majority, 536 out of 647 (82.84%), in the Psi Score-low group were associated with Cluster2 (Additional file 2: Fig. S4). Collectively, these findings indicated a high degree of concordance among the three classification methods.
RNA Ψ modification writer expression is associated with senescence level
According to the results presented above, we speculated that tumors with low expression of Ψ modification writers exhibit a high level of senescence. To quantify senescence levels, we referred to the methods provided in the previous study [35]. We calculated the CS score by utilizing gene sets from the CellAge database [36] and applying the GSVA algorithm [31] to measure the differential activity of CS-related genes. We investigated the effect of RNA Ψ modification writers on senescence levels and revealed both significant correlations and notable interactions among these writers. The results indicated that tumors with low expression levels of Ψ modification writers were consistent with high Psi Score, while high expression levels were associated with low Psi Score (Fig. 3A).
RNA Ψ modification writers associated with senescence level. A Expression level of 13 RNA Ψ modification writers. The names of CRC cohorts were used as sample annotations. Red, high expression of writers; blue, low expression. B Senescence level in Ψ modification writers DKC1 knockdown datasets. NC: negative control; sh: shRNA. C The Venn diagram depicted the intersection of gene sets that construct the two scoring systems (Psi Score and CS Score)
We next constructed WiDr stable cell lines featuring DKC1 knockdown and conducted transcriptome analysis (GSE271639) through RNA-seq. Our results indicated that WiDr stable cell lines with knockdown of the DKC1 gene exhibited higher senescence levels than the control group. Furthermore, analysis of two additional Ψ modification gene knockdown datasets from the Gene Expression Omnibus (Additional file 1: Table S2) revealed a significant increase in senescence levels in tumor samples deficient in these genes. These findings highlight the potential of Ψ modification writers as crucial regulators of senescence and provide further evidence that supports our conclusions. (Fig. 3B and Additional file 2: Fig. S5). Based on these results, we hypothesized that the low Ψ modification writers’ expression levels might induce senescence.
Subsequently, we assessed the relationship between the Psi Score and CS Score to identify common genes. The analysis identified only 11 overlapping genes (Fig. 3C), including CDKN2B (Additional file 2: Fig. S6), which is a recognized marker of cellular senescence [58]. Elevated expression of CDKN2B was observed in the Psi Score-high group, indicating a correlation between Psi Score and levels of senescence.
Overall, the negative correlation of Ψ modification writers’ expression with CS Score reveals that the changes in Ψ modification writers’ expression promote cellular senescence.
Psi Score is associated with cellular senescence in CRC
After exploring the senescence landscape associated with Ψ modification writers, we hypothesized that Ψ writers-mediated genes may regulate senescence. Next, we applied our approach to systematically profile the levels of senescence and Ψ modification in 1293 CRC samples from the GEO database (Additional file 1: Table S3). We revealed the levels of senescence and Ψ modification in different CRC types. Next, we examined the correlation between the Psi Score and the CS Score. A significant positive correlation was observed between Psi Score and CS Score (correlation coefficient = 0.74; Spearman rank correlation, P < 2.2 × 10−16) (Fig. 4A). Additionally, we found that Psi Score had a significantly positive correlation with positively CS-correlated genes (rho = 0.84, p < 2.2e − 16) and a negative correlation with negatively CS-correlated genes (rho = − 0.36, p < 2.2e − 16) (Additional file 2: Fig. S7), suggesting that positively CS-correlated genes largely mediated the association between Ψ modification levels and senescence.
Construction of the Psi Score and Senescence landscape associated with Psi Score in CRC. A The correlation coefficients of the Psi Score and CS Score were calculated by using Pearson’s correlation coefficient formula. The red scatter points represent CRC samples and the blue line indicates the correlation. B, C, and D Comparison of GSVA scores and GSEA for characteristics associated with multiple pathways of senescence between Psi Score-high and -low groups. Cellage Induces: CS-related positive genes in CellAge database; Cellage Inhibits: CS-related negative genes in CellAge database; Stress induced: the CS curated gene sets of stress-induced senescence in CellAge database; Aging Atlas: The CS curated gene sets in Aging Atlas database; GenAge: The CS curated gene sets in GenAge database; SenMayo and SenoRanger: The CS curated gene sets collected from the published literature; Stress induced: Stress-induced senescence genes in CellAge database; Oncogene induced: oncogene-induced senescence genes in CellAge database; Replicative: replicative senescence genes in CellAge database; Pink, Psi Score above median; orange, Psi Score below median. E Heatmap shows the correlation between Psi Score and CS Score in five GEO-CRC datasets and the 50 Hallmark signaling pathways according to the Pearson’s correlation coefficient formula. Orange, high correlation of Psi/CS Score and hallmark pathways; purple, low correlation. * p < 0.05, ** p < 0.01, *** p < 0.001
To further validate the results above, gene set enrichment analysis for CS-curated gene sets collected from the published literature was performed [37,38,39,40,41]. First, we evaluated the pathways enriched by the CS-correlated genes collected from CellAge database [36], which showed that the CS-correlated genes were differentially expressed between the Psi Score-low and Psi Score-high groups (Fig. 4B). Positive CS-correlated genes were significantly enriched in Psi Score-high groups. In contrast, negative CS-correlated genes were significantly enriched in Psi Score-low groups. Among them, stress-induced senescence signaling pathways were notably enriched in Psi Score-high groups (p < 2.2e − 16) (Fig. 4C). Additionally, gene set enrichment analysis (GSEA) showed significant enrichment of senescence-related pathways in Psi Score-high groups (Fig. 4C). Furthermore, gene sets from other databases such as GenAge, SenMayo, SenoRanger, and Aging Atlas were obtained from published literature on the hallmarks of senescence [37,38,39, 41]. As anticipated, it was found that CS genes were significantly enriched in Psi Score-high groups, indicating a significant correlation with RNA Ψ modification levels (Fig. 4D). GO biological processes enrichment analysis also showed similar trends with CS-correlated genes predominantly enriched in Psi Score-high groups (Additional file 2: Fig. S8). The heatmap illustrating the correlation of the Psi Score and CS score with Hallmark signaling pathways indicated a consistent correlation between the two scores, reflecting the high correlation between Ψ writers-specific genes-mediated pathway and senescence-mediated pathway (Fig. 4E).
To evaluate the correlation between Psi Score and CS Score, 1293 CRC samples were divided into four groups based on the median Psi Score and a threshold of 0 for the CS Score (Fig. 5A). Survival analyses revealed that patients with cellular senescence suffered poorer OS, with the worst prognostic outcomes in patients with Psi-related cellular senescence (p = 0.0376) (Fig. 5B). To explore the underlying cause of poor outcomes with senescence groups, we performed GSVA enrichment analysis (Fig. 5C). The signaling pathways activated mostly in the Psi-related cellular senescence groups were WNT, EMT, and KRAS signaling pathways, while the characteristic signaling pathways activated in Psi-independent cellular senescence groups were those related to the p53, DNA repair, G2M checkpoint, MYC target and E2F target signaling pathways. Notably, the signaling pathways involved in cellular senescence groups (Psi-related CS and Psi-independent CS) were mainly enriched in mTOR, TGF-β, TNF-α, and inflammatory response signaling pathways. These pathways are involved in tumor development and inflammation in cellular senescence. Thus, we speculated that RNA Ψ modification may mediate these pathways to regulate senescence. Furthermore, to identify whether age is associated with the Psi Score, CRC samples were grouped by age 50 (EOCRC and LOCRC) (Additional file 2: Fig. S9). No significant association was found between Psi Score and age, and the close correlation between Psi Score and CS Score was independent of age.
Association of Psi Score and Cellular Senescence with prognostic outcomes and biological pathway activation in CRC. A The four-quadrant diagram of the Psi Score and CS Score. Psi-independent CS: Psi-independent cellular senescence. Psi-related CS: Psi-related cellular senescence. Psi-related non-CS: Psi-related non-cellular senescence. Psi-waived non-CS: Psi-waived non-cellular senescence. B Kaplan–Meier curves compare overall survival for the four combined CS and PSi score subgroups in GSE39582. The grouping of CRC samples is shown at the bottom of the chart. p < 0.05 in the two-sided log-rank test was considered statistically significant. C A heatmap visualizing the GSVA enrichment analysis shows the activation states of biological pathways under four subgroups of CS and Psi Score. Red, activated pathways; blue, inhibited pathways
Moreover, these observations were consistent with the current view that abnormal activation of oncogenes and inactivation of tumor suppressors can promote CS [42]. Significant negative correlation between GSVA enrichment scores and Psi Score scores for DNA repair, G2M checkpoint, MYC target, and E2F target signaling pathways, and a significant positive correlation between P53, TNF-α, TGF-β, and BMP pathways and Psi Score (Fig. 4E and Additional file 2: Fig. S10). The above results revealed that Ψ modification levels and cross-talk among Ψ writers may be involved in senescence.
Molecular subtypes and clinical characteristics associated with Psi Score in CRC
Building on our previous findings that the Psi Score is associated with poorer prognosis in CRC patients due to its link with cellular senescence, we next investigated the relationship between the Psi Score and CRC molecular subtypes, as well as their associated clinical characteristics [43, 44].
In our analysis of CRC samples from integrated GEO cohorts and the TCGA Cohort, we observed significant variations of the Psi Score among these CMS subtypes (Additional file 1: Table S3, S5, and S6). Notably, the CMS4 subtype exhibited the highest Psi Score (Fig. 6A–B). TGF-β and VEGF signaling pathways were significantly enriched in the Psi Score-high group. In contrast, the Psi Score-low group considerably enriched signaling pathways related to the cell cycle and ribosomes (Fig. 6C).
Psi Score is associated with molecular phenotypes and clinical features in CRC. A and B Psi Score differences among CMS subtypes of CRC in integrated GEO-CRC datasets (A) and the TCGA-COAD/READ cohort (B). C Gene set enrichment analysis (GSEA) shows the differences in enrichment in the characteristic signaling pathways of CMS subtypes between high- and Psi Score-low groups in integrated GEO-CRC datasets. D Psi Score differences among TNM stages of CRC in GEO-CRC datasets. Wilcoxon test was used to assess the difference. E Differences in the Psi Score between BRAF MT (purple) and BRAF WT (pink) group in integrated GEO-CRC datasets. F Differences in the Psi Score among PDS subtypes of CRC in integrated GEO-CRC datasets
The distinct subtypes of CMS are significantly associated with tumor progression and can be utilized to predict clinical prognosis [45]. Psi Score varied across tumor stages and was higher in advanced CRC (Fig. 6D), suggesting that the parameters that comprise the Psi score are involved in tumor progression. In addition, the Psi Score was also significantly higher in patients with BRAF mutant compared to those with BRAF wild-type (Fig. 6E).
In a recent study, Malla et al. conducted PDS subtyping of CRC using gene ontology and biological activation state information [46]. Similarly, we observed significant differences in Psi scores across PDS subtypes, with the PDS2 subtype having the highest Psi Score (Fig. 6F and Additional file 1: Table S7).
These findings indicated a strong correlation between the Psi Score and CMS subtypes and PDS subtypes, potentially establishing a link between the Psi Score and tumor progression in CRC patients.
Psi Score as a biomarker of cellular senescence in response to anti-aging interventions
Next, we aimed to track the dynamic change of cellular senescence with Psi Score during aging-related interventions. Metformin, which is a widely prescribed medication for type 2 diabetes, has garnered attention for its potential anti-aging properties [47,48,49,50]. It is known for modulating gene expression profiles similar to caloric restriction, a well-established anti-aging regimen [51, 52].
We analyzed datasets with metformin-treated and control groups from GEO databases (Table S2). Compared to control groups, metformin-treated groups had a significantly lower Psi Score (P = 0.025; Fig. 7A). Similarly, we observed a consistent distribution of Psi Score across colorectal cancer cell in GEO datasets (Additional file 1: Table S2) under nutrient-depleted conditions. The results showed that the Psi Score significantly decreased under nutrient-depleted conditions, consistent with the lowered senescence levels observed under caloric restriction (P < 0.001; Fig. 7A) [51]. In contrast, nutrient-rich conditions were associated with a higher Psi Score. The modulatory effect of metformin on nutrient-sensing pathways, which play a crucial role in the anti-aging benefits of caloric restriction, has been extensively researched [52, 53]. Therefore, we hypothesized that the changes in Psi Score with metformin may be related to caloric restriction. Evidence shows that metformin could emulate the benefits associated with dietary restriction through direct inhibition of mTORC1 and sirtuins, with potential effects downstream of caloric restriction.
Psi Score as a biomarker for treatment response in CRC. A Psi Score change in CRC cells treated with 2.5 mM metformin (GSE171478) and under nutrient-rich and nutrient-depleted conditions (GSE245402). Differences between two groups are tested using two-tailed paired t test and considered p < 0.05 as statistical significance. B The correlation between Psi Score and drug sensitivity evaluated by the Spearman analysis. Each column represents a drug. The brightness of the column indicates the significance of the correlation. The height of the column indicates the correlation, indicates that Psi Score related to drug resistance (Rs > 0) or drug sensitive (Rs < 0) to Psi Score. C Signaling pathways targeted by drugs that are resistant or sensitivity to the Psi Score. Drug names are listed on the horizontal axis and the signaling pathway targeted by the drug on the vertical axis. The bar graph on the right shows the number of drugs targeting each signaling pathway. The size of the point indicates the significance of the correlation
In conclusion, our findings suggest that the Psi Score, a marker of cellular senescence, is reduced following treatment with metformin or under caloric restriction conditions. This indicates that the Psi Score could potentially serve as a bioindicator of the effectiveness of interventions such as metformin and CR in combating senescence.
Potential clinical application value of the Psi Score
To further elucidate the impact of the Psi Score on drug response, we analyzed the relationship between the Psi Score and drug sensitivity in cancer cell lines using data from the GDSC database. The GDSC database revealed a significant correlation between 53 drugs and Psi Score (Fig. 7B and Additional file 1: Table S8). We found that the sensitivity to 23 drugs was correlated with the Psi Score. For instance, IGF-1R inhibitor BMS-754807 (Rs = − 0.74, p = 3.65E − 251), mTOR inhibitor AZD8055 (Rs = − 0.51, p = 5.83E − 88), and PKC inhibitor staurosporine (Rs = − 0.43, p = 3.28E − 35) exhibited correlations with the Psi Score indicating drug sensitivity. On the other hand, resistance to 30 drugs correlated with Psi Score, including ERK1/ERK2 inhibitor ulixertinib (Rs = 0.51, p = 1.53E − 92), Bcl-2 inhibitor navitoclax (Rs = 0.42, p = 2.55E − 63), and EGFR inhibitor AZD3759 (Rs = 0.30, p = 3.26E − 32). Further analysis revealed that drugs showing sensitivity associated with high Psi Score, primarily targeting RTK, WNT, and mTOR signaling pathways. In contrast, drugs associated with low Psi scores targeted the EGFR, MAPK, and cell cycle signaling pathways (Fig. 7C).
Notably, 14 of these 23 drugs have been documented for their tumor-inhibitory effects and role in senescence regulation. Among them, AZD8055 has been found to exhibit anti-tumor and anti-senescence properties by targeting the mTOR pathway [54, 55] (Additional file 3: Table S9). These results suggest that RNA Ψ modification patterns are correlated with drug sensitivity. Therefore, the Psi Score may be a potential biomarker for defining appropriate treatment strategies.
Predicting immunotherapy response with the Psi Score and PD-L1 blockade
Building upon our previous investigation into the relationship between the Psi Score and drug sensitivity, we assessed the potential of the Psi Score as a predictive biomarker for immunotherapy response, specifically focusing on the PD-L1 pathway. Additionally, we considered whether the Psi Score, informed by these patterns, could predict the efficacy of immune checkpoint blockade (ICB) therapy.
Correlative studies have suggested that PD-L1 protein expression and tumor mutation burden may be potential biomarkers for predicting response to immune checkpoint inhibitors [56]. Numerous studies have also identified biomarkers that can predict response to immunotherapy [57,58,59]. Therefore, we evaluated the ability of the Psi Score to predict patient responses to ICB therapy, given its possible association with Ψ modification patterns and tumor immune microenvironment (Additional file 2: Fig. S11). According to TIDE analysis, a higher TIDE score indicates an increased risk of immune escape and reduced benefit from immunotherapy. As expected, groups with high Psi Scores showed higher “TIDE scores,” “Exclusion scores,” and “Dysfunction scores” compared to groups with low Psi Scores (Fig. 8A–C).
Association of Psi Score with immune response and clinical outcomes in the IMvigor210 Cohort. A, B, and C TIDE analysis including TIDE score (H), exclusion score (I), and dysfunction score (D). The difference in the Psi Score between distinct clinical outcomes of anti-PD-L1 treatment in the IMvigor210 cohort. SD, stable disease; PD, progressive disease; CR, complete response; PR, partial response. E Differences in TMB between Psi Score-high and -low groups in the IMvigor210 cohort. Wilcoxon test was used to assess the difference. F The proportion of patients in the IMvigor210 cohort with different responses to PD-L1 blockade immunotherapy. The chi-squared test was used to determine the statistical significance of the difference. SD, stable disease; PD, progressive disease; CR, complete response; PR, partial response
In the IMvigor210 cohort, patient responses were categorized as complete response (CR), partial response (PR), stable disease (SD), and progressive disease (PD). A total of 348 patients were analyzed, and results indicated that patients with partial response demonstrated lower Psi Scores (Fig. 8D and Additional file 1: Table S10). Furthermore, analysis of anti-PD-L1 therapy outcomes revealed that patients with low Psi Scores experienced significantly improved treatment outcomes (Fig. 8E; P = 0.01014). Additionally, TMB was significantly higher in the group with low Psi Score than those with high (Fig. 8F and Additional file 1: Table S10).
Overall, our analysis shows that the Psi Score can be used to predict the response to anti-PD-L1 immunotherapy and select drugs for CRC.
Discussion
The RNA Ψ modification seems to be involved in the transcriptional alterations underlying cellular senescence in CRC development. Consistent with previous reports, upregulated expression and mutations of Ψ writers promote activation of oncogenic signaling pathways leading to the development of CRC (Fig. 1C, D). The low expression of Ψ writers reveals higher levels of senescence under the context of DKC1 knockdown (Fig. 3B), consistent with previous literature indicating that downregulation of DKC1 in lung adenocarcinoma (LUAD) cells inhibits cell proliferation and shortens telomere, which ultimately leads to telomere-related cellular senescence [34]. Meanwhile, there is a significant positive correlation between Psi Score based on Ψ writers-mediated genes and senescence levels (Fig. 4A–D). This is likely due to Ψ modification activating mTOR, TGF-β, TNF-α, and inflammatory response signaling pathways to regulate senescence (Fig. 5C). In addition, more than half of the drugs associated with the Psi Score had anti-tumor and anti-senescence effects through searching the literature (Additional file 3: Table S9). Thus, our findings offer a potential biomarker for treatment targets related to senescence and cancer by developing the Psi Score.
The relationship between cancer and aging is bilateral. Epigenetic alteration and cellular senescence are both common hallmarks of cancer and aging. Mounting evidence has demonstrated that the RNA modifications are related to cancer invasion. For example, DKC1 has been shown to promote CRC progression, and high DKC1 expression in CRC is associated with poor prognosis [60]. PUS1 has been shown to contribute to hepatocellular carcinoma and breast cancer through mRNA pseudouridylation [61,62,63]. Moreover, PUS7 has been found to promote glioblastoma progression, and chemical inhibitors of PUS7 reduce the level of Ψ modification in glioblastoma and inhibit tumor growth [64]. Notably, DKC1 downregulation has been shown to promote telomere-related cellular senescence in LUAD, and high DKC1 expression in LUAD is associated with poor prognosis [34]. However, no data are available on the senescence levels of PUS in CRC.
Here we identified DKC1 knockdown in CRC cell lines associated with higher senescence levels (Additional file 1: Table S2). In addition, DKC1 downregulation in cutaneous squamous cell carcinoma cells and PUS1 knockdown in breast cancer are also associated with higher levels of senescence by screening public databases (Additional file 2: Fig. S5). The high expression of PUS contributes to tumor progression and the downregulation of PUS exhibits higher senescence levels, in agreement with the previous notion that senescence has been shown to prevent tumor progression. Importantly, Psi-related CS subgroups are enriched in signaling pathways involved in tumorigenesis and inflammation (Fig. 5C) and also exhibit the worst prognosis (Fig. 5B), suggesting that RNA Ψ modification may regulate senescence through pathways such as mTOR, TGF-β, and TNF-α signaling pathways. In addition, the Psi Score is independent of onset age (Additional file 2: Fig. S9), which may account for the notion that RNA Ψ modification exerts its effects on senescence through these pathways, irrespective of age-related factors.
Our results show that RNA Ψ modification is involved in the regulation of cellular senescence through pathways such as mTOR, TNF-α, and TGF-β pathways, which are involved in inflammatory responses, cellular senescence, and inflammaging [65,66,67,68,69,70]. In this context, RNA Ψ modifications might regulate genes involved in these pathways, modulating the inflammatory response and potentially influencing both cellular senescence and tumorigenesis. Interestingly, we observed significant differences of Psi Score in different PDS subtypes that were generated based on biological annotation and expression profile of genes in Malla et al.’s study (Fig. 6F and Additional file 1: Table S7) [46]. The PDS2 subtype exhibited the highest Psi Score, which has elevated stromal and immune-tumor microenvironment lineages. Malla et al. suggested that the PDS2 subtype was predominantly CMS4 or CMS1 [46], which is consistent with our observation that Psi Score were highest in the CMS4 and PDS2 subtypes. Notably, the PDS2 subtype was enriched in inflammatory and immune signaling pathways, such as transforming growth factor β (TGF-β) activation, which is consistent with the pathways enriched in the Psi Score-high group and Psi-related CS group (Figs. 5E and 6C). Our findings provide new insights into the regulation of RNA Ψ modifications on cellular senescence.
In order to track the Psi Score in aging interventions, we found that the Psi Score is significantly decreased under metformin-treated and nutrient-depleted conditions compared to the control groups (Fig. 7A). This result is consistent with metformin and caloric restriction having anti-aging effects. Additionally, we demonstrated the potential therapeutic effects of Psi Score in CRC, and identified 53 drugs with Psi Score that were associated with drug sensitivity for targeting RTK, mTOR, and WNT signaling pathways or with drug resistance for targeting EGFR, MAPK and cell cycle signaling pathways (Fig. 7C). We then consulted the literature for these drugs regarding the anti-tumor effects and senescence-related effects (Additional file 3: Table S9). There are 27 drugs with tumor-inhibitory effects and senescence-regulating effects. AZD8055 is reported to be anti-tumor and anti-senescence by targeting the mTOR pathway [54, 55], which is also consistent with our analysis. Patients with higher Psi Score may benefit from drugs targeting these signaling pathways, according to these results. Psi Score may serve as a predictor of the clinical efficacy of chemotherapy or targeted therapy in cancer.
Although this study presents some important findings, there are some limitations. Firstly, we defined Ψ modification levels based on a GSVA-based method that relied on DEGs responsible for RNA Ψ modification patterns. The dependency on predefined gene sets, which might exclude relevant genes not accounted for in these collections, thus limiting the comprehensiveness of the analysis. In addition, the method is susceptible to biases arising from the determination of differential expression thresholds, which can introduce inaccuracies in the quantification of Ψ modification levels. Thirdly, we only observed associations of Psi Score with anti-cancer treatment in CRC. These limitations are somewhat inevitable given the current state of bioinformatics tools and the inherent complexity of biological systems. Further research is needed to verify these findings in larger samples and other cancer types. In addition, the Psi Score is based on a high-throughput transcriptomic profile, which may limit its use in clinical settings. Further studies generating a smaller gene panel could facilitate its integration into routine clinical workflows for guiding personalized treatment. Additionally, there is also a need to develop more sophisticated algorithms capable of integrating a wider range of genomic data, including non-coding RNAs and epigenetic modifications. Finally, for a comprehensive understanding of the relationship between RNA Ψ modification levels and senescence levels, more experimental validation is needed.
Conclusions
We performed a systematic and comprehensive analysis of the transcriptomic regulation of RNA Ψ modification in CRC. The gained regulation of RNA Ψ modification is associated with CS, which is mediated by mTOR, TGF-β, and TNF-α signaling pathways. Moreover, Psi Score, based on the quantification of RNA Ψ modification levels, could predict senescence-based anti-cancer treatment outcomes. Overall, our findings highlight the potential of RNA Ψ modification as a biomarker for anti-cancer treatment targeting senescence, which may provide novel insights into CRC therapy.
Methods
Data and resources
The integrated dataset consisted of gene expression data from the same sequencing platform as well as clinical and molecular annotations retrieved from the Gene Expression Omnibus (GEO) and the Cancer Genome Atlas (TCGA) database. The analysis encompassed five CRC cohorts and the TCGA-COAD/READ cohort, as detailed in Table S2. Based on the TCGA database, somatic mutations, and clinical characteristics, such as age, gender, tumor stage, and overall survival (OS) were collected. The mutation data of COAD was further analyzed and visualized using the R package “maftools.” Table S2 provides a summary of the dataset information. To harmonize the data, the ComBat algorithm from the sva package (version 3.50.0) [71] used batch effect correction. This algorithm operates within an empirical Bayesian framework. We used the log-transformed expression data as input for ComBat, and the batch effects were corrected by applying ComBat’s empirical Bayes approach to adjust for systematic variations across batches. This process involved using a full model matrix that excludes batch variables, thereby mitigating the impact of batch effects on the integrated dataset. Data analysis was conducted using R (version 4.3.1).
The dataset on immunotherapy was obtained from the IMvigor210 cohort and is accessible through the IMvigor210 package, which contains data from a previous study [72]. The study of the correlation between Psi Score and immunotherapy prognosis involved analysis of the IMvigor210 cohort.
Clustering expression pattern of 13 RNA Ψ modification writers
A total of 13 proteins, including DKC1, PUS1, PUS3, PUS7, PUS10, PUSL1, PUS7L, TRUB1, TRUB2, RPUSD1, RPUSD2, RPUSD3, and RPUSD4, have been identified as the writers for RNA Ψ modification [73, 74]. The unsupervised clustering analysis was applied to categorize 1293 CRC patients based on the expression of RNA Ψ modification writers. Consensus Cluster Plus (version 1.66.0) was utilized for the above steps based on the k-means algorithm and repeated 1000 times to ensure stability [75, 76].
Gene Set Variation Analysis (GSVA)
To investigate the impact of RNA Ψ modification on biological processes, we conducted Gene Set Variation Analysis (GSVA) using the GSVA R package [31]. We utilized gene sets from the MSigDB repository, which had been annotated for functional profiling with the “clusterProfiler” R package [77].
Construction of Psi Score
-
1)
Identification of DEGs between RNA Ψ modification patterns. The CRC samples were grouped into two Ψ clusters based on the expression levels of Ψ writers. Subsequently, differential gene expression analysis was conducted between the two Ψ clusters with the Limma R package (version 3.58.1) [78]. The voom function [79] was applied to estimate the mean–variance relationship, followed by linear modeling, and the Empirical Bayesian statistics were then used to detect DEGs between the two RNA Ψ modification patterns, as previously described. The DEGs were selected based on the p-value (< 0.05) and absolute value of logFC (> 1).
-
2)
Construction of Psi Score. We utilized GSVA with Ψ-specific DEGs to calculate the Psi Score for evaluation on the Ψ-dependent transcriptomic profile in each individual sample. To investigate the potential biological role and clinical significance of the Psi Score, CRC samples were categorized into two groups according to the median value of the Psi Score.
Cellular Senescence (CS) Score calculation
Cell Senescence genes were obtained from the CellAge database of the Human Ageing Genomic Resources [80]. The CS signature comprised 1259 genes, categorized into 525 positively and 734 negatively associated with cell senescence [36]. GSVA was utilized to determine the functions of two sets of CS-related genes in each sample, separately. CS scores were derived by subtracting negative CS-related activities from positive CS-related activities.
Association of the Psi Score with senescence
The Pearson correlation test assessed the correlation between the Psi Score and CS Score in the 1293 CRC samples. Correlations were categorized as follows: |R|< 0.2, very weak; 0.2 ≤|R|< 0.4, weak; 0.4 ≤|R|< 0.6, moderate; 0.6 ≤|R|< 0.8, strong; |R|≥ 0.8, very strong.
To analyze the association of Psi Score with CS signature, the CS curated gene sets were collected from the published literature: CellAge database of the Human Ageing Genomic Resources [36], GenAge database of the Human Ageing Genomic Resources [41], Aging Atlas database [39], SAUL_SEN_MAYO [37], SenoRanger [38], and GeneCards [81]. The CS signature across different Psi Score groups was evaluated using GSVA and represented as a boxplot for CRC samples. Statistical significance was determined using the Wilcoxon test with Benjamini and Hochberg correction for the false discovery rate (FDR), under a threshold of FDR < 0.05.
To further investigate the changes in Psi Score resulting from anti-aging interventions, we included CRC cells treated with 2.5 mM Metformin under both nutrient-rich and nutrient-depleted conditions for analysis (Table S2). Table S2 provides a summary of the dataset information. The effects of different treatments on Psi Score groups were assessed using GSVA and depicted as boxplots. Differences between the two groups were tested using a two-tailed paired t-test, with statistical significance considered at p < 0.05.
Construction of DKC1 knockdown cell lines
Lentiviral transduction was used to generate DKC1 knockdown in the WiDr cell line. For lentivirus production, shRNA targeting DKC1 (shDKC1) plasmids and their control plasmids were purchased from IGEBIO (IGE Biotechnology Co., Ltd. Guangzhou. China). Next, plasmids were transfected into the 293 T cell line. The lentiviral supernatant was collected after filtration and added to the WiDr cell line along with polybrene. Stable cell lines were selected using puromycin.
Transcriptomic profiling for DKC1-knockdown colorectal cancer cells
Total RNA was extracted using the Trizol reagent kit (Invitrogen, Carlsbad, CA, USA) according to the manufacturer’s instructions. The RNA quality was assessed using an Agilent 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA, USA) and confirmed by RNase-free agarose gel electrophoresis. Following extraction, eukaryotic mRNA was selectively enriched with Oligo(dT) beads. The mRNA was enriched and then fragmented into short sequences using a fragmentation buffer. Subsequently, it was reverse-transcribed into cDNA with the NEBNext Ultra RNA Library Prep Kit for Illumina (NEB #7530, New England Biolabs, Ipswich, MA, USA). The resulting double-stranded cDNA fragments were purified, end-repaired, had an A base added, and were ligated to Illumina sequencing adapters. This ligation product was purified with AMPure XP Beads (1.0X) and amplified using PCR. Finally, the cDNA library was sequenced on the Illumina NovaSeq 6000 (Repugene Technology, Hangzhou, China).
Association analysis of Psi Score and drug sensitivity
The Genomics of Drug Sensitivity in Cancer (GDSC) provides transcription profiles for approximately 1000 cancer cell lines, drug response measurements, and targets/pathways of drugs [82]. The correlation between Psi Score and drug sensitivity was determined using Spearman correlation analysis and with a threshold of |Rs|> 0.3 and P < 0.05 for significance. ICI therapy responses in patients were predicted by using the Tumor Immune Dysfunction and Exclusion (TIDE) algorithm [83].
Calculation of TME cell invasion abundance
We employed the CIBERSORT algorithm to quantify the relative abundance of 22 types of immune cells in colorectal cancer. The parameters used were as follows: our gene expression matrix served as the input mixture matrix, gene signature reference for 22 immune cell types was obtained from Newman et al. [84], permutation test was conducted 100 times, and RNA-seq data was used without quantile normalization, while microarray data was normalized using quantile normalization.
Statistical analysis
We employed Wilcoxon tests for non-normally distributed data and Student’s t-tests for normally distributed data to compare the two groups. The Kruskal–Wallis test was utilized to analyze non-normally distributed data from more than two groups. For statistical analysis comparing non-categorical values among more than two groups, the Cochran-Mantel–Haenszel test was used. Spearman correlation analysis was conducted to assess the correlations between two variables. In survival analysis, the statistical significance of the Kaplan–Meier survival curve was determined using the log-rank test, with P < 0.05 considered statistically significant. All statistical analyses were performed using R (version 4.3.1).
Data availability
All data generated and analyzed during this study are included in this article and supplementary information files. The data that support the findings of this study have been deposited in the Gene Expression Omnibus under accession codes GSE271639 (for RNA-seq data). Source data are provided with this paper. All other data supporting the findings of this study are available from the corresponding authors upon reasonable request.
Abbreviations
- CRC:
-
Colorectal cancer
- CS:
-
Cellular senescence
- Ψ:
-
Pseudouridine
- RS:
-
Replicative senescence
- SIPS:
-
Stress-induced premature senescence
- PUS:
-
Pseudouridine synthases
- LUAD:
-
Lung adenocarcinoma
- CR:
-
Complete response
- PR:
-
Partial response
- SD:
-
Stable disease
- PD:
-
Progressive disease
- DEGs:
-
Differentially expressed genes
- GSVA:
-
Gene Set Variation Analysis
References
Di Micco R, Krizhanovsky V, Baker D. d’Adda di Fagagna F. Cellular senescence in ageing: from mechanisms to therapeutic opportunities. Nat Rev Mol Cell Biol. 2021;22(2):75–95.
Romesser PB, Lowe SW. The potent and paradoxical biology of cellular senescence in cancer. Ann Rev Cancer Biol. 2023;7:207–28.
Schmitt CA, Wang B, Demaria M. Senescence and cancer - role and therapeutic opportunities. Nat Rev Clin Oncol. 2022;19(10):619–36.
Lucas V, Cavadas C, Aveleira CA. Cellular senescence: from mechanisms to current biomarkers and senotherapies. Pharmacol Rev. 2023;75(4):675–713.
Wang R, Sun L, Xia S, Wu H, Ma Y, Zhan S, Zhang G, Zhang X, Shi T, Chen W. B7–H3 suppresses doxorubicin-induced senescence-like growth arrest in colorectal cancer through the AKT/TM4SF1/SIRT1 pathway. Cell Death Dis. 2021;12(5):453.
Montes M, Lubas M, Arendrup FS, Mentz B, Rohatgi N, Tumas S, Harder LM, Skanderup AJ, Andersen JS, Lund AH. The long non-coding RNA MIR31HG regulates the senescence associated secretory phenotype. Nat Commun. 2021;12(1):2459.
López-Otín C, Pietrocola F, Roiz-Valle D, Galluzzi L, Kroemer G. Meta-hallmarks of aging and cancer. Cell Metab. 2023;35(1):12–35.
Liu WW, Zheng SQ, Li T, Fei YF, Wang C, Zhang S, Wang F, Jiang GM, Wang H. RNA modifications in cellular metabolism: implications for metabolism-targeted therapy and immunotherapy. Signal Transduct Target Ther. 2024;9(1):70.
Nachtergaele S, He C. Chemical modifications in the life of an mRNA transcript. Annu Rev Genet. 2018;52:349–72.
Xiong Q, Zhang Y. Small RNA modifications: regulatory molecules and potential applications. J Hematol Oncol. 2023;16(1):64.
Charette M, Gray MW. Pseudouridine in RNA: what, where, how, and why. IUBMB Life. 2000;49(5):341–51.
Ge J, Yu YT. RNA pseudouridylation: new insights into an old modification. Trends Biochem Sci. 2013;38(4):210–8.
Karijolich J, Yi C, Yu YT. Transcriptome-wide dynamics of RNA pseudouridylation. Nat Rev Mol Cell Biol. 2015;16(10):581–5.
Hamma T, Ferré-D’Amaré AR. Pseudouridine synthases. Chem Biol. 2006;13(11):1125–35.
Wu G, Yu AT, Kantartzis A, Yu YT. Functions and mechanisms of spliceosomal small nuclear RNA pseudouridylation. Wiley Interdiscip Rev RNA. 2011;2(4):571–81.
Morais P, Adachi H, Yu YT. The Critical Contribution of Pseudouridine to mRNA COVID-19 Vaccines. Frontiers in cell and developmental biology. 2021;9: 789427.
Cerneckis J, Cui Q, He C, Yi C, Shi Y. Decoding pseudouridine: an emerging target for therapeutic development. Trends Pharmacol Sci. 2022;43(6):522–35.
Nance KD, Meier JL. Modifications in an emergency: the role of N1-Methylpseudouridine in COVID-19 vaccines. ACS Cent Sci. 2021;7(5):748–56.
Corbett KS, Edwards DK, Leist SR, Abiona OM, Boyoglu-Barnum S, Gillespie RA, Himansu S, Schäfer A, Ziwawo CT, DiPiazza AT, et al. SARS-CoV-2 mRNA vaccine design enabled by prototype pathogen preparedness. Nature. 2020;586(7830):567–71.
Nir R, Hoernes TP, Muramatsu H, Faserl K, Karikó K, Erlacher MD, Sas-Chen A, Schwartz S. A systematic dissection of determinants and consequences of snoRNA-guided pseudouridylation of human mRNA. Nucleic Acids Res. 2022;50(9):4900–16.
Ofengand J. Ribosomal RNA pseudouridines and pseudouridine synthases. FEBS Lett. 2002;514(1):17–25.
Han L, Kon Y, Phizicky EM. Functional importance of Ψ38 and Ψ39 in distinct tRNAs, amplified for tRNAGln(UUG) by unexpected temperature sensitivity of the s2U modification in yeast. RNA (New York, NY). 2015;21(2):188–201.
Papaccio F, García-Mico B, Gimeno-Valiente F, Cabeza-Segura M, Gambardella V, Gutiérrez-Bravo MF, Alfaro-Cervelló C, Martinez-Ciarpaglini C, Rentero-Garrido P, Zúñiga-Trejos S, et al. Proteotranscriptomic analysis of advanced colorectal cancer patient derived organoids for drug sensitivity prediction. J Exp Clin Cancer Res : CR. 2023;42(1):8.
Siegel RL, Wagle NS, Cercek A, Smith RA, Jemal A. Colorectal cancer statistics, 2023. CA: Cancer J Clin. 2023;73(3):233–54.
Fraser HC, Kuan V, Johnen R, Zwierzyna M, Hingorani AD, Beyer A, Partridge L. Biological mechanisms of aging predict age-related disease co-occurrence in patients. Aging Cell. 2022;21(4): e13524.
Zhu Y, Qiao L, Zhou Y, Ma N, Wang C, Zhou J. Long non-coding RNA FOXD2-AS1 contributes to colorectal cancer proliferation through its interaction with microRNA-185-5p. Cancer Sci. 2018;109(7):2235–42.
Penzo M, Guerrieri AN, Zacchini F, Treré D, Montanaro L. RNA pseudouridylation in physiology and medicine: for better and for worse. Genes. 2017;8(11):301. https://doiorg.publicaciones.saludcastillayleon.es/10.3390/genes8110301.
Barbieri I, Kouzarides T. Role of RNA modifications in cancer. Nat Rev Cancer. 2020;20(6):303–22.
Conde C, Mark M, Oliver FJ, Huber A, de Murcia G, Ménissier-de Murcia J. Loss of poly(ADP-ribose) polymerase-1 causes increased tumour latency in p53-deficient mice. EMBO J. 2001;20(13):3535–43.
Shi J, Wang L, Yin X, Wang L, Bo L, Liu K, Feng K, Lin S, Xu Y, Ning S, et al. Comprehensive characterization of clonality of driver genes revealing their clinical relevance in colorectal cancer. J Transl Med. 2022;20(1):362.
Hänzelmann S, Castelo R, Guinney J. GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinformatics. 2013;14: 7.
Zundler S, Neurath MF. Integrating immunologic signaling networks: the JAK/STAT pathway in colitis and colitis-associated cancer. Vaccines. 2016;4(1):5. https://doiorg.publicaciones.saludcastillayleon.es/10.3390/vaccines4010005.
Waldner MJ, Neurath MF. Master regulator of intestinal disease: IL-6 in chronic inflammation and cancer development. Semin Immunol. 2014;26(1):75–9.
Kan G, Wang Z, Sheng C, Yao C, Mao Y, Chen S. Inhibition of DKC1 induces telomere-related senescence and apoptosis in lung adenocarcinoma. J Transl Med. 2021;19(1):161.
Wang X, Ma L, Pei X, Wang H, Tang X, Pei JF, Ding YN, Qu S, Wei ZY, Wang HY, et al. Comprehensive assessment of cellular senescence in the tumor microenvironment. Brief Bioinform. 2022;23(3):bbac118.
Avelar RA, Ortega JG, Tacutu R, Tyler EJ, Bennett D, Binetti P, Budovsky A, Chatsirisupachai K, Johnson E, Murray A, et al. A multidimensional systems biology analysis of cellular senescence in aging and disease. Genome Biol. 2020;21(1):91.
Saul D, Kosinsky RL, Atkinson EJ, Doolittle ML, Zhang X, LeBrasseur NK, Pignolo RJ, Robbins PD, Niedernhofer LJ, Ikeno Y, et al. A new gene set identifies senescent cells and predicts senescence-associated pathways across tissues. Nat Commun. 2022;13(1):4827.
Deng EZ, Fleishman RH, Xie Z, Marino GB, Clarke DJB, Ma’ayan A. Computational screen to identify potential targets for immunotherapeutic identification and removal of senescence cells. Aging Cell. 2023;22(6): e13809.
Atlas A. a multi-omics database for aging biology. Nucleic Acids Res. 2021;49(D1):D825-d830.
Chatsirisupachai K, Palmer D, Ferreira S, de Magalhães JP. A human tissue-specific transcriptomic analysis reveals a complex relationship between aging, cancer, and cellular senescence. Aging Cell. 2019;18(6): e13041.
Tacutu R, Craig T, Budovsky A, Wuttke D, Lehmann G, Taranukha D, Costa J, Fraifeld VE, de Magalhães JP: Human ageing genomic resources: integrated databases and tools for the biology and genetics of ageing. Nucleic Acids Res. 2013;41(Database issue):D1027–1033.
Hernandez-Segura A, Nehme J, Demaria M. Hallmarks of cellular senescence. Trends Cell Biol. 2018;28(6):436–53.
Guinney J, Dienstmann R, Wang X, de Reyniès A, Schlicker A, Soneson C, Marisa L, Roepman P, Nyamundanda G, Angelino P, et al. The consensus molecular subtypes of colorectal cancer. Nat Med. 2015;21(11):1350–6.
Bhukdee D, Nuwongsri P, Israsena N, Sriswasdi S. Improved delineation of colorectal cancer molecular subtypes and functional profiles with a 62-gene panel. Mol Cancer Res : MCR. 2023;21(3):240–52.
Dienstmann R, Vermeulen L, Guinney J, Kopetz S, Tejpar S, Tabernero J. Consensus molecular subtypes and the evolution of precision medicine in colorectal cancer. Nat Rev Cancer. 2017;17(4):268.
Malla SB, Byrne RM, Lafarge MW, Corry SM, Fisher NC, Tsantoulis PK, Mills ML, Ridgway RA, Lannagan TRM, Najumudeen AK, et al. Pathway level subtyping identifies a slow-cycling biological phenotype associated with poor clinical outcomes in colorectal cancer. Nat Genet. 2024;56(3):458–72.
Zhou G, Myers R, Li Y, Chen Y, Shen X, Fenyk-Melody J, Wu M, Ventre J, Doebber T, Fujii N, et al. Role of AMP-activated protein kinase in mechanism of metformin action. J Clin Investig. 2001;108(8):1167–74.
Barzilai N, Crandall JP, Kritchevsky SB, Espeland MA. Metformin as a tool to target aging. Cell Metab. 2016;23(6):1060–5.
Thangthaeng N, Rutledge M, Wong JM, Vann PH, Forster MJ, Sumien N. Metformin impairs spatial memory and visual acuity in old male mice. Aging Dis. 2017;8(1):17–30.
Triggle CR, Mohammed I, Bshesh K, Marei I, Ye K, Ding H, MacDonald R, Hollenberg MD, Hill MA. Metformin: is it a drug for all reasons and diseases? Metabolism. 2022;133:155223.
Weiss EP, Fontana L. Caloric restriction: powerful protection for the aging heart and vasculature. Am J Physiol Heart Circ Physiol. 2011;301(4):H1205-1219.
Tatar M, Post S, Yu K. Nutrient control of drosophila longevity. Trends Endocrinol Metab. 2014;25(10):509–17.
Kulkarni AS, Gubbi S, Barzilai N. Benefits of metformin in attenuating the hallmarks of aging. Cell Metab. 2020;32(1):15–30.
Chresta CM, Davies BR, Hickson I, Harding T, Cosulich S, Critchlow SE, Vincent JP, Ellston R, Jones D, Sini P, et al. AZD8055 is a potent, selective, and orally bioavailable ATP-competitive mammalian target of rapamycin kinase inhibitor with in vitro and in vivo antitumor activity. Can Res. 2010;70(1):288–98.
Walters HE, Deneka-Hannemann S, Cox LS. Reversal of phenotypes of cellular senescence by pan-mTOR inhibition. Aging. 2016;8(2):231–44.
Gibney GT, Weiner LM, Atkins MB. Predictive biomarkers for checkpoint inhibitor-based immunotherapy. Lancet Oncol. 2016;17(12):e542–51.
Küçükköse E, Heesters BA, Villaudy J, Verheem A, Cercel M, van Hal S, Boj SF, Borel Rinkes IHM, Punt CJA, Roodhart JML, et al. Modeling resistance of colorectal peritoneal metastases to immune checkpoint blockade in humanized mice. J Immunother Cancer. 2022;10(12):e005345.
Chan TA, Yarchoan M, Jaffee E, Swanton C, Quezada SA, Stenzinger A, Peters S. Development of tumor mutation burden as an immunotherapy biomarker: utility for the oncology clinic. Ann Oncol : Off J Eur Soc Med Oncol. 2019;30(1):44–56.
Vandekerkhove G, Lavoie JM, Annala M, Murtha AJ, Sundahl N, Walz S, Sano T, Taavitsainen S, Ritch E, Fazli L, et al. Plasma ctDNA is a tumor tissue surrogate and enables clinical-genomic stratification of metastatic bladder cancer. Nat Commun. 2021;12(1):184.
Kan G, Wang Z, Sheng C, Chen G, Yao C, Mao Y, Chen S. Dual Inhibition of DKC1 and MEK1/2 synergistically restrains the growth of colorectal cancer cells. Adv Sci (Weinheim, Baden-Wurttemberg, Germany). 2021;8(10):2004344.
Fang Z, Shen HY, Xu Q, Zhou HL, Li L, Yang SY, Zhu Z, Tang JH. PUS1 is a novel biomarker for predicting poor outcomes and triple-negative status in breast cancer. Front Oncol. 2022;12: 1030571.
Lan C, Huang X, Liao X, Zhou X, Peng K, Wei Y, Han C, Peng T, Wang J, Zhu G. PUS1 may be a potential prognostic biomarker and therapeutic target for hepatocellular carcinoma. Pharmacogenomics Personalized Med. 2023;16:337–55.
Hu YX, Diao LT, Hou YR, et al. Pseudouridine synthase 1 promotes hepatocellular carcinoma through mRNA pseudouridylation to enhance the translation of oncogenic mRNAs. Hepatology. 2024;80(5):1058–73. https://doiorg.publicaciones.saludcastillayleon.es/10.1097/HEP.0000000000000702.
Cui Q, Yin K, Zhang X, Ye P, Chen X, Chao J, Meng H, Wei J, Roeth D, Li L, et al. Targeting PUS7 suppresses tRNA pseudouridylation and glioblastoma tumorigenesis. Nat Cancer. 2021;2(9):932–49.
Park EJ, Lee JH, Yu GY, He G, Ali SR, Holzer RG, Osterreicher CH, Takahashi H, Karin M. Dietary and genetic obesity promote liver inflammation and tumorigenesis by enhancing IL-6 and TNF expression. Cell. 2010;140(2):197–208.
Ihara S, Hirata Y, Koike K. TGF-β in inflammatory bowel disease: a key regulator of immune cells, epithelium, and the intestinal microbiota. J Gastroenterol. 2017;52(7):777–87.
Weichhart T, Hengstschläger M, Linke M. Regulation of innate immune cell function by mTOR. Nat Rev Immunol. 2015;15(10):599–614.
Kim YY, Jee HJ, Um JH, Kim YM, Bae SS, Yun J. Cooperation between p21 and Akt is required for p53-dependent cellular senescence. Aging Cell. 2017;16(5):1094–103.
Li P, Gan Y, Xu Y, Song L, Wang L, Ouyang B, Zhang C, Zhou Q. The inflammatory cytokine TNF-α promotes the premature senescence of rat nucleus pulposus cells via the PI3K/Akt signaling pathway. Sci Rep. 2017;7: 42938.
Kawamura H, Nakatsuka R, Matsuoka Y, Sumide K, Fujioka T, Asano H, Iida H, Sonoda Y. TGF-β signaling accelerates senescence of human bone-derived CD271 and SSEA-4 double-positive mesenchymal stromal cells. Stem cell reports. 2018;10(3):920–32.
Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD. The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics (Oxford, England). 2012;28(6):882–3.
Mariathasan S, Turley SJ, Nickles D, Castiglioni A, Yuen K, Wang Y, Kadel EE III, Koeppen H, Astarita JL, Cubas R, et al. TGFβ attenuates tumour response to PD-L1 blockade by contributing to exclusion of T cells. Nature. 2018;554(7693):544–8.
Li X, Ma S, Yi C. Pseudouridine: the fifth RNA nucleotide with renewed interests. Curr Opin Chem Biol. 2016;33:108–16.
Sun H, Li K, Liu C, Yi C. Regulation and functions of non-m(6)A mRNA modifications. Nat Rev Mol Cell Biol. 2023;24(10):714–31.
Wilkerson MD, Hayes DN. ConsensusClusterPlus: a class discovery tool with confidence assessments and item tracking. Bioinformatics (Oxford, England). 2010;26(12):1572–3.
Steinley D. K-means clustering: a half-century synthesis. Br J Math Stat Psychol. 2006;59(Pt 1):1–34.
Yu G, Wang LG, Han Y, He QY. clusterProfiler: an R package for comparing biological themes among gene clusters. OMICS. 2012;16(5):284–7.
Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, Smyth GK. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43(7): e47.
Law CW, Chen Y, Shi W, Smyth GK. voom: Precision weights unlock linear model analysis tools for RNA-seq read counts. Genome Biol. 2014;15(2): R29.
Tacutu R, Thornton D, Johnson E, Budovsky A, Barardo D, Craig T, Diana E, Lehmann G, Toren D, Wang J, et al. Human ageing genomic resources: new and updated databases. Nucleic Acids Res. 2018;46(D1):D1083-d1090.
Rebhan M, Chalifa-Caspi V, Prilusky J, Lancet D. GeneCards: a novel functional genomics compendium with automated data mining and query reformulation support. Bioinformatics (Oxford, England). 1998;14(8):656–64.
Yang W, Soares J, Greninger P, Edelman EJ, Lightfoot H, Forbes S, Bindal N, Beare D, Smith JA, Thompson IR, et al. Genomics of Drug Sensitivity in Cancer (GDSC): a resource for therapeutic biomarker discovery in cancer cells. Nucleic Acids Res. 2013;41(Database issue):D955-961.
Jiang P, Gu S, Pan D, Fu J, Sahu A, Hu X, Li Z, Traugh N, Bu X, Li B, et al. Signatures of T cell dysfunction and exclusion predict cancer immunotherapy response. Nat Med. 2018;24(10):1550–8.
Becht E, Giraldo NA, Lacroix L, Buttard B, Elarouci N, Petitprez F, Selves J, Laurent-Puig P, Sautès-Fridman C, Fridman WH, et al. Estimating the population abundance of tissue-infiltrating immune and stromal cell populations using gene expression. Genome Biol. 2016;17(1):218.
Acknowledgements
None.
Funding
Support for these studies was provided by the National Natural Science Foundation of China (No. 82173067, YL; No. 81972245, YL; No. 82272965, HY; No. 81902877, HY; No. 31900505, YL), the Natural Science Foundation of Guangdong Province (No. 2022A1515012656, HY; No. 2021A1515010639, XL; No. 2021A1515010134, MH; No. 2022A1515012457, JC), Science and Technology Program of Guangzhou (No. 202201011004, HY), the "Five Five" Talent Team Construction Project of the Sixth Affiliated Hospital of Sun Yat-sen University (No. P20150227202010244, JW; No. P20150227202010251, YL), the Excellent Talent Training Project of the Sixth Affiliated Hospital of Sun Yat-sen University (No. R2021217202512965, YL), the Scientific Research Project of the Sixth Affiliated Hospital Of Sun Yat-sen University (No. 2022JBGS07, YL), the Fundamental Research Funds for the Central Universities, Sun Yat-sen University (No. 23ykbj007, HY), the Program of Introducing Talents of Discipline to Universities, and National Key Clinical Discipline (2012).
Author information
Authors and Affiliations
Contributions
Conceptualization: H.C.Y., J.X.L., and Y.G.G.; Investigation, Writing – original draft: Y.G.G., Z.Y., and J.R.W.; Methodology, Writing – review and editing: J.R.W., M.Z.H., T.Y.L., Y.H.W., and X.L.; Data acquisition: Y.G.G., Z.Y., J.R.W., and M.Z.H.; Data analysis and interpretation: Y.G.G., K.X.L., J.Y.H., H.T.L., Z.X.W., Z.M.L., Z.H.C., and J.C.; Study administration: H.C.Y., J.X.L., Y.X.L., and M.J.H.; All authors read and approved the final manuscript.
Corresponding authors
Ethics declarations
Ethics approval and consent to participate
The study protocol was reviewed and approved by the Institutional Review Board of the Sixth Affiliated Hospital of Sun Yat-sen University (no. 2024ZSLYEC-318).
Consent for publication
All authors read and approved the final manuscript.
Competing interests
The authors declare no competing interests.
Additional information
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Supplementary Information
12915_2025_2170_MOESM1_ESM.xlsx
Additional file 1. Figure S1 and S2. GO and KEGG enrichment analysis of the 43 DEGs. The x-axis represents gene counts within each GO term, with column brightness indicating the statistical significance of enrichment. Figure S3. Differences in the Psi Score between RNA Ψ modification patterns and gene clusters in five GEO-CRC cohorts. The Wilcoxon test was used to determine the statistical significance of the difference, and p < 0.05 was considered statistically significant. Figure S4. Venn diagram showing the overlap between high/low Psi Score groups and Cluster 1/2 in CRC, as well as geneCluster A/B. Figure S5. Senescence level in Ψ modification writers knockdown datasets. Figure S6. Heatmap of expression levels for intersecting genes between scoring systems. Figure S7. Correlation coefficients of the Psi Score and OE/UE score were calculated by using Pearson’s correlation coefficient formula. The red scatter points represent CRC samples and the blue line indicates the correlation. Figure S8. Comparison of GSVA scores for characteristics associated with GOBP pathways of senescence between Psi Score-high and -low groups. Figure S9. The difference of Age in Psi Score-high and -low groups in integrated cohort. Figure S10. Correlation coefficients of the Psi Score and two RNA Ψ modification pattern enrichment pathway score were calculated by using Pearson’s correlation coefficient formula. The red scatter points represent CRC samples and the blue line indicates the correlation. Figure S11. The differences in TME infiltration between Psi Score-high and -low groups in integrated cohort.
12915_2025_2170_MOESM2_ESM.zip
Additional file 2. Table S1. The annotation of 13 RNA Ψ modification writers based on the published data. Table S2. Clinical information of CRC cohorts from GEO. Table S3. Group information of samples in integrated GEO-CRC cohorts. Table S4. Enrichment score of Hallmark pathways in integrated GEO-CRC cohorts. Table S5. Sample group information in the 2024 TCGA-COAD cohort. Table S6. Summary of Clinical Features Associated with Psi Score and CS Score in integrated GEO-CRC cohorts and 2024 TCGA-COAD cohorts. Table S7. Group information of samples in integrated GEO-CRC cohorts. Table S8. The association of Psi Score and drugs sensitivity in GDSC database. Table S9. Review of literature on drugs sensitive to Psi Score-high tumors. Table S10. Summary of Clinical Features Associated with Psi Score in IMvigor210 cohorts.
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.
About this article
Cite this article
Gan, Y., Yuan, Z., Weng, J. et al. Transcriptomic profile of RNA pseudouridine modification as a biomarker for cellular senescence associated with survival outcomes in colorectal cancer. BMC Biol 23, 61 (2025). https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s12915-025-02170-6
Received:
Accepted:
Published:
DOI: https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s12915-025-02170-6