A novel prognostic cancer-related lncRNA signature in papillary renal cell carcinoma

Papillary renal cell carcinoma (pRCC) ranks second in renal cell carcinoma and the prognosis of pRCC remains poor. Here, we aimed to screen and identify a novel prognostic cancer-related lncRNA signature in pRCC. The RNA-seq profile and clinical feature of pRCC cases were downloaded from TCGA database. Significant cancer-related lncRNAs were obtained from the Immlnc database. Differentially expressed cancer-related lncRNAs (DECRLs) in pRCC were screened for further analysis. Cox regression report was implemented to identify prognostic cancer-related lncRNAs and establish a prognostic risk model, and ROC curve analysis was used to evaluate its precision. The correlation between RP11-63A11.1 and clinical characteristics was further analyzed. Finally, the expression level and role of RP11-63A11.1 were studied in vitro. A total of 367 DECRLs were finally screened and 26 prognostic cancer-related lncRNAs were identified. Among them, ten lncRNAs (RP11-573D15.8, LINC01317, RNF144A-AS1, TFAP2A-AS1, LINC00702, GAS6-AS1, RP11-400K9.4, LUCAT1, RP11-63A11.1, and RP11-156L14.1) were independently associated with prognosis of pRCC. These ten lncRNAs were incorporated into a prognostic risk model. In accordance with the median value of the riskscore, pRCC cases were separated into high and low risk groups. Survival analysis indicated that there was a significant difference on overall survival (OS) rate between the two groups. The area under curve (AUC) in different years indicated that the model was of high efficiency in prognosis prediction. RP11-63A11.1 was mainly expressed in renal tissues and it correlated with the tumor stage, T, M, N classifications, OS, PFS, and DSS of pRCC patients. Consistent with the expression in pRCC tissue samples, RP11-63A11.1 was also down-regulated in pRCC cells. More importantly, up-regulation of RP11-63A11.1 attenuated cell survival and induced apoptosis. Ten cancer-related lncRNAs were incorporated into a powerful model for prognosis evaluation. RP11-63A11.1 functioned as a cancer suppressor in pRCC and it might be a potential therapeutic target for treating pRCC.

malignancies and involved in the oncogenesis and progression of human cancers, and they have great potential to serve as biomarkers for cancer diagnosis and targets for cancer treatment [5]. LncRNAs play roles in multiple cellular processes of cancer, including cell survival, migration, invasion, metastasis, and epithelial-mesenchymal transition (EMT), etc. [6,7]. Additionally, lncR-NAs are crucial regulators in many aspects, such as target mRNA expression regulation and cancer immunity regulation [8,9]. What's more, a series of lncRNAs have been identified as cancer-related lncRNAs [10][11][12][13]. For more cancer-related lncRNAs identification, different methods were performed and corresponding database was established [14,15]. Therefore, lncRNAs, especially cancer-related lncRNAs are promising biomolecules in human cancers.
More interestingly, RP11-63A11.1 was a renal tissuespecific lncRNA among the ten lncRNAs. It was reported that RP11-63A11.1 was significantly down-regulated in pRCC samples and correlated with clinicopathological features of pRCC patients [17,18], which was consistent with the analysis of ours. However, the role of RP11-63A11.1 is still not completely understood in pRCC. Thus, we further detected the expression of RP11-63A11.1 in pRCC cells and investigated its role on cell survival and apoptosis, which would add more convincing evidence to previous conclusion.

Acquisition and processing of samples data
The KIRP RNA-seq data (TCGA-KIRP) were acquired from UCSC, and clinical data of pRCC was downloaded from TCGA database. A total of 289 tumor samples and 32 matched non-tumor samples were included based on The Cancer Genome Atlas (TCGA) database. The samples with overall survival (OS) less than 1 month were eliminated, since these patients were likely to die of other reasons, such as serious infection [19,20]. A total of 275 pRCC samples were screened out and divided into training set with 136 cases and testing set with 139 cases by the "caret" package of R software at random. The testing set and the combination set (the 275 samples) were used to validate the outcomes acquired from the training set.

