Identification and Validation of a Potent Multi-miRNA Signature for Prediction of Prognosis of Osteosarcoma Patients.

BACKGROUND Osteosarcoma, the most common solid malignancy, has high incidence and mortality rates. We constructed a miRNA-based signature that can be used to assess the prognosis of osteosarcoma patients. MATERIAL AND METHODS The miRNA profile was derived from the Gene Expression Omnibus (GEO) website, with matched clinical records. The miRNA-based overall survival (OS)-predicting signature was established by LASSO Cox regression analysis. Receiver operating characteristic (ROC) curve and Kaplan-Meier (K-M) analyses were performed to examine the stability and discriminatory ability of the OS-predicting signatures. Pathway enrichment analyses were performed to uncover potential mechanisms. RESULTS Three miRNAs (miR-153, miR-212, and miR-591) independently related to the OS were extracted to build a risk score formula. The ROC curve and K-M analyses revealed good discrimination ability of the OS signature for osteosarcoma patients in both the training cohort (P=0.00015, AUC=0.962) and the validation cohort (P=0.0065, AUC=0.793). As shown in multivariate analysis, the classifier showed favorable predictive accuracy similar to the recurrence status to be an independent risk factor for osteosarcoma. Furthermore, the nomogram showed a synergistic effect by combining the clinicopathological features with our classifier. Also, the enrichment analyses of the target genes may contribute to improved treatment of osteosarcoma. CONCLUSIONS The 3-miRNA-based classifier serves as an effective prognosis-predicting signature for osteosarcoma patients.

Background

Osteosarcoma, one of the most common primary malignancies, has the capacity to produce immature osteoid and bones and mostly occurs in children and adolescents [1-3]. Osteosarcoma accounts for 3–5% of newly diagnosed cancers in children [4], and there are approximately 2 million osteosarcoma cases worldwide [5]. Although the 5-year survival of osteosarcoma patients has improved dramatically and is up to 60–70% [6-8], as adjuvant chemotherapy and various surgical methods advance, patients with metastasis, especially to the lungs, have an extremely low survival rate of 19% and high mortality [8,9]. Micro-metastasis is estimated to have occurred in 60% of osteosarcoma patients by the initial time of examination [10]. Furthermore, over 1/3 of patients experience local relapse and distant metastases after surgery and potent chemotherapy, and the mean survival of these patients is less than 1 year [11]. The need to develop early diagnosis and accurate prognosis prediction systems to improve clinical outcomes in patients with osteosarcoma is urgent.

Recently, microRNAs (miRNAs or miRs), which are widely present in eukaryotes [12], have been exploited as potential prognostic biomarkers and therapeutic targets for a variety of cancers. They regulate gene expression at the posttranscriptional level by translational inhibition and degradation [13]. miRNAs can act as either tumor suppressors or oncogenes by regulating the targeted genes, participating in tumor initiation and progression processes [14,15]. To date, multiple studies have shown the association of differential expression of miRNAs with the initiation and etiology of osteosarcoma, and upregulated and downregulated miRNAs were found in human osteosarcoma tissues or osteosarcoma cell lines [16-18]. Upregulated hsa-miR-889-3p expression was found to regulate cell cycle progression and influence osteosarcoma tumor size in vivo [18]. Fujiwara et al. [19] revealed that high serum concentrations of miR-25-3p were detected in osteosarcoma patients with lower OS. In a study by Roberto et al. [20], miR-138-5p was proposed as an intracellular mediator of invasion because upregulation of miR-138-5p is associated with reduced event-free survival and relapse.

Although the prediction of osteosarcoma prognosis by miRNAs has been extensively reported, as indicated above, there are few specific prognostic models, and most of them cannot be used clinically. Thus, identifying practical miRNA-based classifiers is of value and is promising for the diagnosis of osteosarcoma and predicting the prognosis of osteosarcoma patients. Here, we established and validated a reliable miRNA-based-OS predicting signature for osteosarcoma patients.

Material and Methods

Data acquisition

miRNA profile and correlated clinical records were extracted from the Gene Expression Omnibus (GEO) website (GSE39058). Then, the 91 osteosarcoma patients were randomly separated into a discovery group (64 samples) and a validation group (27 samples) for subsequently analyses.

Model establishment, Kaplan-Meier (K-M) analysis, and receiver operating characteristics (ROC) analyses

Univariate Cox regression and K-M analyses were performed to determine the OS-related miRNA candidates. miRNAs with P values less than 0.05 were identified as prognosis-related miRNAs (key miRNAs). Then, the formula was established referring to the hazard ratio (HR), and co-efficient (co-ef) derived from LASSO bagging Cox regression analysis. The osteosarcoma patients were grouped into low- and high-risk subsets, with the cut-off value set as the median risk score calculated by the formula. Then, the ROC curve was drawn to assess the stability of the model.

Stratified survival analysis

Multivariate Cox regression and K-M analyses were performed to determine the prognostic effect of the risk signature and various clinicopathological features, including age, sex, percent necrosis, and recurrence. In addition, to further explore the influence of our classifier on prognosis prediction in different subgroups, we performed subgroup analyses according to prognosis-related clinical features, including recurrence status (yes or no) and percent necrosis (<50% or >50%). The log-rank test was used to determine the significance of differences in survival curves, with a P value < 0.05 indicating significance.

