Development and validation of a novel prognostic signature for osteosarcoma

Introduction

Osteosarcoma (OS) is the most common primary malignant bone tumor, mainly affecting children and adolescents. This malignancy is distinguished by the synthesis of osteoid matrix by neoplastic cells and is frequently associated with aggressive clinical manifestations and a pronounced propensity for metastasis, especially to pulmonary sites. The prognosis for OS remains generally unfavorable, with survival outcomes heavily contingent upon a variety of clinical and molecular determinants.

The prognosis of OS is determined by a multitude of factors, including patient age, tumor size, anatomical location, and the presence of metastases at the time of diagnosis (1,2). Additionally, molecular and genetic determinants significantly impact OS outcomes. For example, elevated expression levels of specific genes, such as CCL5, have been found to be linked to enhanced metastatic potential and poorer prognostic outcomes (3). Immune cell infiltration within the tumor microenvironment significantly influences overall survival prognosis. Specifically, the presence and activity of tumor-associated macrophages have been correlated with tumor growth and metastatic behavior, thereby affecting patient prognosis (4). Recent research has underscored the prognostic relevance of genetic signatures and molecular subtypes in OS. The identification of immune-related glycosylation genes and their categorization into molecular subtypes have been demonstrated to be able to predict prognosis and inform therapeutic strategies, highlighting the significance of molecular profiling in the management of OS (5).

The discovery of cancer-related genes and clustering models increasingly relies on bioinformatics and data integration due to advancements in high-throughput sequencing technologies and data analysis methodologies. Recent studies have confirmed several biomarkers associated with OS progression and clinical diagnosis (6,7). The identification of tumor biomarkers and elucidation of gene functions are pivotal to numerous research efforts (8,9). However, well-validated and clinically applicable gene signatures for guiding therapeutic decisions in OS remain insufficient. Consequently, the identification of tumor markers and the development of feature gene models continue to be pivotal areas of investigation.

This study aims to develop and validate an innovative genetic prognostic signature for OS. Two OS patient datasets from the Gene Expression Omnibus (GEO) database were utilized to identify potential prognostic genes using survival analysis and intersection methods. These genes were utilized to classify OS patients into two distinct clusters using unsupervised clustering methods. Univariate Cox and the least absolute shrinkage and selection operator (LASSO) regression analyses identified ten significant genes. Ultimately, a novel gene signature was developed and validated. We present this article in accordance with the TRIPOD reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-2025-1491/rc).

Methods

Data acquisition

Gene expression data related to OS were sourced from the GEO in the National Center for Biotechnology Information (NCBI) database. The dataset GSE39058, which employed the GPL14951 platform and was deposited in 2013, includes 42 OS samples. The GSE21257 dataset, which was linked to the GPL10295 platform and deposited in 2012, contains 53 OS samples. The study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments.

Survival analysis and intersection analysis

Kaplan-Meier plots were generated using the R packages ’survival’ and ’survminer’. Additionally, Venn diagrams for analysis and mapping can be accessed via the website at http://bioinformatics.psb.ugent.be/webtools/Venn/.

Consensus clustering and unsupervised clustering

The R package “ConsensusClusterPlus” was employed to analyze gene expression patterns for the purpose of differentiating samples into distinct subtypes. Clustering results were validated using principal component analysis (PCA), t-distributed stochastic neighbor embedding (t-SNE), and uniform manifold approximation and projection (UMAP).

Differential expression analysis

Candidate genes with differential expression were identified using the ‘limma’ software, applying criteria of log₂(fold change) >1 and P value <0.05. Subsequently, hierarchical clustering plots of these genes were constructed using the “pheatmap” package.

Immune cell infiltration analysis

We performed single-sample Gene Set Enrichment Analysis (ssGSEA) on OS samples to clarify differences in immune cell composition across clusters. The characteristics of these groups were evaluated using Gene Set Enrichment Analysis (GSEA) and Gene Set Variation Analysis (GSVA) with integrated gene expression data.

Construction of gene signature

