Abstract This study investigated whether m5C-related Long non-coding RNAs (lncRNAs) can predict clear cell renal cell carcinoma (ccRCC) patient prognosis. Co-expression and Cox regression analyses identified 9 prognostic lncRNAs, which were closely associated with tumor immune characteristics and immune escape. The model also predicted the sensitivity of drugs, including Entinostat, SB216763, and Sapitinib. In vitro experiments showed that GNG12-AS1 inhibited ccRCC cell proliferation and migration by reducing the activity of the ERK/GSK-3β/β-catenin pathway. Overall, these findings suggest that the 9 m5C-related lncRNAs can accurately predict ccRCC patient prognosis, providing potential applications for clinical and immunotherapy approaches. GNG12-AS1 emerges as a promising prognostic biomarker for predicting survival outcomes in ccRCC, potentially influencing cell migration through the activation of the ERK/GSK-3β/β-catenin signaling pathway. Keywords: m5C, lncRNA, ccRCC, Prognosis, GNG12-AS1 Introduction In men, kidney cancer accounts for 5% of all malignancies and is the sixth most common cancer. In women, it accounts for 3% of malignancies and is the tenth most common cancer [[34]1]. The most common histological type of RCC is clear cell renal cell carcinoma, which accounts for about 70–90% of all types of RCC [[35]2]. Compared with other subtypes of RCC, ccRCC has the characteristics of highly aggressive and high recurrence rate [[36]3]. At the time of initial diagnosis, approximately 30% of RCC have already metastasized [[37]4]. Renal cell carcinoma is poorly sensitive to chemo- and radiotherapy, thus has a high risk of recurrence after nephrectomy, with a 5-year survival rate of only 60–70% [[38]5]. Therefore, understanding the molecular mechanism of ccRCC and finding sensitive and reliable molecular markers that can predict the occurrence and development of ccRCC are of great significance to improve the early diagnosis and treatment of ccRCC. LncRNA is a kind of RNA with a length of more than 200 nucleotides, which is unable or has weak ability to be translated into proteins [[39]6]. Like most non-coding RNAs, lncRNAs are transcribed by RNA polymerase II and have a 5′ methyl-cytosine cap and 3′ poly(A) tail. Previously, lncRNAs were thought to had no biological function [[40]7]. However, with the deepening of understanding, it has been found that lncRNAs are widely distributed in eukaryotes and participate in various important biological functions, including the regulation of autophagy [[41]8], specific transcriptional and post-transcriptional controls [[42]9], transduction of intermolecular signals [[43]10], and regulation of protein stability after translation [[44]11]. In the research of cancer, lncRNAs are differentially expressed in prostate cancer, breast cancer, colorectal cancer, pancreatic cancer, non-small cell lung cancer, hepatocellular carcinoma and so on, and play an important role in the proliferation, migration and invasion of tumor cells [[45]12–[46]17]. LncRNA CACClnc could be used as a predictor of the efficacy of chemotherapy in patients with colorectal cancer [[47]18]. LncRNA SNHG7 could promote disease progression and drug resistance by inducing autophagy activity in non-small cell lung cancer [[48]19]. However, the role of lncRNAs in ccRCC needs to be further studied. Some studies had shown that RNA modification plays an important role in the development of tumors [[49]20]. N6-methylcytosine (m6A) is the main type of modification in eukaryotic RNA, while 5-methylcytosine (m5C) is another common RNA modification [[50]21]. More and more studies had shown that m5C modification can affect the progression of a variety of cancers. For example, the latest study by Chen Xue et al. pointed out that ALYER-mediated RNA m5C modification plays an important role in the progression of hepatocellular carcinoma [[51]22]. Similarly, Min Yang et al. demonstrated that NSUN2 promoted osteosarcoma progression by inducing m5C modification to stabilize FABP5 [[52]23]. As an epigenetic modification of LncRNA, m5C regulates the function and interaction of lncRNA by changing its own structure and properties [[53]24]. M5C has similar effects with m6A, and the epigenetic features of its related lncRNAs are also considered as potential biomarkers for prognosis [[54]24]. Ruoxin Xu et al. identified five m5C-related lncRNAs that may serve as prognostic classifiers for patients with lung squamous cell carcinoma [[55]25]. However, the relationship between lncRNAs and m5C modification in ccRCC has not been reported. The abnormal expression of lncRNAs and modification of m5C may lead to the occurrence and development of cancers. However, the interactions between lncRNAs and m5C in ccRCC remain unclear. In this study, we aimed to use bioinformatics methods to construct a risk prediction model of m5C-related lncRNAs to predict the prognosis of ccRCC patients. Materials and methods Data collection The transcriptome profiling, clinical data and simple nucleotide variation of ccRCC cases, and normal cases were all downloaded from the TCGA database. A total of 31 m5C-related genes (MRGs) were collected from the summary of the previous reports [[56]26, [57]27]. Identification of m5C-related lncRNAs in ccRCC LncRNAs were obtained in RNA-Seq data through gene annotation of TCGA. These lncRNAs were identified using the "limma" package in R software to determine those significantly correlated with m5C. Correlation was considered significant with a Pearson correlation coefficient > 0.75 and p < 0.05. The co-expression network was visualized using a Sankey diagram. Construction of m5C-related lncRNAs prognostic signature for ccRCC First, the "survival" package in R was used to identify m5C-related lncRNAs associated with overall survival (OS) prognosis through univariate Cox regression analysis (p < 0.05). The hazard ratio (HR) for each m5C-related lncRNA was calculated. Next, the "glmnet" package in R was utilized for LASSO regression analysis to determine the optimal number of m5C-related lncRNAs to include in the model. The samples were then randomly divided into a training group and a testing group in a 1:1 ratio to reduce bias, improve the model’s generalization ability, and enhance its prediction performance for new data. The training group was used to develop the prognostic model, while the testing group was used to validate the model's accuracy. Multifactorial regression analysis was conducted to identify the best combination of m5C-related lncRNAs for the model. The optimal parameters from this analysis were used to construct the final model, and the risk score for each sample was calculated as follows: risk score = (Expression LncRNA1 × Coefficient LncRNA1) + (Expression LncRNA2 × Coefficient LncRNA2) + … + (Expression LncRNAn × Coefficient LncRNAn). Analysis of risk signature Patients in the training and testing groups were divided into high-risk and low-risk groups based on median risk score. The "survival" package in R software was used to identify the differences in OS and progression-free survival (PFS) between patients in two groups. And the chi-square test was used to assess the correlation between the model and clinical features. Based on survival, caret, glmnet, rms, survminer and time ROC, we generated 1-year, 3-year and 5-year ROC curves of each group, and calculated the area under the curves (AUC). The consistency index (C-index) was applied to evaluate the accuracy of the model. Construction of nomogram and validation of clinical subgroups Based on the regplot and rms packages, nomograms were constructed using the sample's age, gender, grade, stage and risk score. The calibration curves of 1 year, 3 years and 5 years were drawn to assess the effect of prognostic model.The c-index curve was used to judge the accuracy of the model.We classified patients into stages I-II and III-IV to observe whether the constructed model is applicable to both early and late stage patients. Principal component analysis (PCA) and functional enrichment analysis The limma package and scatterplot3d package were used to perform PCA analysis of all genes, m5C genes, m5C-related lncRNAs and model lncRNAs, and the distribution of patients with different risk scores was observed. We used limma package to detect the differentially expressed genes (DEGs) between high and low risk groups with the threshold of false discovery rate (FDR) < 0.05 and |log2FC|≥ 1. Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis were performed on these DEGs. Outcoms with P < 0.05 and FDR < 0.05 were considered statistically significant. Immune-related functional analysis and tumor mutation burden (TMB) analysis The immune-related functions of the DEGs between high and low risk groups were analyzed based on limma, GSVA and GSEABase, and the pheatmap package was used to visualize the heatmap. We downloaded the data of TMB of ccRCC patients from the TCGA, and then used the maftools package to visualize the differences of tumor mutation data between high and low risk groups. The TMB was visualized using the ggpubr package. Immunotherapy analysis and pharmaceutical screening The tumor immune dysfunction and exclusion (TIDE) data of ccRCC was downloaded from the TIDE website ([58]http://tide.dfci.harvard.edu/). TIDE is a computational method developed to predict immune checkpoint blocking responses and simulate tumor immune escape [[59]28]. It can be used to evaluate the immune escape potential and the efficacy of immunotherapy in high or low risk groups. OncoPredict, parallel, limma, ggplot2 and ggpubr packages were used to analyze the differences of drug sensitivity between high and low risk groups. The OncoPredict package generated predictive models of drug response, which could be used to screen therapeutic drugs and observe drug sensitivity [[60]29]. Cell transfection The overexpression plasmid of GNG12-AS1 was constructed by Genechem (Shanghai, China), after which the plasmid was transfected using lipofectamine 3000 (Invitrogen, USA), all procedures were performed in strict accordance with the operating instructions. Transfection efficacy was confirmed by qRT-PCR analysis. Quantitative reverse-transcription polymerase chain reaction RNA was isolated from cell lines using Trizol reagent (Invitrogen, USA) and reverse transcribed into complementary DNA for amplification using PrimeScript™ RT reagent Kit with gDNA Eraser (TaKaRa, Japan). The primer design for GNG12-AS1 was as follows: forward: 5 '-GGAACCTGCGGATACAGA-3', reverse: 5 '-ATGACAAACAGGGACCACA-3'. GAPDH was used as an internal control (forward: 5 '-AGAAGGCTGGGGTCATTTG-3', reverse: 5 '-GCAGGAGGCATTGCTGATGAT-3').We used TB Green Premix Ex Taq II (Tli RNaseH Plus)(2 ×)(TaKaRa, Japan) as the marker for RT-qPCR and performed the experiments using the Roche LightCycler480 II PCR apparatus(Roche, Switzerland). After completing the experiment, we examined the amplification and melting curves of each sample to exclude any abnormal data. We recorded the cycle threshold (Ct) values at which each sample reached the fluorescence detection threshold and calculated the relative gene expression using the 2^(-ΔΔCt) method. Finally, we conducted statistical analysis and data visualization with GraphPad Prism software (version 8.0). Wound-healing assay Cells were seeded in 6-well plates at a density of 70 percent and incubated for 24 h. Wounds were generated on the cell surface through a sterile 200 μL pipettor suction tip. The wound healing was observed at 12, 24 and 36 h. Transwell assay A total of 3 × 10^4 over-expressing GNG12-AS1 cells were seeded into the upper chamber of Transwell plates containing serum-free medium. Whereas the lower chamber contained 20%FBS as a chemoattractant for the culture. Cells on the lower surface of the filter were fixed with tissue fixator after 24 h, stained with 0.1% crystal violet, and counted under a microscope. All experiments were performed three times. Cell viability assay A total of 5 × 10^3 cells were inoculated in 96-well plates and cultured for 24, 48, 72, 96, and 120 h, respectively. The cells were then cultured with 10μlCCK-8 solution (Beyotime Biotechnology, China) at 37 ℃ for 2 h. Absorbance was measured by spectrophotometer at a wavelength of 450 nm. All results are expressed as the mean of three independent experiments ± SD. Western blot analysis The cells (786-O, OS-RC-2) were collected for Western blot analysis 48 h after successful plasmid transfection. The cells were lysed with RIPA buffer (P0013B, Beyotime Biotechnology) supplemented with a protease inhibitor (P1005, Beyotime Biotechnology) for 30 min on ice. The concentration of the extracted protein was determined using a colorimetry method. The entire procedure was conducted at low temperature. The collected protein samples were loaded onto a 10% SDS-PAGE gel, with 20 μL of sample in each well, and then separated by electrophoresis. Subsequently, the proteins were transferred from the gel to a nitrocellulose membrane. The membrane was blocked with milk for an additional two hours to prevent non-specific binding.For primary antibody incubation, the following antibodies were used: β-catenin (#8480), E-cadherin (#14472), N-cadherin (#13116), AKT (#4691), Phosphorylation-AKT (#4060), Erk1/2 (#4695), Phosphorylation-Erk1/2 (#9101), GAPDH (#97166) from Cell Signaling Technology, China, GSK-3β (Cat No.22104-1-AP), Phosphorylation-GSK-3β (Cat No.67558-1-lg) from Proteintech. The primary antibodies were incubated overnight at 4 °C while shaking. After primary antibody incubation, the membrane was incubated with secondary antibodies against mouse or rabbit IgG for subsequent detection. Specifically, the anti-mouse IgG secondary antibody (bs-0296G-HRP, Bioss) and anti-rabbit IgG secondary antibody (bs-0295G-HRP, Bioss) were used. Color development was performed using Enhanced Chemiluminescence (ECL) reagent (EK-5008-500 mL, ECOTOP) after completion of the incubation. However, it should be noted that due to the small sample size for each sample, some of our blots were lysed prior to hybridization with the antibody. Nevertheless, within the same cell line, we ensured consistency between the overexpression group and the control group in each experiment. Finally, we developed the protein bands through a G-BOX-mini9 (SYNGENE, USA) multifunctional chemiluminescence imager. We processed the results with GeneSys software, and when exporting the images, we set the DPI parameter of the images to 300, which was the maximum setting available in GeneSys software, and all the images had the same pixel size. Results Identification of prognostic m5C-Related lncRNAs in ccRCC Using |R|> 0.5 and P < 0.001 as filtering criteria, lncRNAs with co-expression relationship were identified from 16,876 lncRNAs and 31 MRGs, and a total of 422 co-expression lncRNAs were obtained. These results were visualized in Sankey diagram (Fig. [61]1A). Using P < 0.05 as the filtering criterion, the training group was firstly subjected to univariate COX analysis, and 36 m5C-related lncRNAs related to prognosis were screened out. Then lasso regression analysis was performed on the results of univariate analysis, and cross-validation was used to screen out the point with the smallest error (Fig. [62]1B–D). Then the COX model was constructed to obtain the optimal model for the 9 m5C-related lncRNAs screened by lasso regression. Correlation heat map was used to show the relationship between MRGs and m5C-related lncRNAs constructed by the model, in which red represented positive correlation and blue represented negative correlation (Fig. [63]1E). Risk score = (−0.732105208538704*AP004609.3) + (1.24333216856082*AC004837.2)  + (1.25215186631223*[64]AL133255.1) + (1.09407448306765*[65]AC044781.1 ) + (−1.72051844326631*GNG12-AS1) + (−2.27930443465364*[66]AP002755.1)  + (−2.47688775750224*[67]AC083805.2) + (−4.36659077891361*ERI3-IT1) + ( 4.50355534403254*[68]AC078795.1). Through the obtained model formula, we could calculate the risk score of each sample in the training and testing groups. According to the median value of risk score, we divided the samples into high and low risk groups. Fig. 1. [69]Fig. 1 [70]Open in a new tab Identification of prognostic m5C-Related lncRNAs in ccRCC. A Sankey diagram showed the lncRNAs with co-expression relationship with 31 MRGs. B Univariate COX regression analysis identified 36 m5C-related lncRNAs associated with prognosis. C Cross-validation curve for lasso regression analysis. D The coefficient of m5C-related lncRNAs in LASSO regression analysis was calculated. E Multivariate Cox regression determined the correlation between 9 m5C-related lncRNAs and 31 MRGs Survival analysis of the signature Through the risk curve, we could observe the differences in the risk score and survival time of patients in high and low risk group (Fig. [71]2A, B). The heat map showed the expression level of 9 m5C-related lncRNAs in high and low risk groups. Among them, AC004837.2, [72]AL133255.1, [73]AC044781.1 and [74]AC078795.1 were high-risk lncRNAs, while AP004609.3, GNG12-AS1, [75]AP002755.1, [76]AC083805.2 and ERI3-IT1 were low-risk lncRNAs (Fig. [77]2C). In order to compare the survival of patients between high and low risk groups, survival analyses were performed on the training group, testing group and all groups. We found that the OS and PFS of patients in high risk group were significantly shorter than that in low risk group (Fig. [78]2D, E). Fig. 2. [79]Fig. 2 [80]Open in a new tab Comparisons of OS and PFS between high risk and low risk groups. A Distribution of risk scores among patients. B Distribution of survival time and survival status. C Heat map of the expression of 9 m5C-related lncRNAs. Comparisons of OS (D) and PFS (E) in the training group, testing group and all groups Independent prognostic analysis Univariate and multivariate independent prognostic analyses were performed to evaluate whether the prognostic model could be used as an independent prognostic factor characteristics. The results showed that our model was relevant to the prognosis of patients in both univariate and multivariate analyses. In univariate COX regression analysis, HR (hazard ratio) (95% CI) of risk score was 1.054 (1.041–1.067) with P < 0.05 (Fig. [81]3A). In multivariate COX regression analysis, the HR (95% CI) of risk score was 1.034 (1.020–1.048) with P < 0.05 (Fig. [82]3B). Then, we constructed ROC curve to judge the accuracy of our model in predicting patients’ survival. And the results showed that the AUC value of risk score (0.769) was significantly higher than that of age (0.652), gender (0.505) and grade (721) (Fig. [83]3D). In addition, the AUCs of one-year, three-year and five-year risk scores were 0.769, 0.718 and 0.724, respectively, indicating that the model had reliable accuracy in predicting the survival of ccRCC patients (Fig. [84]3C). Fig. 3. [85]Fig. 3 [86]Open in a new tab The reliability and accuracy of the prognostic model were evaluated. A Univariate independent prognostic analysis. B Multivariate independent prognostic analysis. C Area under the ROC curve at 1 year, 3 years, and 5 years. D The area under the ROC curve of the constructed model compared with that of other clinical characteristics Construction of nomogram and PCA analysis We used age, sex, grade, stage, T-stage and risk score to construct a nomogram for predicting 1 -, 3 -, and 5-year survival (Fig. [87]4A). For example, if the patient has a score of 133, the probability of survival over one year is 0.896, over three years is 0.725, and over five years is 0.578. The calibration curve showed that the nomogram had good reliability in predicting the survival of patients at 1, 3, and 5 years (Fig. [88]4B). In addition, we constructed the C-index curve and found that the C-index value of the prognostic model was significant (Fig. [89]4C), indicating that the accuracy of the prognostic model was reliable. Patients were divided into subgroups according to clinical characteristics (gender, grade, stage, and T-stage), and then survival analyses were conducted on these subgroups based on risk scores. It was found that the model was applicable to patients in all subgroups (Fig. [90]4D–G). Finally, we performed PCA on the model lncRNAs, m5C-related lncRNAs, MRGs, and all genes (Fig. [91]5). We observed that the lncRNAs in our model could well distinguish patients with high and low risk groups, indicating the prognostic model was reliable. Fig. 4. [92]Fig. 4 [93]Open in a new tab Nomogram and validation of clinical grouping model. A Prognostic nomogram of OS in ccRCC patients. B Calibration curves at 1, 3, and 5 years. C Risk curves for clinical characteristics and prognostic models. D Validation of the model in Stage I-II and Stage III–IV patients. E Validation of the model in T1-T2 and T3-T4 patients. F Validation of the model in male and female patients. G Validation of the model in G1-G2 and G3-G4 patients Fig. 5. [94]Fig. 5 [95]Open in a new tab PCA analyses. A Spatial distribution of model lncRNAs. B Spatial distribution of m5C-related lncRNAs. C Spatial distribution of MRGs. D Spatial distribution of all genes Enrichment analysis and immune-related function analysis GO enrichment analysis of DEGs between high and low risk groups showed that the biological process (BP), cellular component (CC) and molecular function (MF) were mostly enriched (Fig. [96]6A). The main components enriched in BP were phagocytosis and recognition, complement activation, and humoral immune response mediated by circulating immunoglobulin; the main enriched components in CC were immunoglobulin complex, circulation of immunoglobulin complex and external side of plasma membrane; and the main components enriched in MF were antigen binding, immunoglobulin receptor binding and endopeptidase inhibitor activity (Fig. [97]6A). KEGG analysis of these DEGs also suggested that these m5C-related lncRNAs might be associated with complement and coagulation cascades, viral protein interaction with cytokine receptor, and NF − kappa B signaling pathway (Fig. [98]6B). These results suggested that these lncRNAs might be involved in tumor development. Further analysis of immune-related functions showed that, except for major histocompatibility complex class I molecules (MHC_class_I), other types of immune responses were significantly higher in the high risk group than low risk group (Fig. [99]6C). This might explain why the high risk patients had worse survival compared with low risk patients. Fig. 6. [100]Fig. 6 [101]Open in a new tab Enrichment analyses and immune-related function analysis. A GO enrichment analysis of DEGs. B KEGG enrichment analysis of DEGs. C Immune-related functional analysis of DEGs TMB analysis and TIDE analysis By comparing the mutation frequency between high and low risk groups, it could be observed that the top ten genes with high mutation frequency in high risk group were VHL, PBRM1, TTN, SETD2, BAP1, MUC16, MTOR, DNAH9, KDM5C, DST, HMCN1, LRP2, CSMD3, FBN2, KMT2C. TTN, SETD2, BAP1, MTOR, HMCN1 and CSMD3 had a higher mutation frequency in the high-risk group (Fig. [102]7A), while VHL, PBRM1, DNAH9, KDM5C, DST and KMT2C had high mutation frequency in low risk group (Fig. [103]7B). In general, there was no significant difference in TMB between the two groups (Fig. [104]7C). In addition, we performed survival analyses on patients with high and low TMB, and observed that the OS of patients with high TMB was significantly longer than patients with low TMB. And this result was also true in the case of a more detailed division of the groups (Fig. [105]7D, E). We then scored TIDE for each sample expressing the input file through the online website ([106]http://tide.dfci.harvard.edu/), and observed that patients in high risk group had a greater potential for immune escape than low risk group, indicating that they had a worse response to immunotherapy (Fig. [107]7F). Fig. 7. [108]Fig. 7 [109]Open in a new tab TMB analyses and TIDE analysis. A Waterfall plot of the high risk group. B Waterfall plot of the low risk group. C Analysis of differences between high- and low- TMB groups. D Survival analysis between the high- and low-TMB groups. E Survival analysis between the four groups(H-TMB + high risk, H-TMB + low risk, L-TMB + high risk, L-TMB + low risk). F TIDE algorithm analysis for high and low risk groups Drug sensitivity analysis According to the expression file and drug susceptibility file of the database, we scored each sample of ccRCC for drug susceptibility. A lower score indicated more susceptible to this drug. Then we could obtain sensitivity scores of each sample for 198 drugs. Using P < 0.01 as a filter, we found that the IC50 (half maximal inhibitory concentration) values of Entinostat and SB216763 were higher in high risk group compared with low risk group (Fig. [110]8A, B), indicating patients in low risk group were more sensitive to these two drugs. In addition, the higher IC50 value of Sapitinib in low risk group suggested that patients in high risk group were more sensitive to Sapitinib (Fig. [111]8C). Fig. 8. [112]Fig. 8 [113]Open in a new tab Drug sensitivity analyses. A Drug sensitivity of Entinostat. B Drug sensitivity of SB216763. C Drug sensitivity of Sapitinib Role of lncRNA GNG12-AS1 in KIRC cells Since 9 related lncRNAs were identified, the expression profiles of these 9 lncRNAs in normal and tumor tissues were analyzed. Six lncRNAs (AC004837.2, [114]AL133255.1, [115]AC044781.1, [116]AC078795.1, GNG12-AS1, ERI3-IT1) were found to have statistically significant expression differences, and only GNG12-AS1 expression was down-regulated in tumor tissues (Fig. [117]9A). Therefore, we hypothesized that GNG12-AS1 may play an important role in the occurrence and development of KIRC. We then performed survival analysis of GNG12-AS1 expression and patient survival using the online database GEPIA2, and found that among KIRC patients, those with high GNG12-AS1 expression had better overall survival (Fig. [118]9B). GNG12-AS1 may act as a tumor suppressor gene and play an important role in the occurrence and development of KIRC. We performed a series of experiments to further explore how GNG12-AS1 regulates ccRCC progression in vitro. The plasmid overexpressing GNG12-AS1 was transfected into 786-O, 769-P, OS-RC-2 cells. The qPCR results showed good transfection efficiency (Fig. [119]9C). The results of wound healing assay showed that the scratch healing rates of 786-O, 769-P and OS-RC-2 cell lines overexpressing GNG12-AS1 were significantly lower than those of the control group. Transwell assays also showed that higher GNG12-AS1 expression reduced the migratory capacity of 786-O, 769-P, and OS-RC-2 cell lines. In the cell viability assay, the OD value at 450 nm of the GNG12-AS1 overexpressing cell line was lower than that of the control group at all time points (Fig. [120]10A–C). Fig. 9. [121]Fig. 9 [122]Open in a new tab Overexpression of GNG12-AS1. A Differential analysis of model lncRNAs expression. B Kaplan–Meier analysis of overall survival in patients with low and high GNG12-AS1 levels. C Overexpression of GNG12-AS1 in 786-O, 769-P, and OS-RC-2 cell line Fig. 10. [123]Fig. 10 [124]Open in a new tab Cell function Experiments. A Wound-Healing Assay. B Transwell Assay. C Cell Viability Assay GNG12-AS1 may activate ERK/GSK-3β/β-catenin signaling pathway The aforementioned results indicate that GNG12-AS1 has a significant impact on ccRCC cell migration. After reviewing relevant literature, it was found that GNG12-AS1 can activate the AKT/GSK-3β/β-catenin signaling pathway in glioma cells, thereby promoting the migration and metastasis of these cells [[125]30]. Therefore, we investigated related biomarkers such as N-cadherin, E-cadherin, β-catenin, GSK-3β, AKT, etc. The results demonstrated a significant decrease in β-catenin expression when GNG12-AS1 was overexpressed. However, the total protein expression of N-cadherin, E-cadherin, AKT, and GSK-3β did not show significant changes. Furthermore, the expression of Phosphorylation-AKT increased while Phosphorylation-GSK-3β decreased (Fig. [126]11). These findings indicate that GNG12-AS1 regulates cell migration through the GSK-3β/β-catenin signaling pathway in ccRCC. However, there may be another upstream signaling pathway that influences GSK-3β. The ERK/GSK-3β/Snail1 signaling pathway plays a crucial role in the formation and development of renal fibrosis [[127]31]. Based on this, we hypothesized that GNG12-AS1 may affect ccRCC cell migration through the ERK/GSK-3β/β-catenin signaling pathway. To test this hypothesis, we conducted Western blot experiments, which revealed that the total protein expression of ERK1/2 did not significantly change upon GNG12-AS1 overexpression (Fig. [128]11). However, the expression of Phosphorylation-ERK1/2 was significantly decreased. In conclusion, our study suggests that GNG12-AS1 may have a significant impact on ccRCC cell migration by activating the ERK/GSK-3β/β-catenin signaling pathway. Fig. 11. Fig. 11 [129]Open in a new tab Upregulation of GNG12-AS1 inhibits ERK/GSK-3β/β-catenin signaling pathway in renal clear cell carcinoma Discussion Although the basic treatment of ccRCC is surgery, due to the lack of sensitivity of renal cancer to radiotherapy and traditional chemotherapy, the five-year survival rate of patients will decrease rapidly once cancer cells occur invasion and metastasis [[130]32, [131]33]. Therefore, it is important to identify the prognosis of ccRCC with molecular heterogeneity. M5C is involved in various important biological processes in vivo [[132]34]. It is well known that lncRNAs serve as therapeutic targets and innovative biomarkers for a variety of cancers due to their tissue-specific expression characteristics and genome-wide expression patterns [[133]35]. Through TCGA data analysis, our research identified nine m5C-related lncRNAs (AC004837.2, [134]AL133255.1, [135]AC044781.1, [136]AC078795.1, AP004609.3, GNG12-AS1, [137]AP002755.1, [138]AC083805.2,ERI3-IT). We established relevant prognostic models based on their risk scores, and their predictive efficiency was well verified in both training and testing group. By reviewing the relevant literature, we found that only similar [139]AP004609.1 and GNG12-AS1 have been reported in the relevant literature. [140]AP004609.1 may serve as a potential biomarker for the prognosis of rectal cancer and colon adenocarcinoma [[141]36, [142]37]. GNG12-AS1 is a stable lncRNA located in the nucleus, which is involved in the regulation of cell cycle and trans-function of cell migration [[143]38]. Moreover, studies had shown that it can regulate the proliferation and migration of glioma cells through AKT/GSK-3β/ β-catenin signaling pathway [[144]30]. Survival analysis, ROC curve, nomogram and heat map all showed that our prognostic model based on 9 m5C-related lncRNAs had good performance in predicting the prognosis of patients. The model could also accurately distinguish the early and late stage patients. Moreover, both univariate and multivariate analyses showed that the risk score of our model could be used as an independent prognostic factor. Pathway enrichment analysis showed that complement and coagulation cascades, cytokine − cytokine receptor interaction, viral protein interaction with cytokine and cytokine receptor and seven other related pathways were activated. Recent studies had shown that complement C1QC plays an important role in the prognosis and immune infiltration of ccRCC tumors [[145]39]. The cytokine-cytokine receptor interaction plays an important regulatory role in the development of human B lymphocytes, thereby affecting human humoral immunity [[146]40]. Studies had shown that SARS-CoV-2 can down-regulate the expression of angiotensin converting enzyme 2, thereby affecting the prognosis of KIRP(kidney renal papillary cell carcinoma) patients through metabolism and immune regulation [[147]41]. NF − kappa B is involved in the regulation of hypoxia in the tumor microenvironment and plays an important role in the development and maintenance of tumors [[148]42]. As a multifunctional lymphokine, more and more studies had shown that IL-17 has a dual role in tumor development, but its main role in ccRCC seemed to be promoting tumor, which needs further explore [[149]43]. CYP27A1, a key enzyme that regulates cholesterol homeostasis by activating LXR/ABCA1, plays an important role in the proliferation and migration of ccRCC [[150]44]. Immune-related function analysis showed that APC_co_inhibition, APC_co_stimulation, CCR, Check-point, Cytolytic_activity, HLA, Inflammation-promoting, Parainflammation, T_cell_co-inhibition, T_cell_co-stimulation and Type_I_IFN_Reponse were more significantly enriched in the high risk group. However, the enrichment of Type_II_IFN_Reponse in the low risk group was more significant. It has been reported that I_IFN can play a role in anti-viral infection by inducing a variety of antiviral genes, and can participate in the treatment of hepatitis C [[151]45]. In addition to its role in the immune response against tumors, II_IFN is also involved in maintaining the homeostasis of the immune environment in the human body [[152]46]. Survival analysis showed that patients with high risk TMB had a better prognosis, and patients in the high risk group had higher TTN and SETD2 expression. In colorectal cancer, cutaneous melanoma and thyroid cancer, TTN gene mutation is often used as a predictive biomarker to evaluate the immune microenvironment and prognosis [[153]47–[154]49]. SETD2 is a tumor suppressor whose mutation is prevalent in 7% of lung adenocarcinomas [[155]50]. As a mutated chromatin modifier gene, it has the highest mutation rate (13%) in ccRCC [[156]51]. Our study suggests that SETD2 and TTN can be used as targets for cancer immunotherapy. It is important to note that our conclusions are based on the TCGA public database and lack validation of clinical patient specimens and external datasets. We found that GNG12-AS1 may play an important role in ccRCC by activating ERK/GSK-3β/β-catenin signaling pathway, but we have not further verified this result in reversal experiments and animal experiments, and we will further pay attention to and improve these issues in the future. In addition, the expression of DEGs and related lncRNAs needs to be verified by experiments and other methods. Finally, how these lncRNAs interact with m5C and the underlying biological mechanisms by which M5C-related lncRNAs participate in ccRCC remain unclear. Conclusion Our study investigated the role of m5C-related lncRNAs in ccRCC. A prognostic model based on 9 m5C-related lncRNAs was constructed, which could effectively predict the prognosis and immunotherapy response of ccRCC patients, and GNG12-AS1 was a promising prognostic biomarker for predicting survival outcome of ccRCC. GNG12-AS1 may have a significant impact on ccRCC cell migration by activating the ERK/GSK-3β/β-catenin signaling pathway. Acknowledgements