Novel Multiple miRNA-Based Signatures for Predicting Overall Survival and Recurrence-Free Survival of Colorectal Cancer Patients.

BACKGROUND Colorectal cancer (CRC) has become a heavy health burden around the world, accounting for about 10% of newly diagnosed cancer cases. In the present study, we aimed to establish the miRNA-based prediction signature to assess the prognosis of CRC patients. MATERIAL AND METHODS A total of 451 CRC patients' expression profiles and clinical information were download from the TCGA database. LASSO Cox regression was conducted to construct the overall survival (OS)- and recurrence-free survival (RFS)-associated prediction signatures, by which CRC patients were divided into low- and high-risk groups. Kaplan-Meier (K-M) curve and receiver operating characteristic (ROC) curves were used to explore the discriminatory ability and stability of the signatures. Functional enrichment analyses were performed to identify the probable mechanisms. RESULTS miRNA-216a, miRNA-887, miRNA-376b, and miRNA-891a were used to build the prediction formula associated with OS, while miR-1343, miR-149, miR-181a-1, miR-217, miR-3130-1, miR-378a, miR-542, miR-6716, miR-7-3, miR-7702, miR-677, and miR-891a were obtained to construct the formula related to RFS. K-M curve and ROC curve revealed the good discrimination and efficiency of OS in the training (P<0.001, AUC=0.712) and validation cohorts (P=0.019, AUC=0.657), as well as the results of RFS in the training (P<0.001, AUC=0.714) and validation cohorts (P=0.042, AUC=0.651). The function annotations for the targeted genes of these miRNAs show the potential mechanisms of CRC. CONCLUSIONS We established 2 novel miRNA-based prediction signatures of OS and RFS, which are reliable tools to assess the prognosis of CRC patients.

Background

Colorectal cancer (CRC) is a serious health problem worldwide. It is the third most prevalence cancer, accounting for about 10% of new cancer cases, and is the second most common cause of cancer-related deaths [1]. In China, CRC is the third most common cancer and is the fifth most common cause of cancer-related deaths. There are an estimated 376.3 newly diagnosed CRC cases per 100 000 population every year [2]. As in other Asian countries, CRC is third most prevalent cancer in males and the second most prevalent cancer in females in Japan [3]. In Europe, CRC is the second most common cancer, and the mortality rate is decreasing both in men (−6.7%) and women (−7.5%) from 2012 to 2018 [4]. In the USA, CRC is the third most common cancer, and the 5-year relative survival rate is 80.1% to 88.1% in the patients with stage I-II, and is only 12.6% in stage IV patients [5]. Although the number of new CRC diagnoses and deaths has sharply fallen in the USA, more and more CRC patients younger than 50 years old were diagnosed in past decades, and there has been an increase in early death rates since 2000 [6]. Thus, for CRC patients, an efficient prediction biomarker is essential to predict prognosis.

MicroRNAs (miRNAs) are a class of small non-coding RNAs that regulate the protein-encoded gene after transcription by binding to the 3′-UTR of mRNA, and about half of mRNAs are thus regulated by miRNAs [7,8]. Recently, more and more miRNAs are found to have pivotal roles in cancer cell proliferation, differentiation, or apoptosis [9-11]. Thus, miRNAs could be involved in tumorigenesis or recurrence as oncogenes or suppressors by altering cell signaling pathways [12,13]. With the development of gene sequencing, miRNAs are seen as potential biomarkers, not only in early-stage detection of cancer, but also in the prediction of prognosis [14-17].

In the current study, we established the miRNA-based prediction signature to evaluate the death risk and recurrence risk for CRC patients, based on the miRNA matrix and clinical information of TCGA CRC patients.

Material and Methods

CRC patients’ database download

We obtained CRC datasets from the TCGA database [18], which contains miRNA expression data of 451 CRC patients, as well as overall survival (OS) and recurrence-free survival (RFS) information. We randomly divided these 451 CRC patients into a training cohort and a validation cohort.

MiRNAs candidate screening and establishment of prediction signature