Genes linked to prognosis were identified through a univariate Cox proportional hazards analysis utilizing the ’survival’ package. A forest plot was created to display the significant genes. Key genes were further refined using LASSO regression analysis with the “glmnet” package. Prognostic genes were chosen to formulate a risk score, calculated as: Risk Score = β_gene1 × Expression_gene1 + β_gene2 × Expression_gene2 + ... + β_gene(n) × Expression gene(n), where β denotes the model-derived regression coefficient.

Validation and evaluation of the gene signature

In the GSE21257 cohort, OS patients were divided into high-risk and low-risk groups using the median risk score as a threshold. Kaplan-Meier survival curves and the log-rank test were used to analyze and compare the survival outcomes of these groups. Receiver operating characteristic (ROC) curve analysis was performed to assess the predictive accuracy of the gene risk scores, emphasizing sensitivity and specificity metrics. A P value below 0.05 was used to establish statistical significance. The gene signature was validated by applying the same risk formula to a different cohort from the GEO. A predictive nomogram for OS patients was created by combining risk assessment with clinical characteristics using the ‘rms’ package in R, and its predictive accuracy was evaluated through calibration plots. Decision curve analysis (DCA) was performed to assess the clinical utility of the 10-gene nomogram by measuring net benefits across different threshold probabilities.

Evaluation of the tumor microenvironment

This study utilized the CIBERSORT algorithm to examine immune cell infiltration differences between high- and low-risk groups using gene signatures. Subsequently, the Spearman method was applied to assess correlations among immune cells. Additionally, the ESTIMATE method was employed to quantify the levels of immune cells and immunologically relevant molecules in each GEO sample. Gene expression profiles of both groups were analyzed to calculate and evaluate the stromal, immune, and ESTIMATE scores, highlighting the variations.

Statistical analyses

Statistical analyses were performed using R statistical software (version 4.5.1). KM analysis was applied to estimate overall survival, while univariate analysis was carried out using Cox proportional hazards regression. All statistical P values were two-sided, and P<0.05 was considered statistically significant.

Results

Establishment of candidate genes and clusters

Survival analysis on samples from two GEO datasets identified 123 candidate genes through data intersection (Figure 1A). Subsequently, an unsupervised clustering technique was employed to classify the OS specimens within the GSE21257 cohort, facilitating a more comprehensive examination of the expression profiles of the candidate genes. The integration of these analyses indicated that the optimal number of groups for differentiation was two (k=2), as evidenced by the minimal overlap observed in the consensus matrix (Figure 1B). Consequently, the OS specimens were divided into two subgroups, with Figure 1C depicting the differential outcomes between them. Figure 1D-1F demonstrate the precise clustering results within the GSE21257 cohort, as validated by principal component analysis (PCA), UMAP, and t-SNE. A comparative analysis of the 123 candidate genes between clusters A and B in the GSE21257 dataset showed significant differences (P<0.05) (Figure 2A). The heatmap in Figure 2B illustrates the expression patterns of the 123 genes alongside clinical characteristics within each cluster. Moreover, the two groups displayed distinctly different immune characteristics. The ssGSEA algorithm revealed that 3 out of 23 immune cell types showed upregulation from cluster A to cluster B, as illustrated in Figure 2C. These findings suggest a potential correlation between these gene subtypes and immune status, which may inform treatment strategies for patients with OS. To explore the distinctions between the two clusters, we conducted GSVA and GSEA on each cluster. According to the GSVA results, cluster A showed increased activity in pathways like glycerolipid metabolism relative to cluster B (Figure 2D). GSEA analysis revealed that cluster A exhibited heightened activity in the intestinal immune network for IgA production, natural killer cell-mediated cytotoxicity, and the Toll-like receptor signaling pathway (Figure 2E).

Figure 1 Results of intersection analysis and unsupervised clustering analysis. (A) The identification of 123 candidate genes derived from survival analysis across two datasets. (B) The unsupervised clustering of these 123 candidate genes, highlighting the optimal consistency matrix with k=2 within the GSE21257 dataset. (C) Survival rates corresponding to the two subtypes identified in GSE21257. (D-F) The differentiation between the two subtypes using PCA (D), t-SNE (E), and UMAP (F). PCA, principal component analysis; t-SNE, t-distributed stochastic neighbor embedding; UMAP, uniform manifold approximation and projection.

