As one of the most common tumors of urinary system, bladder cancer carried a poor prognosis with hundreds of thousands of new cases and deaths every year[2]. Currently, clinician usually take the TNM staging system as the most important standard for prognostic predictions for BC patients, which was developed by the American Joint Committee on Cancer [18]. However, the TNM system distinguishes different cancer stages by anatomical disease progression. It is difficult to explain why patients with same TNM stage exhibit variable responses to the similar therapy. This phenomenon prompts us to pay attention on the genomic heterogeneity of malignant lesions. Therefore, we aimed to screen several effective biomarkers for earlier BC diagnosis and prognosis prediction.
According to the accumulating evidence, lncRNAs are involved in tumorigenesis and cancer progression. LncRNAs participate in cell proliferation, apoptosis, invasion, metastasis and immune response [19–21]. Without function of protein coding, lncRNAs are attracting as miRNA sponges or ceRNAs to regulate post-transcriptional processes or gene expression indirectly [22]. In recent years, various ceRNA networks have been reported in which lncRNAs regulate mRNA expression by interacting with miRNAs. As a miRNA sponge distributed mostly in the cytoplasm, the up-regulated expression of lncRNA DANCR is detected in BC. It positively regulates the expression of MSI2 by adsorbing miR-149, and thus promotes BC cell proliferation, migration, invasion and tumorigenicity [23]. HOTAIR is a well-known lncRNA involved in various tumors. In breast cancer, the high expression of HOTAIR inhibits cell apoptosis, facilitates cell growth and metastasis through miR-20a-5p/HMGA2 axis [24]. What’s more, due to the initial status in ceRNA network, lncRNAs can be taken as potential diagnostic or prognostic biomarkers for cancer patients. In the hepatocellular carcinoma(HCC), a ceRNA network with 37 lncRNAs, 10 miRNAs, and 26 mRNAs included was constructed based on TCGA data. And then by means of multivariate Cox regression and lasso analysis, a risk score system consisting of 13 lncRNAs was established for prognositic prediction of patients with HCC [25]. Nevertheless, the ceRNA network of BC progression and the prognostic values of lncRNA signature in the ceRNA network have been rarely investigated. Therefore, the study of RNA interactions in BC progression will help us to a better understand the biological processes in tumor.
In the present study, RNA sequencing data of BC was firstly extracted from TGGA database. We randomly devided the BC patients into two group, training set and validation set. Differentially expressed lncRNAs, miRNAs, and mRNAs were screened by comparing early tumor stages (I/II) tissues with normal tissues, late tumor stages (III/IV) with normal tissues respectively according to data of training set. Then the overlapping 130 lncRNAs, 159 miRNAs, and 2,048 mRNAs in were identified. We screened the lncRNA-miRNA-mRNA chain as its expression pattern is up-down-up-regulation or down-up-down-regulation in BC tissues after predicting the interactions between lncRNA and miRNA, miRNA and mRNA. Then a BC-associated ceRNA regulatory network including 70 DElncRNAs, 30 DEmiRNAs, and 62 DEmRNAs was constructed. Nextly, the functions and signaling pathways in which these dysregulated mRNAs were involved were determined by means of GO and KEGG analyses. We performed K-M survival analysis to dig out the RNAs in this network which were significantly associated with OS. Subsequently, we established a 7-lncRNA-based risk prediction model (AC004803.1, AC009690.2, AC055713.1, AC105942.1, BX890604.1, MIR200CHG, RNF139-AS1) via univariate Cox regression and LASSO regression analysis. According to the risk score, the patients were devided into low-risk and high-risk groups in the K-M survival analysis. For the time-dependent ROC curve, the AUC values of one year, three years and five years were 0.665, 0.675 and 0.661 respectively, showing moderate efficacy in prognostic prediction. Then we verified the prognostic value of the mutiple lncRNAs signature in the validation set. We also conducted univariate and multivariate Cox regression analyses to combine the risk score system and clinical characteristics, and the result suggested a good predictive power of the risk score, which can be an independent factor for prognosis. Finally, based on ceRNA network theory, the correlation between lncRNAs and mRNAs was explored.
Through GO and pathway analyses, we can explore the functions of dysregulated mRNAs in the ceRNA network. Several GO terms of the differentially expressed mRNAs belonged predominantly to the following categories: “epithelial cell proliferation”, “mesenchymal cell proliferation”, “cell cycle G1/S phase transition”, “regulation of cell cycle G1/S phase transition” and “lymphocyte differentiation”, suggesting that the alteration of cell cycle, cell proliferation and immune response may play crucial roles in BC progression. The result of KEGG pathway analysis showed that mRNAs were most enriched in “microRNA in cancer”, “breast cancer”, “gastric cancer”, “hepatitis B”, “proteoglycans in cancer”, “signaling pathways regulating pluripotency of stem cells”, “cell senescene”, “melanoma”, “colorectal cancer” and “prostate cancer”. This indicated that different cancer types may share the common abnormal signaling pathways.
Among the 7 lncRNAs involved in the signature model, limited studies were performed on the role of MIR200CHG. Previously, an up-regulated expression of MIR200CHG was demonstrated in colon adenocarcinoma. MIR200CHG was identified as an immune-related lncRNA and independent biomarker for prognositic prediction of patients with colon adenocarcinoma[26]. MIR200CHG was also included in a prognostic four-lncRNA based Risk Score evaluation model for urothelial bladder carcinoma. The up-regulated expression of MIR200CHG was observed in tumor than normal tissues in bladder cancer. However, the Risk Score system showed that MIR200CHG acts as a protective factor in BC prognosis [27]. The result may be contributed to the various functional mechanisms of lncRNAs, which regulated expression of multiple mRNAs by competitively bind miRNAs. In our study, the upregulated expression of MIR200CHG and its protective action in BC progression was consistent with the study above. This indicates that different regulatory function of MIR200CHG may be involved in BC tumorigenesis and progression.
In another study, Hang et al.[28] selected AC105942.1 as an EMT-related lncRNA of which the expression is negatively related to BC prognosis. Consistently, AC105942.1 is a risky factor in BC prognosis in our study. According to the correlation test, AC105942.1 showed high positive correlations with the RECK mRNA via miR-182-5p, and STK40 mRNA via miR-31-5p. No reports that address the roles of AC004803.1, AC009690.2, AC055713.1, BX890604.1 and RNF139-AS1, but these lncRNAs were speculated to be involved in BC oncogenesis and progression.
Several limitations of our study remain to be solved. First of all, the sequencing data of a small amount of tumor-free normal bladder tissues from TCGA may cause some error for DEGs identification. More samples of control group should be included in the study to pledge the accuracy of analysis. Second, external validation was not performed. The valuable of the risk model could be proved with more data available in other databases. Third, the functions of most lncRNAs in the risk model remain unclear because of the lack of relevant literature. And thus, the functions of unreported lncRNAs should be evaluated via a series of exploratory experiments.