Univariate Cox analysis of each single miRNA with OS and RFS was evaluated in the training cohort to determine all the potential miRNAs. The formula of miRNA-based prediction signatures was establishment based on LASSO Cox regression analyses, according to the hazards ratio (HR) and co-efficient (co-ef). R (v3.3.1) was used to complete all the analysis, and we used the “glmnet” package (2.0–10) for LASSO analysis.

Kaplan-Meier (K-M) curve and receiver operating characteristic (ROC) curve

We obtained the risk scores of all the patients in the training cohort and validation cohort, with the pre-established OS- and RFS-related miRNA prediction signatures, respectively. In each cohort, CRC patients were also divided into low-risk and high-risk groups to complete the following evaluation, based on the risk score <0 or >0. We measured the difference in OS between the high-risk and low-risk patients stratified by the OS-related miRNA prediction model, and also analyzed the difference in RFS between the high-risk and low-risk patients stratified by the RFS-related miRNA prediction model. ROC curves and the area under the ROC curve (AUC) were also calculated to assess the discriminatory ability and stability of the OS- and RFS-related miRNAs prediction signatures.

Target prediction and function annotation

To further explore the function of miRNAs enrolled in the CRC prognosis prediction signatures, we predicted the downstream correlated genes with the web-interactive prediction tool, TargetScan [19]. For all the enrolled downstream genes, we managed the pathway annotation to screen the hypothetical biological pathways involved in the OS- or RFS- related positive prediction miRNAs. GO ontology and KEGG pathway analyses were conducted by DAVID (P value 0.05). Visualization of these miRNAs and genes association was performed using Cytoscape software (Cytoscape Consortium, San Diego, CA, USA).

Results

Establishment of the OS- and RFS- prediction signature for CRC patients

Univariable Cox analysis and multiple LASSO Cox regression analysis were conducted to choose the most OS- and RFS-related miRNAs among CRC patients. With the miRNA data matrix of the TCGA database, we screened all the miRNAs that were positively associated with the OS and RFS of CRC patients, respectively. Finally, there were 27 miRNAs obtained for establishing the prediction signature of CRC patients’ OS (Supplementary Figure 1A). Simultaneously, another 19 miRNAs reflecting the process of RFS were also acquired (Supplementary Figure 1B).

Using LASSO Cox regression analysis, we drew out the prediction signature-involved miRNAs to evaluate OS and RFS. For predicting the OS of CRC patients, we found that miR-216a (HR=1.620, 95% CI=1.190–2.190, P=0.002, co-ef=0.480), miR-887 (HR=1.990, 95% CI=1.330–2.970, P<0.001, co-ef=0.690), miR-376b (HR=0.520, 95% CI=0.370–0.740, P<0.001, co-ef=−0.660), miR-891a (HR=1.200, 95% CI=0.960–1.500, P=0.108, co-ef=0.180) (Figure 1A), and the risk score formula of the OS for each CRC patient was risk score=0.480×miR-216a+0.690×miR-887−0.660×miR-376b+0.180×miR-891a (Supplementary Table 1). For predicting the RFS of CRC patients, we extracted miR-1343 (HR=1.64, 95% CI=0.88–3.1, P=0.119, co-ef=0.495), miR-149 (HR=1.06, 95% CI=0.79–1.4, P=0.703, co-ef=0.060), miR-181a-1 (HR=1, 95% CI=0.65–1.5, P=0.995, co-ef=0.001), miR-217 (HR=1.16, 95% CI=0.86–1.6, P=0.329, co-ef=0.150), miR-3130-1 (HR=0.92, 95% CI=0.67–1.3, P=0.595, co-ef=−0.087), miR-378a (HR=0.81, 95% CI=0.55–1.2, P=0.285, co-ef=−0.210), miR-542 (HR=1.22, 95% CI=0.82–1.8, P=0.337, co-ef=0.194), miR-6716 (HR=1.42, 95% CI=0.74–2.7, P=0.295, co-ef=0.359), miR-7-3 (HR=0.83, 95% CI=0.6–1.1, P=0.239, co-ef=−0.192), miR-7702 (HR=0.77, 95% CI=0.52–1.1, P=0.194, co-ef=−0.258), miR-677 (HR=1.26, 95% CI=0.89–1.8, P=0.193, co-ef=0.233), miR-891a (HR=1.44, 95% CI=1.09–1.9, P=0.01, co-ef=0.362), (Figure 1B). The risk score formula of the OS for each CRC patient is risk score=0.495×miR-1343+0.06×miR-149+0.001×miR-181a-1+0.15×miR-217−0.087×miR-3130-1−0.21×miR-378a+0.194×miR-542+0.359×miR-6716−0.192×miR-7-3−0.258×miR-7702+0.233×miR-677+0.362×miR-891a (Supplementary Table 2). Along with the signature prediction formulas, the risk scores of OS and RFS for each patient in the training cohort and validation cohort were assessed.