Figure 2 Subtypes of OS associated with 123 candidate genes identified in the GSE21257 dataset. (A) Differential expression of these 123 genes across the two identified clusters. (B) Heat map depicting the expression of the candidate genes, categorized by clinical information and the two subtypes. (C) Distinctive expression patterns of immune cells corresponding to each subtype. (D) Results of GSVA for the two subtypes, while (E) the GSEA specific to subtype A. *, P<0.05; **, P<0.01; ***, P<0.001. GSEA, Gene Set Enrichment Analysis; GSVA, Gene Set Variation Analysis; OS, osteosarcoma.

Development and validation of the gene signature

A univariate Cox regression analysis was conducted to assess the prognostic significance of 123 candidate genes in individuals with OS. Figure 3A illustrates a significant correlation between 46 genes and overall survival. Our prognostic model for OS was developed using LASSO regression analysis, as depicted in Figure 3B,3C. The prognostic risk score was determined using the formula: (3.4065 × FAM45B) + (1.5741 × KIR2DS1) + (2.4418 × KRTAP1-3) + (0.2860 × MCAM) − (3.4604 × MEST) + (3.1166 × MRPL4) + (5.8418 × PDILT) − (0.9179 × TGFB2) + (0.3215 × TMEM88) − (0.6040 × UBE2L3). Patients in the GSE21257 dataset were categorized into low-risk and high-risk groups according to the median risk score. Survival rate analysis showed that high-risk patients had significantly shorter overall survival than low-risk patients (P<0.05; see Figure 3D). The model’s prognostic accuracy was assessed using the ROC curve. Figure 3E indicates that the gene signature is a reliable predictor of overall survival in OS patients, with a P value of 0.006 in the GSE39058 dataset. Figure 3F shows the 3-year area under the curve (AUC) as 0.945. Figure 4A-4F illustrate the analysis of Z-score transformed patient risk scores and survival status in the GSE21257 dataset and GSE39058 dataset. The findings indicate that OS patients in the high-risk category exhibit a notably higher mortality rate than those in the low-risk category. To evaluate the robustness of this model, the same coefficients were applied to the validation cohort. Patients from the GSE39058 dataset were classified into low-risk and high-risk categories using a standardized risk formula, which demonstrated significantly different overall survival outcomes.

Figure 3 The results of regression analysis, survival analysis, and ROC curve assessment. (A) The forest plot derived from univariate Cox regression analysis identifies 46 key genes significantly associated with overall survival. (B,C) LASSO regression results demonstrate that all 10 key genes are crucial for the model’s construction. (D,E) Kaplan-Meier survival curves for the training set (GSE21257 dataset) and the validation set (GSE39058) depict overall survival stratified by low and high risk groups, as determined by the median risk score. (F) The ROC curves for 1-, 3-, and 5-year intervals are shown for the training groups based on the gene signature. AUC, area under the curve; LASSO, least absolute shrinkage and selection operator; ROC, receiver operating characteristic curve.

Figure 4 The risk scores associated with the gene signature across two GEO cohorts. Heatmaps illustrating the expression levels of the 10 key genes within the GSE21257 cohort (A) and the GSE39058 cohort (B). (C,D) The distribution of patients stratified by varying risk scores within these two GEO cohorts. (E,F) The survival status of patients categorized by different risk scores in each of the GEO cohorts. GEO, Gene Expression Omnibus.

Creating and assessing a nomogram for OS

To derive a prognostic signature with potential clinical applicability, a nomogram was developed (Figure 5A). The calibration curves for the nomogram-projected overall survival demonstrated satisfactory concordance between the observed and predicted survival probabilities (Figure 5B). The DCA demonstrated that the nomogram, derived from the 10-gene signature model, improved clinical utility for predicting overall survival prognosis (Figure 5C).

Figure 5 Prognostic nomogram and the association between the gene signature and immune status in patients with OS. (A) Nomogram plot; (B) calibration plot for the nomogram; (C) DCA for the nomogram and risk score concerning 3-year overall survival in OS; (D) immune cell expression levels between high- and low-risk groups. (E) violin plots representing immune cell expression; (F) correlation analysis of immune cells. *, P<0.05; ***, P<0.001. DCA, decision curve analysis; OS, osteosarcoma.