Network construction and pathway enrichment analysis

Cytoscape was used to visualize the miRNA-mRNA network, while pathway enrichment analysis was performed to verify the biological functions of these targeted genes. R software was used in these statistical analyses (Version 3.4.0).

Results

Construction of the miRNA-based OS-predicting signature

miRNA array profiles and related clinical records for osteosarcoma patients were extracted from the GEO database. Then, 64 osteosarcoma patients were assigned to the training cohort, while 17 patients were assigned to the validation cohort. By univariate Cox analyses, we identified 39 miRNAs with the potential ability to predict the prognosis and survival statuses of osteosarcoma patients (P<0.05, Kaplan-Meier method, Supplementary Figure 1). Then, we rebuilt 1000 resample OS matrices based on the training set. RMIP of the selected miRNAs is displayed in Figure 1A and Table 1. A risk score generated from LASSO Cox regression analysis was estimated from the following formula: 4.93×miR-153–1.817×miR-591+2.38×miR-212 (Figure 1B, Table 2).

Figure 1

The establishment of the miRNA-based OS/RFS-classifier. (A) Distribution of 1000-times resampled results in OS-related miRNAs in the training cohort. (B) Hazard ratios of the enrolled OS-related miRNAs performed by LASSO Cox regression analysis.

Table 1

Clinical parameters of the training and validation datasets.

Parameters	Training set (%) N=64	Validation set (%) N=27	P value	SMD
Age (year)			0.678	0.175
≤18	53 (82.8)	24 (88.9)
>18	11 (17.2)	3 (11.1)
Sex			1.000	0.017
Male	29 (45.3)	12 (44.4)
Female	35 (54.7)	15 (55.6)
Recurrence			0.197	0.361
Yes	30 (46.9)	8 (29.6)
No	34 (53.1)	19 (70.4)
Death			0.864	0.101
Yes	17 (26.4)	6 (22.2)
No	47 (73.4)	21 (77.8)

SMD – Std mean difference

Table 2

Cox regression analysis was conducted to calculate the co-efficient of the OS-related miRNAs.

Gene_ID	Co-ef	Exp (co-ef)	Se (co-ef)	z	Pr (>|z|)
hsa-miR-153	4.938073	139.5012	1.237753	3.989546	6.62E-05
hsa-miR-591	−1.81787	0.162372	0.621404	−2.92542	0.00344
hsa-miR-212	2.383968	10.84787	0.575825	4.140094	3.47E-05

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

To test the discrimination ability of our miRNA-based OS-predicting signature, K-M analysis was used. The results revealed that most patients with high-risk scores have a significantly poorer OS compared to patients with low-risk scores (Figure 2A, 2C). We also used the ROC analysis to assess the stability of the OS signature, and the AUC value was 0.952 (95%CI: 0.891–1.00) in the training cohort and 0.793 (95%CI: 0.532–1.00) in the validation cohort (Figure 2B, 2D).

Figure 2

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

Comparison between miRNA-based predicting signature and clinicopathological features

Previous studies have shown age, sex, the percent of necrosis, and recurrence status had relatively significant impacts on prognosis. Therefore, we performed a multivariate Cox regression analysis to assess the relevant variables in predicting OS, including age, sex, percent necrosis, recurrence, and the 3-miRNA classifier (Table 3). The results identified recurrence factor and the miRNA-based signature with an HR of 12.28 as the most significant contributors to poor prognosis in the overall set (Figure 3A and Supplementary Figure 2). Compared to the recurrence factor, risk score had a similar prognostic accuracy as assessed by ROC curve (AUC=0.834 for recurrence status, AUC=0.800 for classifier) (Figure 3B). Clinicopathologic features (age, sex, and percent necrosis) failed to exhibit a consistently independent role for predicting OS outcomes in the overall set. In addition, we attempted to establish a nomogram integrating age, sex, percent necrosis, recurrence, and the 3-miRNA classifier to provide better overall estimation of the prognosis for osteosarcoma patients, with a predictive performance similar to that of recurrence status (Figure 3B).

Table 3

Cox regression analyses of OS-related miRNA signature and clinical features were used to evaluate the co-efficient.

Parameters	Co-ef	Exp (co-ef)	Se (co-ef)	z	Pr (>|z|)
Age>10	0.21771	1.243227	0.487284	0.446784	0.655031
Sex (male)	−0.15089	0.859944	0.465414	−0.3242	0.745785
Percent necrosis > 50%	−0.04133	0.95951	0.454333	−0.09097	0.927513
Recurrence (Y)	2.732481	15.37097	0.764471	3.57434	0.000351
Classifier (high-risk)	1.627931	5.093327	0.569527	2.858393	0.004258

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

Figure 3

The differences between the OS-related classifier and clinicopathological features. (A) Hazard ratios of the enrolled OS-related miRNAs established by LASSO Cox regression analysis. (B) ROC curve identified differences between the miRNA-based OS classifier and clinicopathological features in the overall cohort.

