Exploring osteosarcoma based on the tumor microenvironment

Abstract

Osteosarcoma is a cancerous bone tumor that develops from mesenchymal cells and is characterized by early metastasis, easy drug resistance, high disability, and mortality. Immunological characteristics of the tumor microenvironment (TME) have attracted attention for the prognosis and treatment of osteosarcoma, and there is a need to explore a signature with high sensitivity for prognosis. In the present study, a total of 84 samples of osteosarcoma were acquired from the UCSC Xena database, analyzed for immune infiltration and classified into two categories depending on their immune properties, and then screened for DEGs between the two groups and analyzed for enrichment, with the majority of DEGs enriched in the immune domain. To further analyze their immune characteristics, the immune-related genes were obtained from the TIMER database. We performed an intersection analysis to identify immune-related differentially expressed genes (IR-DEGs), which were analyzed using a univariate COX regression, and LASSO analysis was used to obtain the ideal genes to construct the risk model, and to uncover the prognostic distinctions between high-risk scoring group and low-risk scoring group, a survival analysis was conducted. The risk assessment model developed in this study revealed a notable variation in survival analysis outcomes between the high-risk and low-risk scoring groups, and the conclusions reached by the model are consistent with the findings of previous scholars. They also yield meaningful results when analyzing immune checkpoints. The risk assessment model developed in this study is precise and dependable for forecasting outcomes and analyzing characteristics of osteosarcoma.

Keywords: immune-related genes; immunization checkpoints; immunotherapy; osteosarcoma; tumor microenvironment.

Figure 1

The ESTIMATE algorithm was utilized to conduct survival analysis on patient groups categorized based on high and low immune scores, and the results were visually represented through Kaplan-Meier (K-M) curves. (A) ESTIMATE score, (B) Tumor purity score, (C) Stromal score, and (D) Immune score.

Figure 2

Immune subtype identification and comparative analysis. The symbols * represent p-values less than 0.05, ** represent p-values less than 0.01, *** represent p-values less than 0.001, **** represent p-values less than 0.0001. (A) ssGSEA analyses of 84 samples, divided into two groups based on 30 ssGSEA scores per sample. (B) Unsupervised clustering of the samples based on their immunological characteristics, where the number of clusters with the highest average within-group agreement is 2. (C–F) shows, in order, the differences in ESTIMATEScore, ImmuneScore, StromalScore, and TumorPurity between the high and low immunity groups. (G) A box plot is utilized to display the levels of immune cell infiltration in groups categorized as either having high or low immunity. In this visualization, the red boxes correspond to the high immunity group, while the blue boxes correspond to the low immunity group. (H) PCA analysis of the two immune subtypes, with purple and yellow dots representing the immunity–high and —low groups, respectively. "ns" stands for no statistical difference.

Figure 3

Mutations between high and low immunity groups. The symbols represent p-values less than 0.01. (A) Mutation status of genes in the high and low immunity groups. (B) TMB distribution of all samples. (C) Bar graph showing the difference in TMB between the high and low immunity groups. ** represent p-values less than 0.01.

Figure 4

Enrichment analysis was conducted on the DEGs identified in the two immune subtypes. (A) Volcano diagram showing the regulation of DEG expression, with green, grey, and red dots representing down-regulation, unregulation, and up-regulation, respectively. (B) Bubble diagram showing the top ten pathways according to Gene Set Enrichment Analysis (GSEA). (C) Bubble plots showing the top 10 enriched GO BP, CC, and MF. (D) Bubble plots showing the top 20 enriched terms of the KEGG pathway, with the size of the dots representing the number of enrichments.

Figure 5

Identification and enrichment analysis of differentially expressed genes (DEGs) associated with the immune system. (A) Venn diagram showing 221 immune-associated DEGs overlapping 836 DEGs and 1811 IRGs. (B) Protein interaction network diagram of DEGs (C) Bubble plots showing the top 10 enriched GO BP, CC, and MF. (D) Bubble plots showing the top 20 enriched terms of the KEGG pathway, with the size of the dots representing the number of enrichments. (E) Bubble plots showing the top 20 enriched terms of the KEGG pathway, based on GSEA analyses of the top 10 pathways with the highest gene enrichment.

Figure 6

Analysis of immune cell infiltration was conducted on the ten hub genes, visualized using lollipop charts.From left to right, the order is CCR5, TNF, IL10, IL6, CD8A, CD4, IL1B, CCR7, CCL5, and CCL2.

Figure 7

Ten hub genes were subjected to survival analysis.

Figure 8

Construction of the total risk profile. The symbols * represent p-values less than 0.05, *** represent p-values less than 0.001. (A, B) Employing the LASSO method for the identification of key candidate genes. (C–E) Survival analysis K-M curves for the training cohort, validation cohort, and initial combined cohort, respectively. (F–H) ROC curves were generated to assess the prognostic significance of risk features in the training, validation, and original merged groups. (I) Differences in PDK1 and PPARG gene expression between high and low-risk score groups. (J) Differences in TMB between high and low-risk score groups. (K) Differences in TIDEscore in high and low-risk score groups. (L) Difference in risk scores between high and low TMB groups.

Figure 9

Role of risk models for immune checkpoints. The symbols * represent p-values less than 0.05, ** represent p-values less than 0.01, and *** represent p-values less than 0.001. (A) Comparison of the expression levels of immune checkpoint-related genes in the group with high-risk scores and the group with low-risk scores. (B, C) The expression levels of immune checkpoint-related genes between the two groups by dividing the gene expression of PDK1 and PPARG into high and low groups according to the median. (D) CIBERSORT immune infiltration analysis of all samples. (E, F) Infiltration abundance levels of some immune cells by dividing gene expression of PDK1 and PPARG into high and low groups by median. (G) MicroenvironmentScore in the high and low-risk scoring groups.