Figure 1

Establishment of the multiple miRNA prediction signatures associated with OS and RFS of colorectal cancer. Hazard ratio of the enrolled OS- (A) and RFS- (B) related miRNA conducted by LASSO Cox regression analysis.

Discrimination of the miRNAs prediction signature in the training and validation cohorts

To evaluate the discrimination of the miRNA-based prediction signatures, we managed the K-M curve to compare the OS or RFS in different risk groups. For OS prediction signature, patients with low-risk scores presented much better OS time compare to the high-risk group of the training cohort (P<0.001) (Figure 2A), as well as in the validation group (P=0.019) (Figure 2C). These results indicated that the prediction signatures based on miR-216a, miR-887, miR-376b, and miR-819a could reveal the prognosis of OS time for CRC patients. For the RFS prediction signature, patients with low risk scores had much better recurrence-free survival times compare to the high-risk group of the training cohort (P<0.001) (Figure 3A), as well as in the validation group (P=0.0011) (Figure 3C). These results indicated that the 12-miRNA-based prediction signature could predict the RFS time for CRC patients.

Figure 2

OS-related miRNA prediction signature performance in colorectal cancer patients. Kaplan-Meier curves of the low- and high-risk groups divided by the OS prediction signature in the training cohort (A) and validation cohort (C); ROC curves of the low- and high-risk groups divided by the OS prediction signature in the training cohort (B) and validation cohort (D).

Figure 3

RFS-related miRNA prediction signature performance in colorectal cancer patients. Kaplan-Meier curves of the low- and high-risk groups divided by the RFS-related miRNA prediction signature in the training cohort (A), and validation cohort (C); ROC curves of the low- and high-risk groups divided by 4 RFS-related miRNA prediction signatures in the training cohort (B) and validation cohort (D).

Efficiency of the miRNAs prediction signature in the training and validation cohorts

ROC curve and AUC were used to assess the discriminatory ability and stability for OS or RFS prediction signature in training and validation cohorts. For the OS prediction signature, the dependent variable was whether the CRC patient was alive or not. The ROC curve in the training cohort shown an AUC of 0.712, with 95% CI of 0.637–0.788 (Figure 2B), while the AUC was 0.657 with 95% CI of 0.462–0.853 (Figure 2D) for the validation cohort. These results disclosed that the 4 miRNAs-based prediction signatures precisely and steadily assessed the OS risk for CRC patients. For RFS prediction signature, the dependent variable was whether the CRC was recurrent or not. The ROC curve in the training cohort showed the AUC was 0.714 with 95% CI of 0.611–0.818 (Figure 3B), while it was 0.651 with 95% CI of 0.502–0.800 (Figure 3D) for the validation cohort. These results disclosed that the 4 miRNAs-based prediction signatures precisely assessed judge the RFS of CRC patients.

Function annotation of the downstream genes for OS- and RFS-related miRNAs