Identification of differentially expressed mRNAs and cancer-related lncRNAs
The differentially expressed mRNAs (DEMs) between pRCC and its corresponding normal samples were analyzed by the "limma" package of R software. The significant cancer-related lncRNAs in pRCC were obtained from the Immlnc database. And then, their expression data in TCGA database was extracted and the "limma" package of R software was performed to analyze the differentially expressed cancer-related lncRNAs (DECRLs). |log2FC| > 2 and adjust P value < 0.01 were considered as the threshold. The DECRLs were removed if their expression values were equal to 0 in more than 135 samples. Finally, 367 DECRLs were screened out for further investigation.

Identification of prognostic cancer-related lncRNAs and establishment of a signature for prognosis prediction
We firstly implemented the univariate Cox regression analysis on the 367 DECRLs and clinical survival data in training set to determine prognostic cancer-related lncRNAs. And then we further analyzed the significant results using multivariate Cox regression analysis to generate a risk model. Each patient obtained a risk score with this formula: riskscore = β1 × exp(lnc1) + β2 × exp (lnc2)+···+βi × exp(lnci). Explnc is the expression value of DECRLs in the patients, while β refers to the regression co-efficient of cancer-related lncRNAs generated from the risk model in training set. Based on the medial riskscore, patients in all three sets were divided into high and low risk groups. Subsequently, the survival rate of the two groups was evaluated. The area under the ROC curve (AUC) was utilized to assess the precision of the risk model. The testing set and combination set were utilized to validate the above results.

Independence analysis of the cancer-related lncRNA signature in pRCC
Both the univariate and multivariate Cox regression were implemented to analyze the independence of this cancer-related lncRNA signature from the clinical features of pRCC patients, including age, gender, tumor stage and T classification. The data of M and N classifications were excluded in the independence analysis because the number of samples with undetermined classification was far beyond half of the all samples.

Weighted gene co-expression analysis (WGCNA) and gene enrichment analysis
In order to explore correlation between lncRNAs signature and biologic functions of pRCC, we performed WGCNA to construct the gene co-expression modules among differentially expressed mRNAs as previous study mentioned [21]. The modules with the maximal absolute module significance correlated to lncRNAs signature were screened out. Then, gene enrichment analysis of genes in the most lncRNAs signature related module was carried out using DAVID [22]. Gene Ontology (GO) terms with false discovery rate (FDR) less than 0.05 were considered as significant biologic functions.

Aberrant expressed mRNAs correlated with RP11-63A11.1
To explore the potential mechanism of RP11-63A11.1 involved in the progression of pRCC, we analyzed the aberrant expressed mRNAs which were correlated with RP11-63A11.1 expression. We listed aberrant expressed mRNAs in the combination set using Pearson method. The positively related mRNAs were screened out with correlation coefficient over 0.3.