Evaluation of the immune microenvironment

CIBERSORTx was employed to assess the presence of 22 immune cell infiltrates within the OS tissue of patients from the GSE21257 dataset. Immune cell infiltration patterns were identified by ranking OS patients according to their risk scores, from highest to lowest (Figure 5D,5E). The analysis indicated that the low-risk group exhibited notably higher levels of CD8⁺ T cells and regulatory T cells (Tregs), along with lower levels of resting memory CD4⁺ T cells and M0 macrophages, in comparison to the high-risk group. Spearman’s correlation analysis demonstrated a robust inverse relationship, as well as a positive correlation, between gene signature risk scores and the levels of memory B cells, CD8⁺ T cells, Tregs, M0 macrophages, and resting memory CD4⁺ T cells (Figure 5F, Figure 6A-6E). A heatmap illustrated the association between 22 immune cells and a 10-gene signature (Figure 6F). We also analyzed the immune microenvironment variations between high-risk and low-risk groups in the OS dataset. Low-risk patients demonstrated elevated stromal, immune, and ESTIMATE scores relative to high-risk patients (Figure 6G).

Figure 6 Analysis of correlation and immune microenvironment scores. (A-E) The correlation analysis between risk scores and the abundance of five distinct immune cell types. (F) The correlation analysis between immune cells and the gene signature. (G) The variation in immune microenvironment scores between high-risk and low-risk groups. *, P<0.05; **, P<0.01; ***, P<0.001. TME, tumor microenvironment.

Discussion

The stagnation in clinical outcomes for OS patients represents a critical issue, underscored by the enduring absence of innovative therapies and biomarkers. OS, as the most prevalent primary malignant bone tumor, has experienced minimal progress in therapeutic advancements over the past several decades. This lack of progress persists despite the disease’s aggressive nature and the associated high mortality rate. The prevailing standard treatment, which involves surgical resection in conjunction with chemotherapy, has remained largely unchanged since the 1970s. Consequently, survival rates have plateaued, particularly among patients with metastatic or recurrent disease, where prognosis continues to be poor (10,11). Therefore, the identification of OS tumor markers and the development of feature gene models remain central to ongoing research efforts.

The study utilized two OS datasets from the GEO database to identify 123 candidate genes. An unsupervised clustering algorithm categorized the OS samples into two clusters, Cluster A and Cluster B, based on candidate gene expression profiles. The immune infiltration analysis revealed a significantly higher prevalence of diverse immune cell types in Cluster A than in Cluster B, indicating a meaningful correlation between gene subtypes and immune status. These findings highlight the potential of these subtypes to classify OS patients into distinct tumor categories, offering a new framework for molecular classification and immunotherapy strategies in OS. Using univariate Cox regression and LASSO analysis, ten genes significantly associated with patient prognosis (P<0.05) were identified. A predictive signature was created and validated for effectively differentiating between high-risk and low-risk OS patient groups. Significant differences in overall survival were observed among these risk categories, with comparable survival discrepancies noted in validation cohorts. ROC analysis demonstrated that the gene signature effectively predicted the prognosis of patients with OS across the specified cohorts. CIBERSORTx analysis indicated that the low-risk group had significantly elevated levels of CD8⁺ T cells and Tregs, and decreased levels of resting memory CD4⁺ T cells and M0 macrophages, in comparison to the high-risk group. Additionally, a significant inverse relationship was found between gene signature risk scores and the abundance of memory B cells, CD8⁺ T cells, and Tregs. A strong positive correlation was found between the risk scores and the levels of M0 macrophages and resting memory CD4⁺ T cells. The analysis of stromal, immune, and ESTIMATE scores revealed notable differences in the immune microenvironment between high- and low-risk groups, underscoring the unique features of the tumor microenvironment based on patient risk stratification.