We predicted the downstream correlated genes using the TargetScan web interactive prediction tool [19]. The networks of downstream genes for the 4 OS-related miRNAs are displayed in Figure 4A visualized with Cytoscape, and the downstream genes of the RFS related miRNAs are shown in Figure 4B. Pathway annotation was performed to discover the biological pathways involved in the progression of CRC (Figure 5, Supplementary Figure 2). For the target genes predicted by OS-related miRNA signature (Figure 5), GO-BP-Enrich indicated the most significantly pathway, including regulation of cell morphogenesis and protein dephosphorylation; the GO-CC-Enrich items were transcriptional repressor complex, ubiquitin ligase complex, and cell-cell junction, while GO-MF-Erich is involved in ubiquitin-associated enzyme activity, phosphoric ester hydrolase activity, and transcriptional activator activity. KEGG pathway enrichment was shown to participate in ErbB signaling pathway, Proteoglycans in cancer, MAPK signaling pathway, Ras signaling pathway, AMPK signaling pathway, Focal adhesion, Regulation of actin cytoskeleton, and EGFR tyrosine kinase inhibitor resistance. For the downstream genes linked to RFS-related miRNAs, similar results were obtained (Supplementary Figure 2). All these results show that the downstream target genes in the prediction bioscience pathways are related to CRC initiation, progression, and drug resistance.

Figure 4

miRNA-mRNA interaction network. (A) The downstream gene prediction of OS related miRNAs; (B) downstream gene prediction of RFS-related miRNAs.

Figure 5

Functional enrichment analysis depicted the biological pathways and processes associated with OS-correlated genes. The results of GO-BP biological process enrichment (A), GO-CC biological process enrichment (B), GO-MF biological process enrichment (C), Hallmark biological process enrichment (D), and KEGG signaling pathways analysis (E), Reactome biological process enrichment (F).

The combined prognosis value of the miRNA-based classifier and clinical factors

We assessed the combined prognostic values of the miRNA classifiers and clinical factors. For OS, in the multiple nomogram analysis of HR of the factors, we found that age over 60 years old (P=0.0093) and the OS classifier (P=0.00026) are the independent factors predicting risk of death (Supplementary Table 3). The ROC curve revealed that the multiple nomogram (AUC=0.680, 95% CI=0.623–0.736), pathologic stage (AUC=0.722, 95% CI=0.650–0.794), and the miRNA-based classifier (AUC=0.644, 95% CI=0.577–0.711) were all good tools to predict the survival status of CRC patients (Figure 6A).

Figure 6

The comparison of the 2 classifiers and clinicopathological features. ROC curve identified the differences between the miRNA-based OS (A) and RFS (B) classifiers and clinicopathological classifiers in the overall cohort.

For RFS, 12-based predicting miRNAs-based classifier precisely predicted the recurrence status. In the multiple nomogram analysis of hazard ratio of the factors, we found that male sex (P=0.004779), Stage IV (P=0.04232), and the miRNA classifier (P<0.001) are independent factors increasing the risk of death (Supplementary Table 4). The ROC curve revealed that the multiple nomogram (AUC=0.686, 95% CI=0.616–0.755), the miRNA classifier (AUC=0.664, 95% CI=0.592–0.736), and pathologic stage (AUC=0.670, 95% CI=0.589–0.752) are all good tools to predict the recurrence status of HCC patients (Figure 6B). Moreover, details of the clinicopathological features of these patients were shown in Supplementary Table 5.

Discussion

CRC is a common cancer with the high incidence and can cause a high rate of cancer-related death, of which the 5-year survival rate is about 64–67% [20]. Surgery is still the criterion standard treatment for early or even advanced CRC [21,22], and 30–50% of recurrence of CRC after surgery occurs within the first 2 years [23]. It is essential to find efficient prediction tools to predict the prognosis for each individual patient, aiming to provide timely and precise therapy.