Prediction of subcellular location
As the function of lncRNA is associated with its subcellular localization, the sequence of the RP11-63A11.1 was obtained from the LNCipedia website (https:// lncip edia. org/) [23]. The obtained sequence was input into the lncLocator website to predict the subcellular localization of the lncRNA [24]. The results were obtained as a score for each potential subcellular location for RP11-63A11.1, including the cytoplasm, nucleus, ribosome, cytosol and exosome.

The association of RP11-63A11.1 expression and clinicopathological characteristics and prognostic indicators of pRCC
Among the ten lncRNAs, RP11-63A11.1 was primarily expressed in renal tissues and it was the most significant one with the lowest P value for OS in training set. Therefore, we further analyzed the association of RP11-63A11.1 expression and clinicopathological features (including gender, age, tumor stage and TNM classifications) as well as prognostic indicators of pRCC. Furthermore, we also investigated the expression and functions of RP11-63A11.1 in vitro.

Cell culture
Human pRCC SK-RC-39 cell line and normal renal cell line HK-2 were purchased from Cellcook Biotech (Guangzhou, China). The cell mediums (including RPM1-DMEM and RPM1-1640) and Fetal Bovine Serum (FBS) were purchased from Fcmacs Biotech (Nanjing, China). Antibiotic mixture (penicillin and streptomycin) and l-glutamine were purchased from Sangon Biotech (Shanghai, China). 57 mL FBS and 5.7 mL antibiotic mixture as well as 5.7 mL l-glutamine were added into cell mediums to generate complete mediums. SK-RC-39 and HK-2 were fed in complete 1640 and DMEM medium, respectively. Both the two kinds of cell lines were cultured in the 5% CO 2 incubator at 37 °C.

Cell proliferation assay
CCK-8 (Vazyme, Nanjing, China) assay was applied for the evaluation of cell proliferation ability. In brief, after counting with the cell count plate, 1 × 10 4 cells were seeded into each experimental well of a 96-well plate and then cultured at 37 °C for 24 h, 48 h and 72 h. After adding 10% CCK-8 reagent into the corresponding wells, the absorbance value at 450 nm was measured.

Cell apoptosis assay
For cell apoptosis detection, we prepared and adjusted the cell concentration to 1 × 10 6 /mL. Then Annexin V-Alexa Fluor 488/PI Apoptosis detection kit (Fcmacs Biotech, Nanjing, China) was used for cell staining and flow cytometer was used for the assessment of cell apoptosis.

Western Blot
The total protein was extracted using radioimmunoprecipitation assay solution (Beyotime, Shanghai, China) before polypropylene gel electrophoresis. Then the protein was further transferred to polyvinylidene difluoride membranes (Millipore, USA). After blocked by 3% bovine serum albumin (Fcmacs Biotech, Nanjing, China) in room termperature, the membrane was incubated with primary antibody overnight. Sources of primary antibodies (cleaved PARP, cleaved Caspase 3, Bcl-2, Bax, Cytochrome C and GAPDH) were from Cell Signaling (Danvers, MA, USA). RASSF10 antibody was from Abcam (Cambridge, England, UK). After the membrane was incubated with secondary antibody, the proteins were visualized by the gel-imaging system.

Statistical analysis
R software 4.0.2 and SPSS 25 were used to implement statistical analysis in this study. Univariate and multivariate Cox regression analysis were applied to generate a cancer-related lncRNA signature and identify independent prognostic factors for pRCC. The survival rate of pRCC patients in two risk-groups was appraised by Kaplan-Meier method. The P value < 0.05 was generally considered statistically significant.

Differentially expressed mRNAs and cancer-related lncRNAs in pRCC
A total of 1872 cancer-related lncRNAs with expression data in pRCC were extracted (Additional file 1: Table S1). Through the "limma" package of R software, 2030 DEMs were identified, including 219 up-regulated DEMs and 1811 down-regulated DEMs (Fig. 1A, B), while 579 differentially expressed cancer-related lncRNAs (DECRLs) were found, including 84 up-regulated DECRLs and 495 down-regulated DECRLs (Fig. 1C, D). The DECRLs were removed if their expression values were equal to 0 in more than 135 samples. 367 DECRLs were finally determined for further analysis.
According to the riskscore formula, the riskscore of each pRCC patient in all three sets was calculated. On the basis of median riskscore value, pRCC cases were segmented into high and low risk groups. As shown in Fig. 3A, the mortality rate of pRCC patients was progressively increased with the rising of risk score in training-set, testing-set and combination-set. The levels of the ten lncRNAs in two groups were shown in Fig. 3B. The survival curves displayed that patients in high risk group had a shorter survival time than those in low risk group (Fig. 3C). The precision of the 10 cancer-related lncRNA signature in predicting prognosis of pRCC were analyzed through the ROC curves. The results demonstrated that all of the AUC values at 1, 3, and 5 years were over than 0.75 in the training-set, testing-set and the combination-set (Fig. 3D), indicating the ten cancer-related lncRNAs signature had a good sensitivity and specificity in prognosis prediction.

The cancer-related lncRNA signature was an independent factor in prognosis prediction
We further evaluated whether the identified cancerrelated lncRNA signature was an independent factor for prognosis prediction in pRCC. Co-variables including age, gender, tumor stage and the risk model in trainingset, testing-set and combination-set were analyzed by Cox regression analysis. The outcomes from univariate Cox regression revealed that tumor stage and the risk model were significantly relevant to OS in all three sets. Multivariate Cox regression analysis suggested that tumor stage and risk model were independent factors in the combination set. The results were listed in Table 2. Taken together, the cancer-related lncRNA signature was an independent factor for prognosis prediction in pRCC.

Stratification analysis
By performing stratification analysis on the clinicopathological features of pRCC patients in the combination cohort, we revealed that the cancer-related lncRNA signature was a good predictor. As shown in Fig. 4A, B, both older (> 60 years) and younger (≤ 60 years) pRCC cases in high risk group exhibited a worse prognosis compared with those cases in low risk group. Of note, similar results of the cancer-related lncRNA signature were found in different genders (Fig. 4C, D) and tumor stages (Fig. 4E, F).

The cancer-related lncRNA signature was associated with biological processes
The cancer-related lncRNA signature was vital in prognosis of pRCC patients; therefore this signature should be closely associated with biological processes of pRCC. 2030 DEMs were used to construct 8 similar expression modules by average linkage clustering. The relevance with P value between each module and cancer-related lncRNA signature was listed in every module (Fig. 5A). The most negative related module (yellow module, R = − 0.13, P = 0.04) and positive related module (blue module, R = 0.27, P = 5e−06) were selected out for gene enrichment analysis. The FDR of GO term of genes in yellow module was over than 0.05. The genes in blue module were mainly involved in cell division, cell proliferation, and cell cycle process, such as G1/S transition of mitotic cell cycle, sister chromatid cohesion (Fig. 5B), etc.   we found that RP11-63A11.1 was mainly expressed in kidney, suggesting RP11-63A11.1 may be a renal tissue-specific lncRNA (Fig. 6A). It was indicated that the level of RP11-63A11.1 was significantly higher in lower tumor stages and TNM classifications, while it did not show significant difference in age and gender ( Fig. 6B-G). Of note, higher level of RP11-63A11.1 indicates favorable OS (Fig. 7A), PFS (Fig. 7B), as well as DSS (Fig. 7C). Although there is no significant difference of DFS between high expression group and low

Analysis of the potential mechanism of RP11-63A11.1 in pRCC
How the RP11-63A11.1 regulates the apoptosis progression in pRCC cells is still unknown. To further investigate the possible mechanism, we first predict the subcellular location. The sequence information of RP11-63A11.1 was obtained from LNCipedia (Additional file 4: Table. S2). It is indicated in the LncLocator that RP11-63A11.1 is probably located in the cytoplasm ( Table 3), suggesting that RP11-63A11.1 may serve as a ceRNA (competing endogenous RNA). Then, we showed all the genes which were positively correlated with RP11-63A11.1 (Additional file 5: Table S3). Among all the positively correlated genes, we found that RASSF10 was attractive. Finally, we analyzed the expression of RASSF10. Elevated RP11-63A11.1 significantly increased the expression of RASSF10 (Additional file 6: Fig.S3). Thus, we believe that introduction of RP11-63A11.1, as a ceRNA, elevated RASSF10, inducing apoptosis in pRCC cells.

Discussion
Similar to ccRCC, pRCC patients are usually diagnosed at an advanced stage or metastasis. The major treatment for local pRCC is surgical operation. However, approximately 40% patients reoccur after surgical resection, resulting in a worse prognosis [25]. In recent years, increasing researches of cancers focused on molecular characteristic in early diagnosis and prognosis improvement. Of note, lncRNAs have been revealed to function as main regulators in various biological processes [26,27]. Notably, it was believed that cancer-related lncRNAs may be helpful for researchers to understand the mechanism of cancer progression and develop efficient measures to improve the prognosis [28]. And cancer-related lncRNAs have been described in some cancer types, such as esophageal adenocarcinoma [29] and prostate cancer [30]. In the present study, we obtained 1872 caner-related lncR-NAs through integrating Immlnc and TCGA databases and screened 367 significant DECRLs in pRCC. Cox regression analysis suggested that 26 significant DECRLs were related to the survival of patients with pRCC, Based on this signature, survival curves suggested that patients in high risk group had a poorer survival rate than patients in low risk group. Additionally, ROC curves exhibited that the precision of this signature was more than 75% for survival prediction at different years. These results were validated in testing set and combination set and similar outcomes were observed, suggesting the signature was a forceful tool for predicting prognosis. What's more, this cancer-related lncRNA signature was also a forceful tool for prognosis prediction in different classes of some clinicopathological features, including age, gender, and tumor stage.
Other prognostic signatures in pRCC were previously reported. For instance, Wang et al. [31] established a prognostic signature with 15 immune-related genes. Although prognosis of patients could be predicted by their signature, the accuracy of their prognostic signature in predicting prognosis in 1, 3 and 5 years was lower than that of ours (Table 4). In addition, Gao et al. [32] also grouped pRCC patients into training set and testing set. A five mRNAs signature was identified in training set and the accuracy of their signature was 0.82, which was lower than the accuracy in our training set. Although the result of testing set showed that patients could be divided into high and low risk groups by the five mRNAs signature, the accuracy of their signature in testing set was not assessed. Zhang et al. [33] constructed a risk model by 17 mutant genes. The AUC of their signature was 0.907 in 3 years, which was a considerable result. However, the result was not validated in their study and it was unknown whether the mutant-gene signature was efficient for predicting prognosis in 5 years. Similarly, a signature of four lncRNAs was constructed to predict prognosis [34]. In spite of significant accuracy of this signature displayed, a validation result was not provided in their study. Taken together, our cancer-related lncRNA signature might be more beneficial than the previous signatures mainly because of higher accuracy and validation of the findings.
Among the lncRNAs in our cancer-related signature, RNF144A-AS1, now called GRASLND, was reported to serve as an important regulator in stem cell chondrogenesis [35]. RNF144A-AS1 could also facilitate the migration and invasion of bladder cancer cells [36]. LINC00702 was a newly identified lncRNA and was involved in the progression of several malignancies via tumorigenesisassociated pathways such as Wnt/β-catenin pathway and PTEN/PI3K-AKT pathway [37][38][39][40]. LUCAT1 was the only lncRNA widely investigated in human cancers among the ten lncRNAs. LUCAT1 could promote tumorigenesis and development in various cancers. Besides, it could promote anti-tumor drug resistance in some tumor types such as NSCLC and osteosarcoma [41,42]. These results suggest that these lncRNAs are reasonable and of importance in human cancers. However, to our knowledge, they were firstly uncovered to be novel prognostic biomarkers in pRCC in this study. In consistent with previous researches [17,18] in which RP11-63A11.1 was reported, we found that RP11-63A11.1 correlated with clinicopathological features and prognosis of pRCC patients. However, its role in pRCC cells was incompletely understood. We revealed that RP11-63A11.1 was decreased in pRCC cells. Furthermore, increased RP11-63A11.1 inhibited the proliferation and induced apoptosis of pRCC cells, indicating RP11-63A11.1 served as a tumor suppressor in pRCC. Subcellular localization of lncRNA transcripts is of great importance when we investigate its mechanism. Therefore, the subcellular location of the RP11-63A11.1 transcript was predicted  using lncLocator based on its sequence acquired from LNCipedia. RP11-63A11.1 was likely to be located in the cytoplasm. LncRNA located in the cytoplasm or cytosol probably serves as a ceRNA, and regulates mRNA stability or translation. Regarding that elevated RP11-63A11.1 significantly increased the expression of RASSF10, we hypothesized that RP11-63A11.1 acted as ceRNA, and regulated the stability or translation of RASSF10. Further studies are indicated to investigate more specific mechanisms of RP11-63A11.1 in pRCC. The roles of other lncRNAs including RP11-573D15.8, LINC01317, RP11-400K9.4, and RP11-156L14.1 were not reported in previous literatures. These lncRNAs were firstly reported in the present study. When combined with the above-mentioned lncRNAs, these lncRNAs generate a powerful tool for survival prediction in pRCC. To explore the biological functions of the cancer-related lncRNA signature in pRCC, we performed WGCNA to seek gene modules associated with this signature. Two gene modules were finally identified. Gene enrichment analysis showed genes in the blue module (positively correlated) were mainly enriched in biologic processes of cell division, proliferation and cell cycle, indicating that pRCC with high risk score has stronger proliferation capability than pRCC with low risk score.

Conclusions
In summary, we successfully figured out a novel cancerrelated lncRNA signature with powerful predictive function for pRCC prognosis. These results of this study may facilitate understanding of molecular mechanisms in pRCC development, and the cancer-related lncRNA signature may become a promising biomarker and therapeutic target in pRCC. Of note, RP11-63A11.1 may serve as a ceRNA and regulate RASSF10, facilitating the understanding of the development in pRCC.