A total of ten genes were identified within the gene signature. Recent investigations have underscored the significance of melanoma cell adhesion molecule (MCAM) in the progression and metastasis of OS. Evidence indicates that MCAM is upregulated in OS patients exhibiting metastases, thereby suggesting its potential utility as a prognostic biomarker and therapeutic target (12). The increased expression of MCAM in tip-like endothelial cells suggests its potential as a therapeutic target to reduce OS metastasis (13). Transforming growth factor beta 2 (TGFB2), part of the TGF-β family, is involved in multiple aspects of OS progression. TGFB2 is crucial in the epithelial-mesenchymal transition (EMT), a process that increases cancer cells’ metastatic potential. In the context of OS, TGFB2 expression correlates with increased invasiveness and metastatic capability, as it facilitates the transition of cells to a more mesenchymal phenotype, which is essential for metastasis (14,15). Furthermore, recent studies have elucidated the role of ubiquitin-conjugating enzyme E2 L3 (UBE2L3) in the progression of this malignancy. UBE2L3 is a crucial regulator of oxidative stress and necroptosis, both vital in OS (16,17).

Recent studies highlight the crucial role of CD8⁺ T cells in OS prognosis and treatment, emphasizing their potential as therapeutic targets. One study examined the influence of specific genes on CD8⁺ T cells within the context of OS, identifying ecotropic viral integration site 2B (EVI2B) as a protective immune-related gene. The elevated expression of EVI2B was associated with an increased proportion of CD8⁺ T cells, which in turn promoted the expression of granzyme A and K, thereby enhancing the cytotoxic effects on tumor cells (18). Another study concentrated on the prognostic implications of CD8⁺ T cell-associated gene signatures in OS. Researchers developed a risk score model based on CD8⁺ T cell gene markers by analyzing mRNA expression data. The model demonstrated a strong correlation with patient prognosis and immunotherapy response, highlighting the crucial role of CD8⁺ T cells in OS treatment (19). Moreover, interleukin-35 (IL-35), an immunosuppressive cytokine, has been shown to suppress the antitumor function of CD8⁺ T cells in OS patients (20). Furthermore, Tregs are recognized for their ability to suppress the activity of effector T cells, thereby promoting tumor immune evasion and contributing to unfavorable prognoses in OS patients. Research has indicated that the interaction between osteoclasts and regulatory CD4⁺ T cells can substantially modify the tumor microenvironment, resulting in varied survival outcomes for individuals with OS (21). Additionally, studies have revealed that M0 macrophages, along with M2 macrophages, are among the predominant immune cells infiltrating OS tissues, and their presence is associated with poor survival outcomes (22).

This research is subject to several limitations. Primarily, the gene signature was developed using a dataset with a limited sample size, which may result in potential inaccuracies. Furthermore, the gene signature was exclusively derived through bioinformatics approaches. Therefore, it is advisable for subsequent studies to utilize larger datasets and to incorporate experimental, in vivo validation of the findings. The prognostic value of the new gene signature requires further validation in OS patients.

Conclusions

In conclusion, we have discovered a ten-gene signature that serves as a new prognostic predictor for OS patients. This study offers new insights and introduces novel molecular biomarkers pertinent to OS prognosis. Furthermore, the findings have the potential to facilitate the discovery of new therapeutic targets with clinical applicability.

Acknowledgments

The authors are grateful to all study participants.

Footnote

Reporting Checklist: The authors have completed the TRIPOD reporting checklist. Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2025-1491/rc

Peer Review File: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2025-1491/prf

Funding: None.

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tcr.amegroups.com/article/view/10.21037/tcr-2025-1491/coif). The authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. The study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments.

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.

References