Several miRNAs have been reported to be associated with the CRC progression or recurrence. In a clinical study, Baltruskevicienee et al. [24] revealed that miR-148a and miR-625-30 were downregulated in CRC, while patients with high expression of miR-148a had shorter RFS times. Takahashi et al. [25] observed a similar phenomenon among advanced CRC patients, showing that low miR-148a expression was correlated with a remarkably shorter disease-free survival time and indicated a poor OS. Hibino et al. [26] found that miR-148a could promote the invasion of CRC through MMP7. Ashizawa et al. [27] revealed that miR-148a-3p could negatively regulate the expression of PD-L1 in colorectal cancer cells, and further the immunosuppressive tumor microenvironment. Zhuang et al. [28] demonstrated that miR-106b-5p is a suppressor of CRC through the MALAT1/miR-106b-5p/SLAIN2 signaling pathway. Chen et al. [29] discovered that miR-203a-3p promotes the proliferation and migration of CRC through its target gene, PDE4D.

In the current study, we established the miRNA-based prediction signature of OS and RFS for CRC patients. For the OS-related miRNA prediction signature, miRNA-216a, miRNA-887, miRNA-376b, and miRNA-891a were used to build the prediction formula, and the predicted OS of each patient was calculated with the formula, then we divided the patients into a low risk of death group and a high risk of death group among the 2 cohorts (training and validation cohorts). For RFS-related miRNA prediction signature, we used miR-1343, miR-149, miR-181a-1, miR-217, miR-3130-1, miR-378a, miR-542, miR-6716, miR-7-3, miR-7702, miR-677, and miR-891a to construct the formula, and divided the patients with low and high recurrence risk in the training and validation cohorts. The results showed that the novel miRNA prediction signature of OS could precisely and reliably predict the OS of CRC patients, and the RFS prediction formula had the similar results. Several previous publications support it effective prediction of CRC prognosis. Zhang et al. [30] revealed that miR-216a could suppress the function of KIAA1199, and subsequently decreased invasion in vitro and metastasis in vivo. Qiu et al. [31] discovered that the incidence of miR-376b variant is higher in tumor tissue than in adjacent normal tissue, and this variant could influence the target genes of miR-376b. No previous studies have assessed how the other miRNAs are involved in the development and prognosis of CRC. The functional annotation of OS- and RFS-related prediction miRNAs and their downstream genes showed the potential mechanisms of CRC. For example, the AMPK signaling pathway and Ras signaling pathway were associated with the OS of CRC reveled by KEGG analysis, and these results have been verified by other researchers [32,33].

Conclusions

Overall, we established 2 novel miRNA prediction signatures of OS and RFS for CRC patients, which successfully classify the CRC patients into low- and high-risk groups, as well as reveal the risk of death and recurrence for CRC patients. These 2 novel miRNA signatures are reliable tools for use in assessing the prognosis of CRC patients.

Kaplan-Meier curves of survival-associated miRNA detected with univariable Cox regression analysis. (A) Overall survival-related miRNAs; (B) Recurrence-free survival-related miRNAs.

Functional enrichment analysis depicted the biological pathways and processes associated with RFS-correlated genes. (A) The results of GO-BP biological process enrichment. (B) GO-MF biological process enrichment.

Supplementary Table 1.

Cox regression analyses were performed on the training data to determine the coefficient of the OS-related miRNAs.

miR-ID	Co-ef	Exp (co-ef)	Se (co-ef)	z	Pr (> |z|)
hsa-mir-216a	0.481361978	1.618276959	0.155460269	3.096366555	0.00195908
hsa-mir-887	0.686677453	1.987102312	0.205054711	3.348752382	0.000811763
hsa-mir-376b	−0.656208945	0.518814465	0.177832824	−3.690032757	0.000224225
hsa-mir-891a	0.182747919	1.200511744	0.113621665	1.608389737	0.107749848

Co-ef – co-efficient; Exp (co-ef) – expected (co-ef); Se (co-ef) – standard error (co-ef).

Supplementary Table 2.

Cox regression analyses were performed on the training data to determine the coefficient of the RFS-related miRNAs.