Functional enrichment and network visualization

To evaluate the potential biological processes and pathways that the identified miRNA might be involved in, GO, Reactome, hallmark, and KEGG enrichment analyses were conducted for the risky and protective miRNA. For the miRNA, the GO analyses showed target genes took part in cancer-related biological processes, such as the regulation of nervous system formation, Ras protein signal transduction pathway, and ubiquitin-protein ligase activity (Figure 4A–4C). The hallmark and KEGG analyses identified that these predicted genes were enriched in UV response of DNA, Hippo signaling, and Rap1 signaling pathways, which are directly or indirectly correlated with carcinogenesis (Figure 4D, 4E). Significant enrichment of cancer-related pathways was also observed in Reactome analysis (Figure 4F), such as the EGFR signaling pathway, PI3K/AKT pathway, and Receptor Tyrosine Kinases signaling pathway. The network of downstream genes for the 3 related miRNAs is displayed in Figure 5.

Figure 4

The pathway enrichment analyses for the miRNA targeted genes in the miRNA-based OS predicting signature. (A) GO-biological process analysis. (B) GO-cellular component analysis. (C) GO-molecular function enrichment. (D) Hallmark enrichment. (E) KEGG pathway enrichment. (F) Reactome analysis.

Figure 5

The miRNA and targeted gene interaction network.

Discussion

At present, the primary treatment approach for osteosarcoma is multimodality management. In the 1970s, Jaffe announced the first success of methotrexate in treating advanced disease, leading to the increased use of chemotherapy for osteosarcoma [21]. For localized disease, multidrug perioperative treatment combined with curative operation remains a standard treatment option. Many regimens are used for the treatment of osteosarcoma. MAP (high-dose methotrexate, anthracycline, and cisplatin) is used as a reference protocol [22,23]. Curative surgery includes amputation and limb-sparing techniques. Amputation surgery was used widely in the past. Nowadays, limb-sparing surgery is increasingly used. The role of radiotherapy is minimal in osteosarcoma; it can be used to treat residual disease in cases where limitedly resection surgery was performed due to the anatomical locations. In metastatic disease, chemotherapy remains the principal modality. Isolated lung metastasis can be resected in combination with chemotherapy. Chemotherapy has remained relatively unchanged and there has been little progress in the medical treatment of osteosarcoma. Resistance to drugs is troublesome. The above factors show the need to be able to predict the prognosis of the disease.

Growing evidence has revealed that miRNAs are involved in human diseases such as cancer, diabetes, hypertension, and cardiovascular disease [24]. Since a single miRNA can target several mRNAs, and being components of a larger regulatory network [25], miRNAs can affect cellular function, differentiation, proliferation, and apoptosis, as well as tumor formation, metastasis and drug resistance in carcinogenesis [26,27]. Here, we established and validated a 3-miRNA-based classifier to increase the accuracy of predicting the prognosis of patients with OS. The results showed that the classifier could precisely discriminate between the high-risk and low-risk sets, and the high-risk individuals were more likely to have worse OS than patients with a lower risk score. In addition, the signature was proved to serve as an independent prognostic factor for OS, and showed a higher predictive value than the available clinical records, except for recurrence status. The advantage of our signature over the recurrence status is that we obtained the risk score immediately after the surgery, much earlier than recurrence.

For the 3 miRNAs involved in our classifier, all 3 genes (miR-212, miR-591, and miR-153) in our model were confirmed to be downregulated in cancer tissues. Among them, the miR-212 gene has been reported to inhibit cell proliferation and tumor-promoting properties in gastric cancer [28], hepatocellular carcinoma [29], and lung cancer [30]. Luo et al. revealed that miR-212 was significantly downregulated in human osteosarcoma tissues, and ectopic miR-212 expression suppressed the proliferation and invasion of osteosarcoma cells [31], and this was subsequently validated by studies by Liu et al. [32] and Li et al. [33]. The miR-591 gene was initially identified as a tumor-suppressor gene due to its downregulation in breast cancer tissues, and low miR-591 expression was found to be involved in lymph node metastasis and advanced TNM stage in patients with breast cancer [34]. Further, Niu et al. [35] found that miR-153 was downregulated in osteosarcoma tissues compared with normal controls, potentially acting as a tumor suppressor through negatively regulating TGF-b2 expression, which agrees with the results of Wang et al. [36]. The above evidence is consistent with and confirms our findings of the association between the 3-gene model and OS of osteosarcoma patients. Moreover, we performed pathway enrichment analysis to provide new insight into the biological function of these 3 genes. The pathway analyses suggested that the predicted genes of these 3 miRNAs were mostly involved in tumor-related pathways, providing potential therapeutic targets for the treatment of osteosarcoma.

Conclusions

In conclusion, we constructed a 3-miRNA-based signature for predicting prognosis of patients with OS, which may improve clinical decision-making and advance the development of personalized medicine.

Kaplan-Meier curves of survival associated miRNA detected with uni-variable Cox regression analysis.

Kaplan-Meier curves for the high- and low-risk patients in different subgroups.