1. Hagleitner MM, Hoogerbrugge PM, van der Graaf WT, et al. Age as prognostic factor in patients with osteosarcoma. Bone 2011;49:1173-7. [Crossref] [PubMed]
2. Fu P, Shi Y, Chen G, et al. Prognostic Factors in Patients With Osteosarcoma With the Surveillance, Epidemiology, and End Results Database. Technol Cancer Res Treat 2020;19:1533033820947701. [Crossref] [PubMed]
3. Sun K, Gong C, Peng H, et al. High CCL5 expression is associated with osteosarcoma metastasis and poor prognosis of patients with osteosarcoma. Mol Med Rep 2017;16:6953-7. [Crossref] [PubMed]
4. Xiao KW, Liu ZB, Zeng ZH, et al. Construction and Validation of a Macrophage-Associated Risk Model for Predicting the Prognosis of Osteosarcoma. J Oncol 2021;2021:9967954. [Crossref] [PubMed]
5. Wang W, Jiao Y, Du X, et al. Immune-related glycosylation genes based classification predicts prognosis and therapy options of osteosarcoma. Gene 2025;933:148985. [Crossref] [PubMed]
6. Sun W, Ma X, Shen J, et al. Bioinformatics analysis of differentially expressed pathways related to the metastatic characteristics of osteosarcoma. Int J Mol Med 2016;38:466-74. [Crossref] [PubMed]
7. Zhao L, Zhang J, Tan H, et al. Gene function analysis in osteosarcoma based on microarray gene expression profiling. Int J Clin Exp Med 2015;8:10401-10. [PubMed]
8. Wang DW, Yu SY, Cao Y, et al. Identification of CD20, ECM, and ITGA as Biomarkers for Osteosarcoma by Integrating Transcriptome Analysis. Med Sci Monit 2016;22:2075-85. [Crossref] [PubMed]
9. Xiong Y, Wu S, Du Q, et al. Integrated analysis of gene expression and genomic aberration data in osteosarcoma (OS). Cancer Gene Ther 2015;22:524-9. [Crossref] [PubMed]
10. Zhao X, Wu Q, Gong X, et al. Osteosarcoma: a review of current and future therapeutic approaches. Biomed Eng Online 2021;20:24. [Crossref] [PubMed]
11. Harrison DJ, Geller DS, Gill JD, et al. Current and future therapeutic approaches for osteosarcoma. Expert Rev Anticancer Ther 2018;18:39-50. [Crossref] [PubMed]
12. Du X, Zhang Q, Wang S, et al. MCAM is associated with metastasis and poor prognosis in osteosarcoma by modulating tumor cell migration. J Clin Lab Anal 2022;36:e24214. [Crossref] [PubMed]
13. Tang H, Cai Y, Yang M, et al. Single-cell and spatial transcriptomics reveals the key role of MCAM(+) tip-like endothelial cells in osteosarcoma metastasis. NPJ Precis Oncol 2025;9:104. [Crossref] [PubMed]
14. Ruh M, Stemmler MP, Frisch I, et al. The EMT transcription factor ZEB1 blocks osteoblastic differentiation in bone development and osteosarcoma. J Pathol 2021;254:199-211. [Crossref] [PubMed]
15. Trivedi T, Pagnotti GM, Guise TA, et al. The Role of TGF-β in Bone Metastases. Biomolecules 2021;11:1643. [Crossref] [PubMed]
16. Wei N, Chao-Yang G, Wen-Ming Z, et al. A ubiquitin-related gene signature for predicting prognosis and constructing molecular subtypes in osteosarcoma. Front Pharmacol 2022;13:904448. [Crossref] [PubMed]
17. Zhao X, Shan G, Xing D, et al. UBE2L3 Suppresses Oxidative Stress-regulated Necroptosis to Accelerate Osteosarcoma Progression. Recent Pat Anticancer Drug Discov 2025;20:102-12. [Crossref] [PubMed]
18. Xie T, Feng W, He M, et al. Analysis of scRNA-seq and bulk RNA-seq demonstrates the effects of EVI2B or CD361 on CD8(+) T cells in osteosarcoma. Exp Biol Med (Maywood) 2023;248:130-45. [Crossref] [PubMed]
19. Lu Y, Cao N, Zhao M, et al. Importance of CD8 Tex cell-associated gene signatures in the prognosis and immunology of osteosarcoma. Sci Rep 2024;14:9769. [Crossref] [PubMed]
20. Liu MX, Liu QY, Liu Y, et al. Interleukin-35 suppresses antitumor activity of circulating CD8(+) T cells in osteosarcoma patients. Connect Tissue Res 2019;60:367-75. [Crossref] [PubMed]
21. Li F, Tang H, Luo X, et al. Interaction gene set between osteoclasts and regulatory CD4(+) T cells can accurately predict the prognosis of patients with osteosarcoma. Cancer Sci 2023;114:3014-26. [Crossref] [PubMed]
22. Cao H, Quan S, Zhang L, et al. BMPR2 expression level is correlated with low immune infiltration and predicts metastasis and poor survival in osteosarcoma. Oncol Lett 2021;21:391. [Crossref] [PubMed]