miR-ID	Co-ef	Exp (co-ef)	Se (co-ef)	z	Pr (> |z|)
hsa-mir-1343	0.495979821	1.642106422	0.317875979	1.560293493	0.11869054
hsa-mir-149	0.056549768	1.058179277	0.148141377	0.381728382	0.702662845
hsa-mir-181a-1	0.001436688	1.001437721	0.220095953	0.006527553	0.994791803
hsa-mir-217	0.150055104	1.161898266	0.153541186	0.977295455	0.328422902
hsa-mir-3130-1	−0.086633684	0.917012951	0.162530614	−0.533029941	0.594012854
hsa-mir-378a	−0.210284456	0.810353704	0.196674069	−1.069202751	0.284978319
hsa-mir-542	0.194905828	1.215196543	0.202832126	0.960921877	0.33659145
hsa-mir-6716	0.348947533	1.417574812	0.332736581	1.048720075	0.294306972
hsa-mir-7-3	−0.192239893	0.825108906	0.163030248	−1.179167029	0.238331672
hsa-mir-7702	−0.25760838	0.772897855	0.197932617	−1.301495349	0.193088956
hsa-mir-877	0.232934145	1.262298348	0.178909745	1.301964545	0.192928506
hsa-mir-891a	0.361565683	1.435575312	0.140087772	2.580993896	0.009851632

Co-ef – co-efficient; Exp (co-ef) – expected (co-ef); Se (co-ef) – standard error (co-ef).

Supplementary Table 3.

Cox regression analyses of OS-related miRNA signature and clinical features was performed to evaluate the coefficient.

		Co-ef	Exp (co-ef)	Se (co-ef)	z	Pr (> |z|)
Sex	Female	Reference	–	–	–	–
	Male	−0.05395	0.947477	0.216541	−0.24916	0.803238
Age	<60	Reference	–	–	–	–
	≥60	0.678291	1.970507	0.260971	2.599106	0.009347
Stage	I	Reference	–	–	–	–
	II	0.338812	1.403279	1.158715	0.292403	0.769979
	III	1.05351	2.8677	1.123976	0.937307	0.348601
	IV	1.926002	6.862021	1.154961	1.66759	0.095397
T stage	T1+T2	Reference	–	–	–	–
	T3	0.277453	1.319764	1.027374	0.270061	0.787114
	T4	1.054186	2.869638	1.048056	1.005848	0.314489
Classifier	Low risk	Reference	–	–	–	–
	High risk	0.813539	2.255877	0.222891	3.649941	0.000262

Co-ef – co-efficient; Exp (co-ef) – expected (co-ef); Se (co-ef) – standard error (co-ef).

Supplementary Table 4.

Cox regression analyses of RFS-related miRNA signature and clinical features was performed to evaluate the coefficient.

		Co-ef	Exp (co-ef)	Se (co-ef)	z	Pr (> |z|)
Sex	Female	Reference	–	–	–	–
	Male	0.727198	2.069274	0.257727	2.821585	0.004779
Age	<60	Reference	–	–	–	–
	≥60	−0.19094	0.826182	0.265596	−0.71891	0.472195
Stage	I	Reference	–	–	–	–
	II	1.200736	3.32256	0.888132	1.351978	0.176382
	III	1.524229	4.5916	0.840919	1.812574	0.069898
	IV	1.823177	6.191498	0.897959	2.030356	0.04232
T stage	T1+T2	Reference	–	–	–	–
	T3	−0.67801	0.507626	0.746377	−0.9084	0.363666
	T4	−0.03927	0.961495	0.785657	−0.04998	0.960139
Classifier	Low risk	Reference	–	–	–	–
	High risk	1.382055	3.983077	0.288713	4.786954	<0.001

Co-ef – co-efficient; Exp (co-ef) – expected (co-ef); Se (co-ef) – standard error (co-ef).

Supplementary Table 5.

The clinical features of all enrolled CRC patients.

	Male (N=224)	Female (N=200)
Age (years)	67.4 (31–90)	65.6 (34–90)
Alive
Yes	170	155
No	51	44
Recurrence
Yes	145	30
No	49	139
Pathologic stage
I	39	32
II	85	76
III	61	61
IV	33	27
Pathologic T stage
T1	4	8
T2	41	31
T3	152	139
T4	27	22

Four patients lacked data on whether they were still alive, 61 patients lacked data on whether there was recurrence, and 10 patients lacked the data on pathologic stage.