Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 14-Day Trial for You or Your Team.

Learn More →

Unveiling and validation of a disulfidptosis determined prognostic model for osteosarcoma: new insights from prognosis to immunotherapy and chemotherapy

Unveiling and validation of a disulfidptosis determined prognostic model for osteosarcoma: new... IntroductionOsteosarcoma represents one of the most prevalent primary bone malignancies, exhibiting elevated metastasis and recurrence rates, and predominantly afflicts children, adolescents, and young adults [1], [2], [3]. In the United States alone, the incidence of osteosarcoma surpasses 1,000 novel cases, annually [4]. The 5-year event-free survival rate for patients afflicted with nonmetastatic osteosarcoma is approximately 60 %, whereas for patients presenting with a primary tumor site, multiple pulmonary nodules, or other metastatic lesions, this falls below 20 % [5, 6]. Modifications in TP53 and RB1 tumor suppressor genes, the bone microenvironment, and osteosarcoma predisposition have been posited to be linked to the tumorigenesis and metastasis of osteosarcoma [7]. The optimal, standard-of-care therapeutic intervention for osteosarcoma employs high-dose methotrexate, doxorubicin, and cisplatin (MAP) [3, 7]. Additional treatments, including surgical intervention and radiotherapy, have been shown to provide localized modality control for metastatic sites [6]. Nevertheless, the scarcity of osteosarcoma has resulted in a paucity of dependable markers, thereby impeding the establishment of systematic population-based osteosarcoma screening programs as a component of preventative strategies [6]. Consequently, the discernment of pivotal pathways and genes assumes a crucial role in osteosarcoma diagnosis and the exploration of innovative therapeutics.Disulfides denote organic sulfur compounds encompassing disulfide linkages, predominantly occurring in cystine within living organisms [8]. Elevated intracellular disulfide formation can augment the levels of reactive oxygen species (ROS) within the mitochondria of cancer cells, thereby potentially inducing mitochondrial dysfunction and subsequently activating cellular demise [9]. It has previously been proposed that the overexpression of solute carrier family 7 member 11 (SLC7A11) contributes to cellular apoptosis. SLC7A11 assists cancer cells in the uptake of cysteine, thereby sustaining intracellular glutathione levels and safeguarding cells from oxidative stress-induced cell death [10]. SLC7A11 is commonly overexpressed in cancer and induces cell death by converting cystine to the more soluble cysteine, which necessitates a substantial amount of NADPH and renders cancer cells to be reliant on the pentose phosphate pathway (PPP) [10], [11], [12]. Restricting glucose availability to tumor cells augments cystine accumulation and exacerbates the depletion of NADPH and GSH [10, 13]. Pharmaceutical combinations that simultaneously target GLUT1 and GSH synthesis may induce tumor cell death, thus offering a novel potential therapeutic strategy for osteosarcoma treatment [10, 13].Recently, Liu et al. reported disulfidptosis, which is different from apoptosis and ferroptosis [14]. In this study, genes were procured that encode proteins associated with augmented disulfide bonds under conditions of NADPH depletion. Disulfidptosis-associated genes were applied to facilitate the subtyping of osteosarcoma, yielding three subgroups with a distinct prognosis. Subsequently, differentially-expressed genes were identified between best survival and worst survival groups, reflecting tumor alterations post-disulfidptosis. Subsequently, genes indicative of osteosarcoma patient prognosis were screened and prognostic models were established.MethodsBaseline information of cohortsIn this study, two osteosarcoma assemblages participated. For the Therapeutically Applicable Research to Generate Effective Treatments Osteosarcoma (TARGET-OS) cohort, the gene expression profiles of 84 patients were acquired from the UCSC Xena platform (http://xena.ucsc.edu/), which were accompanied by the patient’s clinical information. Furthermore, the GSE21257 collection, encompassing 53 patients, was obtained from the Gene Expression Omnibus platform (GEO, http://www.ncbi.nlm.nih.gov/geo/). The baseline information for patients in both cohorts is presented in Table 1.Table 1:Population demographics and prognostic outcomes of enrolled patients.GSE21257 (n=53)TARGET-OS (n=84)p-ValueAge, years Mean (SD)18.7 (12.2)15.0 (4.83)0.0379a Median [Min, Max]16.7 [3.08, 79.0]14.4 [3.56, 32.4]Gender Female19 (35.8 %)37 (44.0 %)0.376 Male34 (64.2 %)47 (56.0 %)Tumor site Arm/hand8 (15.1 %)6 (7.1 %)0.145 Leg/foot44 (83.0 %)76 (90.5 %) Pelvis0 (0 %)2 (2.4 %) Unknown1 (1.9 %)0 (0 %)Metastatic at diagnosis No39 (73.6 %)63 (75.0 %)0.844 Yes14 (26.4 %)21 (25.0 %)Alive status Alive30 (56.6 %)57 (67.9 %)0.205 Dead23 (43.4 %)27 (32.1 %)Overall survival time, months Mean (SD)68.5 (59.3)48.9 (34.0)0.0316a Median [Min, Max]45.0 [4.00, 246]47.9 [2.43, 143]ap-Value<0.05.Collection of disulfidptosis-associated genesLiu et al. reported a novel type of disulfidptosis. Cells exhibiting elevated SLC7A11 levels after being subjected to glucose scarcity showed augmented cystine uptake. In line with an inadequate NADPH provision, this culminates in NADPH depletion, aberrant disulfide bonding within actin cytoskeleton proteins, actin lattice disintegration, and consequently cellular demise [14]. Ninety cysteine sites were documented with a 1.5-fold augmentation in disulfide bonds subsequent to glucose deprivation, encoded by 77 genes, the ID of the genes were collected from Supplementary Table 2 of Liu et al.’s study.Identifying distinct disulfidptosis subtypesThe consensus clustering and molecular subtype discernment were conducted with the R package named “ConsensusClusterPlus” [15]. Fifty iterations were implemented to reveal the cluster k-means, of which at least 80 % of samples were enrolled each time. The score of cumulative distribution function (CDF), which means the area that under the CDF curve, was calculated for the selection of the optimal cluster number. The principal component analysis (PCA) was harnessed to substantiate the credibility of the consensus clusters, as well as the t-distributed stochastic neighbor embedding (t-SNE).Gene set enrichment analysisDifferentially-expressed genes (DEGs) are important for identification of different mechanisms among subgroups and were collected by the ‘limma’ package in R. An absolute fold-change exceeding 1.5 and an adjusted p-value below 0.01 were established as threshold values for DEGs filtration. Packages ‘org.Hs.e.g.db’ and ‘msigdbr’ were utilized to conduct Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG), and HALLMARK enrichment analyses [16, 17], while ‘ClusterProfiler’ was employed to explicate representative signaling pathways [18].Constructing the prognostic signatureUnivariate Cox examination was executed to investigate prognostic genes in both TARGET-OS and GSE21257 cohorts. Subsequently, a Venn diagram was amalgamated, depicting concordant prognostic genes for ensuing analysis. The least absolute shrinkage and selection operator (LASSO) regression analysis, executed by the ‘glmnet’ package in R [19], was harnessed to refine and regularize all input factors, ultimately yielding a statistical model with exceptional predictive precision and interpretability. LASSO, a regularization technique, can avert overfitting, conduct feature selection more proficiently, and optimize the structure of the regression model more effectively compared to traditional ordinary least squares subset selection methods. It is particularly suitable for situations involving numerous potential predictors where only a few are presumed pertinent. Genes selected by LASSO analysis were employed to calculate the risk score of each patient, thereby utilizing gene expression values and coefficients to devise the following disulfidptosis-related (DSPR) prognostic signature:DSPR=∑i=1nCi×Ei$$\mathrm{D}\mathrm{S}\mathrm{P}\mathrm{R}=\sum\limits _{i=1}^{n}{C}_{i}\times {E}_{i}$$where Ci${C}_{i}$denotes the coefficient of a gene, and Ei${E}_{i}$represents the expression value of the corresponding gene.Establishing a predictive nomogramIn this study, we embarked on a comprehensive multivariate Cox regression examination as the initial endeavor, with the objective of pinpointing autonomous prognostic determinants. Subsequently, the DSPR signature was integrated into the analytical framework. Next, a nomogram was created showing the autonomous prognostic constituents to amalgamate the prognosticated values. The “rms” and “regplot” R packages were employed for this endeavor. To ascertain the accuracy and reliability of our study, a calibration diagram, decision curve scrutiny (DCA), and clinical impact curve assessment were executed.Personalized chemotherapy predictionBy leveraging the pharmacological susceptibility and phenotypic data from GDSC 2016 (https://www.cancerrxgene.org/), the R package “pRRophetic” was employed to predict the chemotherapeutic responsiveness for each patient [20]. Estimated IC50 values (lower IC50 indicates enhanced therapy susceptibility) for each specimen exposed to a specific chemotherapeutic agent were obtained via ridge regression [21], and the prediction accuracy was ascertained through a tenfold cross-validated process [22].StatisticsThe T-test or Mann-Whitney U test was employed to evaluate disparities between subgroups with continuous variables. Survival curves were created utilizing Kaplan-Meier survival analysis and log-rank analysis, while the hazard ratio (HR) and 95 % confidence interval (CI) for candidate genes were ascertained via a univariate Cox regression model. Multivariate Cox analyses were employed to determine the independent prognostic influence of the risk score, subsequent to adjustments for assorted clinical characteristics. Adopting p<0.05 as the criterion for significance, two-tailed statistical evaluation was conducted. R software (4.1.2) was employed for all analytical analyses.ResultsUnveiling the disulfidptosis-defined subcategoriesAs mentioned in the methods section, 79 disulfidptosis associated genes were collected, among which 40 could be identified in both the TARGET-OS cohort and GSE21257 cohort. Therefore, the gene expression profiles of the 40 genes were collected as the input matrix for clustering. The consensus matrices for k=2, 3, 4, 5, 6 were tested (Figure 1A). With the cumulative distribution function (CDF) curve, it was observed that while the area under the CDF did not change significantly with the number of k after k=3 (Figure 1B), and that the relative change of the CDF between k=3 and k=4 showed the biggest alteration. Therefore, it was recognized that k=3 was the best cluster number for subsequent analyses (Figure 1C).Figure 1:Identification the optimal number of clusters by ConsensusClusterPlus. (A) Heatmaps of the consensus matrices for k=2, 3, 4, 5, 6; (B) the cumulative distribution functions of the consensus matrix for each k; (C) the relative change in area under the CDF curve comparing k and k − 1.In this study, the overall survival (OS) prognosis of patients with distinct subtypes was evaluated. Subtype C1 comprised 29 individuals, among whom 15 patients succumbed, while subtype C2 included 41 individuals with 12 fatalities. Subtype C3 encompassed 15 individuals, all of whom survived until the conclusion of the follow-up period. Consequently, among the three subtypes delineated by disulfidptosis-associated genes, subtype C3 exhibited the most favorable prognosis, while subtype C1 demonstrated the least favorable outcome, with a marked disparity in OS prognosis among the three patient cohorts (overall: p=0.0048; C1 vs. C2: p=0.0550; C1 vs. C3: p=0.0072; C2 vs. C3: p=0.0550 Figure 2A). Patients within the three subtypes were both mutually independent and intersected among other subtypes (Figure 2B). Furthermore, an increased prevalence of tumors situated in the arm/hand region among subtype C2 patients was observed, while all tumors in subtype C3 were located in the leg/foot. Regarding disulfidptosis-associated genes, it was observed that subtype C2 exhibited elevated levels of FLNC and MYH14. High expression of MYH4, MYH6, MYH8, MYH13, CHCHD3, and PRDX1 in subtype C1 could potentially contribute to an unfavorable prognosis. In the context of subtype C3, heightened expression of TLN1, ZHX2, MYH9, and MYH10 was observed (Figure 2C). In addition, immune infiltration among three subtypes was evaluated, and it was observed that the best prognosis C3 subtype presented the highest infiltration of immunocytes, especially for CD56 bright natural killer cells, central memory CD8 T cells, eosinophil,s gamma delta T cells, macrophages, and myeloid-derived suppressor. Cells (MDSC). Regarding the C1 subtype with the worst prognosis, immunocyte infiltration was at an intermediate level (Figure 2D).Figure 2:Revealing three disulfidptosis determined molecular subtypes. (A) Diverse overall survival outcome of three subtypes; (B) separation and intersection of the distribution of three subtypes; (C) expression levels of 40 disulfidptosis genes across different patients; (D) differential immune cell infiltration among three subtypes.Determine the disulfidptosis-influenced pivotal genes and pathwaysBecause of apparent differences in OS prognosis between subtypes C1 and C3, DEGs within these two subcategories were compared (Figure 3A). Employing a screening threshold of an absolute log2 (fold change) greater than 1.5 and a p-value less than 0.05, 56 DEGs were identified; the expression patterns of which in distinct patients are depicted in Figure 3B. Furthermore, to investigate the potential functional mechanisms of these 56 DEGs, pathway enrichment analysis was conducted for these genes, and it was revealed that disulfidptosis might influence the signaling of protein digestion and absorption, the Wnt signaling pathway, ECM-receptor interaction, focal adhesion and the TGF-beta signaling pathway as based on KEGG terms. Moreover, based on the HALLMARK terms, the enrichment of epithelial mesenchymal transition, myogenesis and angiogenesis pathways was also evaluated (Figure 3C).Figure 3:Disulfidptosis impacted key genes and pathways. (A) Different expressed genes between C1 and C3 subtype, which representing the best and worst prognostic; (B) specific expression profiles of differentially expressed genes in C1 and C3 subtypes; (C) activated signaling pathways of the identified genes.Ascertain disulfidptosis-associated genes correlated with patient prognosisThe prognostic significance of all genes within the TARGET-OS and GSE21257 cohorts were evaluated, respectively. A total of 356 risk-associated genes (HR>1, p<0.05) and 704 protective genes (HR<1, p<0.05) were identified in TARGET-OS patients (Figure 4A), and 512 risk-associated genes (HR>1, p<0.05) and 426 protective genes (HR<1, p<0.05) were uncovered in GSE21257 patients (Figure 4B). Following the screening for genes exhibiting similar trends in both two cohorts, nine osteosarcoma prognosis-associated risk genes (Figure 4C) and 16 protective genes were identified (Figure 4D). Ultimately, it was discovered that six of these 25 genes were associated with cellular disulfidptosis (Figure 4E).Figure 4:Identifying prognostic genes associated with disulfidptosis. (A) Isolating risky and protective genes correlated with osteosarcoma patient overall survival in TARGET-OS cohort; (B) isolating risky and protective genes correlated with osteosarcoma patient overall survival in GSE21257 cohort; (C) Venn diagram depicting the coexisting risky genes within TARGET-OS and GSE21257 cohorts; (D) Venn diagram depicting the coexisting protective genes within TARGET-OS and GSE21257 cohorts; (E) illustrating risky genes, protective genes, disulfidptosis-related genes, and their intersections.Construction of a prognostic signature based disulfidptosis-associated genesThe six identified genes were incorporated into a LASSO regression to further screen osteosarcoma prognosis-associated genes and establish a prognostic prediction signature. In the LASSO path plot (Figure 5A), each line represents a gene in the model. The X-axis represents the value of the penalty parameter (lambda) used in the LASSO algorithm, and the Y-axis represents the value of the standardized regression coefficient for each variable. The plot shows how the magnitude of the regression coefficients changes as the penalty parameter increases. When the minimal lambda value is 0.0332, the partial likelihood deviance meets the minimal value, which reflects the optimal formula and prevents overfitting (Figure 5B). Ultimately, of the six inputted genes, five were selected by the LASSO regression model, culminating in the construction of a prognostic signature; a disulfidptosis-related (DSPR) prognostic signature. DSPR=0.47887788 × EI24 + 0.22775366 × FAIM + 0.41722950 × FAM81A − 0.63659820 × RAB8B − 0.05951149 × GNG12. Subsequently, risk scores were computed based on the signature for each patient in the TARGET-OS cohort and were presented alongside their prognostic survival status and the expression patterns of each included gene (Figure 5C).Figure 5:Developing the disulfidptosis-related (DSPR) prognostic signature. (A) The LASSO regression model, coupled with the cross-validation technique, was employed for gene screening; (B) the gene regression coefficient map within the LASSO model; (C) risk distribution diagram depicting the risk scores, alive status, and gene expression profiles for each patient within the TARGET-OS cohort.Patients were subsequently stratified into low and high score groups, and it was observed that those with high-DSPR scores encountered adverse overall survival prognosis, exhibiting a hazard ratio of 4.37 compared to the low-DSPR score cohort, with a 95 % confidence interval ranging from 1.837 to 10.388 (Figure 6A). A subsequent evaluation of the overall clinical outcome predictive value was conducted using ROC curves, revealing AUC values of 0.816 at 1-year, 0.740 at 2-year, and 0.767 at 3-year follow-up (Figure 6B). Upon excluding the potential influence of other clinical parameters, such as patient age, gender, metastatic presence at diagnosis, and primary tumor location, it was determined that the DSPR signature serves as an independent prognostic indicator (HR=4.407, 95 % CI: 1.797–10.810, p=0.0012), in addition to metastatic status at diagnosis (HR=5.479, 95 % CI: 2.515–12.01, p<0.001) (Figure 6C).Figure 6:Evaluating the prognostic predictive utility of the DSPR signature for osteosarcoma patients. (A). Kaplan–Meier survival plot illustrating the prognostic disparities among patients with high or low risk scores; (B) ROC curve demonstrating the DSPR signature’s predictive capacity for overall survival at 1-year, 3-year, and 5-year intervals; (C) forest plot showcasing the prognostic prediction ability of the DSPR signature for osteosarcoma patients after adjusting for clinical and pathological parameters; (D) developing a prognostic prediction nomogram based on clinical indicators and the DSPR signature.Nomogram reflects OS with clinical features and DSPR signatureSubsequently, a nomogram was constructed to delineate the prognosis of osteosarcoma patients, encompassing factors, such as age, gender, metastatic presence at diagnosis, primary tumor location, and DSPR signature (Figure 6D). For instance, a 17-year-old female osteosarcoma patient diagnosed with a tumor in the leg/foot, presenting with metastasis would receive a cumulative score of 270, which corresponds to a 1-year mortality risk of 0.345, a 3-year mortality risk of 0.975, and a 5-year mortality risk of 0.994 (Figure 6D). After calculating the scores of all patients using the nomogram, patients were divided into high-point and low-point groups based on the median value. It was discovered that patients with higher points exhibited a poorer prognosis (HR=4.370, 95 % CI: 1.837–10.388, p=0.0012, Figure 7A). The predictive accuracy of the nomogram was further enhanced, as evidenced by the AUC values of 0.952 at 1-year, 0.890 at 2-year, and 0.873 at 3-year (Figure 7B). Calibration curves also demonstrated a high degree of concordance between the prognostic predictions of the nomogram and actual outcomes, with Hosmer-Lemeshow p values of 0.623, 0.887, and 0.593 at 1-year, 3-year, and 5-year, respectively, thus corroborating the nomogram’s commendable ability to predict clinical outcomes (Figure 7C). DCA revealed that when the threshold probability exceeded 16 %, employing the DSPR signature for predicting mortality risk groups resulted in equivalent, if not superior, benefits compared to treating all or none of the patients (Figure 7D). The clinical impact curve demonstrated that administering treatment to all patients with a predicted probability of clinical death greater than 70 % had the most significant clinical impact (Figure 7E). Moreover, to further assess the prognostic value the C-index was calculated, and it was observed that the nomogram presented a C-index high to 0.830, while the other features were less than 0.7, and the metastatic status at diagnosis (C-index: 0.696), thereby adding more to the overall index after DSPR signature (C-index: 0.719) (Figure 7F). A total of 50 patients were randomly selected from the TARGET-OS cohort and the predictive performance of the nomogram was evaluated using the C-index. The data showed that the newly constructed nomogram demonstrated a stable prognostic predictive ability (Figure 7G).Figure 7:Assessing the clinical applicability of the nomogram. (A) Kaplan–Meier survival plot illustrating the overall survival disparities between high-point and low-point patient groups determined by the nomogram; (B) ROC curve demonstrating the nomogram’s predictive capacity for overall survival at 1-year, 3-year, and 5-year intervals; (C) calibration curves depict the concordance between the nomogram’s predicted probability of mortality and the actual probability at the end of 1-year, 3-year and 5-year; (D) decision curve analysis for the nomogram; (E) clinical impact curve for the nomogram; (F) C-index showing the prognostic value of nomogram and other features; (G) C-index for the randomly selected 50 patients from the TARGET-OS cohort.The DSPR signature meets a preferable prognostic value in the external GSE21257 cohortUtilizing the DSPR signature formula derived from LASSO analysis in the TARGET-OS cohort, the risk scores were computed for each patient in the GSE21257 cohort, the results are presented in Figure 8A, and the expression profiles of the five genes are depicted in Figure 8B. Patients were divided into low- and high-score groups, with the high-score group encountering less favorable overall survival (OS) outcomes compared to the low-score group (HR=4.000, 95 % CI: 1.517–10.183, p=0.004, Figure 8C). Additionally, the AUC values demonstrated satisfactory results, registering 0.824 at 1-year, 0.818 at 2-year, and 0.871 at 3-year (Figure 8D). Upon conducting multifactorial Cox regression analysis to eliminate the influence of other factors on patient prognosis, it was observed that, after adjusting for gender, age, HUVOS grade, metastatic presence at diagnosis, and primary tumor location, the DSPR signature served as an independent prognostic indicator in the GSE21257 cohort as well (HR=4.512, 95 % CI: 1.308–15.56, p=0.0171, Figure 8E).Figure 8:Validating the prognostic prediction performance of the DSPR signature in an external cohort. (A and B) Risk distribution diagram depicting the risk scores and gene expression profiles for each patient within the GSE21257 cohort; (C) Kaplan–Meier survival plot illustrating the prognostic disparities among patients with high or low risk scores; (D) ROC curve demonstrating the DSPR signature’s predictive capacity for overall survival at 1-year, 3-year, and 5-year intervals; (E) Forest plot showcasing the prognostic prediction ability of the DSPR signature for osteosarcoma patients after adjusting for clinical and pathological parameters.Elucidating the potential molecular mechanisms underlying distinct point subgroups of patientsTo comprehend the potential molecular characteristics of patients classified into high-point and low-point subgroups, or in other words, with disparate clinical outcomes, clusterProfiler was utilized to enrich the signaling pathways for each subgroup, respectively. It was observed that the high-point group met the activation of biomineral biomineralization of bone ossification, extracellular encapsulating, structure organization, cartilage neuron fate commitment and the response to copper ion regulation (Figure 9A). Regarding the other patients, the favorable prognosis might be regulated by cell-substrate adhesion system development, adhesion-dependent cellular to nutrient vitamin, immune activation, and complement humoral immunity (Figure 9B).Figure 9:Investigating potential signaling pathways influencing the prognosis of osteosarcoma patients, and predicting the effective chemotherapeutic agents. (A) Enrichment results of highly expressed genes in the high points group within the background of GO terms; (B) enrichment results of highly expressed genes in the low points group within the background of GO terms; (C) predicting analysis of efficacious chemotherapeutic agents conducted in TARGET-OS cohort; (D) predicting analysis of efficacious chemotherapeutic agents conducted in GSE21257 cohort.Identifying chemotherapy drugs suitable for the treatment of patients in distinct subgroupsTo provide guidance for the clinical treatment of osteosarcoma, the ‘pRRophetic’ R package was employed to analyze the therapeutic disparities of 86 potential chemotherapeutic agents across various subgroups. It was discovered that, within the TARGET-OS cohort, patients in the low-point group benefited more from LFM-A13 treatment, while those in the high-point group were better suited for BMS-509744, BMS-536924, Crizotinib, Cyclopamine, GNF-2, Imatinib, KIN001-135, PHA-665752, Salubrinal, and Sorafenib therapies (all p<0.05, Figure 9C). Analogous findings were observed in the GSE21257 cohort (all p<0.05, Figure 9D).Patients in the low-risk group exhibit greater immune cell infiltration and may derive more substantial benefits from immunotherapyIn further studies, we evaluated the differences in immune cell infiltration levels between high- and low-risk groups. Our data showed that in the TARGET-OS cohort, low-risk patients exhibited greater immune cell infiltration, including CD56 bright natural killer cells, central memory CD8 T cells, immature B cells, macrophages, and natural killer cells (Figure 10A). Similarly, in the GSE21257 cohort, it was observed that low-risk patients showed an increase in macrophage infiltration (Figure 10B). Consequently, we attempted to assess the response of patients in different risk groups to immunotherapy. Based on submap analysis, it was found that low-risk patients could benefit more from anti-PD-1 immunotherapy, a result that was validated in both cohorts (Figure 10C, D).Figure 10:Differential immune cell infiltration and response to immunotherapy across subgroups. (A) More immunocytes infiltration presenting in low risk group of TARGET-OS cohort; (B) more immunocytes infiltration presenting in low risk group of GSE21257 cohort; (C) patients in low-risk group of TARGET-OS cohort benefiting more from anti-PD-1 therapy; (D) patients in low-risk group of GSE21257 cohort benefiting more from anti-PD-1 therapy.DiscussionProgrammed cell death pathways, including apoptosis, autophagy, necrosis, and ferroptosis, have received extensive attention for their critical roles in regulating cancer cell growth [23]. It is well-documented that the Bcl-2 family is a major participant in the intrinsic apoptotic pathway, and pro-apoptotic proteins are released from the mitochondrial intermembrane space into the cytoplasm due to mitochondrial membrane disruption [24]. Various molecules, such as miR-9, Lnc-MAP6-1:3, and Siglec-15, have been reported to affect cell apoptosis and promote osteosarcoma progression by regulating Bcl-2 expression [25], [26], [27]. Furthermore, ferroptosis, a non-apoptotic type of cell death, has gained attention, the biosynthesis of glutathione or glutathione-dependent antioxidant enzyme glutathione peroxidase 4 is inhibited upon activation of ferroptosis [28]. FANCD2 has been reported to inhibit ferroptosis in osteosarcoma by regulating the JAK2/STAT3 pathway [29]. On the other hand, miRNA-1287-5p has been shown to directly bind to the 3′-untranslated region of GPX4 to suppress its protein level and activity, leading to ferroptosis induction and tumor suppression [30]. Several death-associated signaling cascades may be activated within neoplastic cells, and the interplay between programmed cellular demise pathways can dictate the destiny of tumoral cells and the clinical consequence of therapeutic endeavors. Recently, Liu et al. reported a novel form of cell death called disulfidptosis. Partifularly, this form of cell death is initiated when cells exhibiting elevated concentrations of SLC7A11 protein undergo glucose deprivation. In preclinical models, administration of a glucose inhibitor can provoke disulfidptosis in neoplastic cells overexpressing SLC7A11, thereby efficaciously impeding tumor proliferation without eliciting substantial toxicity in healthy tissues.Tumor subtyping and prognostic signature are wildly identified for the understanding of underlying tumor heterogeneity and the prediction of diverse clinical outcome for several types of tumors [31, 32]. This study aimed to identify disulfidptosis-associated genes and their potential impact on osteosarcoma prognosis. Using gene expression profiles, patients were clustered into subtypes C1, C2, and C3 based on 40 disulfidptosis-associated genes. The data demonstrated subtype C3 had the most favorable prognosis, while subtype C1 had the least favorable outcome. DEGs between subtypes C1 and C3 were compared, identifying 56 genes with potential functional mechanisms of action. It was observed that the higher expression of disulfidptosis genes PCBP1, DHX9, and IPO7 significantly associated with the poor prognosis of C3 subtype. Huang et al. [33] documented that PCBP1 extensively modulates the alternative splicing and expression of genes abundant in oncogenic pathways, encompassing the extracellular matrix, cellular adhesion, diminutive molecule metabolic processes, and apoptosis. DHX9 is a type of RNA-binding protein, which has been reported to regulate the growth of glioma and impact the infiltration of tumor-associated macrophages in the microenvironment [34]. Overexpression of IPO7 facilitated the malignant phenotypes of pancreatic cancer cells, which might act through the axis of the IPO7/p53/LncRNA MALAT1/MiR-129-5p pathway [35]. Pathway enrichment analysis revealed the influence of disulfidptosis on various signaling pathways and cellular processes. Further analysis identified nine osteosarcoma prognosis-associated risk genes and 16 protective genes; six of these 25 genes being associated with disulfidptosis.Subsequently, a prognostic signature, termed a DSPR prognostic signature, was constructed using LASSO regression, and its predictive value was assessed through risk scores and survival outcomes. The DSPR signature is composed by five genes, including EI24, FAIM, FAM81A, RAB8B, and GNG12. High levels of FAIM, EI24 and FAM81A were observed in high-risk patients, while a low level of GNG12 and RAB8B links to poor prognosis. Yuan et al. [36] likewise conveyed that GNG12 serves as a potential biomarker for osteosarcoma prognosis that is suppressed in metastatic osteosarcoma relative to non-metastatic osteosarcoma. FAIM facilitated the tetrameric assembly of GAC by augmenting PRKCE/PKCε-mediated phosphorylation and concurrently stabilize GAC by sequestering it from ClpXP protease-mediated degradation. These impacts enhanced the generation of α-ketoglutarate, leading to the activation of MTOR [37]. Methylation of the 3′ UTR in FAM81A has been reported to have a close association with the stemness of pancreatic cancer and contributes to tumoral immune infiltration [38]. The DSPR signature was found to be an independent prognostic indicator in TARGET-OS cohort, as well as in the external cohort GSE21257.In recent days, along with the thorough study of tumor mechanisms, more and more precision therapy strategies were applied for the clinical treatment, such as some molecular drug chemotherapy and targeted immunotherapy. To understand the potential molecular characteristics of patients with different clinical outcomes, signaling pathways were enriched for each subgroup. The results suggested that the high-point group showed activation of bone ossification and other processes, while the low-point group had a favorable prognosis regulated by cell-substrate adhesion, immune activation, and other factors. Suitable chemotherapy drugs are essential for clinical treatment based on the target, due to the increase rate of drug resistance to typical medicine [39, 40]. In the current study, potential chemotherapy drugs suitable for distinct subgroups were identified, which indicated that patients in the low-point group benefited more from LFM-A13 treatment, while those in the high-point group were better suited for a range of other therapies. The infiltration of immunocytes, especially the activation status of the immune microenvironment, reflects the response to anti-PD-1 immunotherapy [41], while the immunosuppressive microenvironment has shown to lead to poor prognosis of tumor [42, 43]. In the current study, it was observed that DSPR determined low-risk patients are more suitable for anti-PD-1 immunotherapy.While the study described above provides valuable insights into the potential impact of disulfidptosis-associated genes on osteosarcoma prognosis, there are several limitations to consider. First, the sample size may not be large enough to generalize the findings to the broader population. Although an external cohort was used to validate the results, further studies with larger sample sizes are needed to confirm the findings. Second, this study relies heavily on computational and statistical methods, which may introduce potential bias. For instance, the identification of disulfidptosis-associated genes could be influenced by the choice of the clustering method, gene selection criteria, and the statistical threshold. Third, in this study, potential chemotherapy drugs suitable for distinct subgroups were identified, but these findings are based on gene expression profiling and pathway analysis. Therefore, further studies are needed to validate these findings in preclinical and clinical settings.ConclusionsIn conclusion, in this study, we identified disulfidptosis-associated genes and studied their roles in osteosarcoma prognosis, constructed a prognostic signature, and provided guidance for personalized treatment strategies based on the identified subgroups. http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png ONCOLOGIE de Gruyter

Unveiling and validation of a disulfidptosis determined prognostic model for osteosarcoma: new insights from prognosis to immunotherapy and chemotherapy

ONCOLOGIE , Volume 25 (4): 17 – Jul 1, 2023

Loading next page...
 
/lp/de-gruyter/unveiling-and-validation-of-a-disulfidptosis-determined-prognostic-fR0yGmOTwR

References (45)

Publisher
de Gruyter
Copyright
© 2023 the author(s), published by De Gruyter, Berlin/Boston
eISSN
1765-2839
DOI
10.1515/oncologie-2023-0129
Publisher site
See Article on Publisher Site

Abstract

IntroductionOsteosarcoma represents one of the most prevalent primary bone malignancies, exhibiting elevated metastasis and recurrence rates, and predominantly afflicts children, adolescents, and young adults [1], [2], [3]. In the United States alone, the incidence of osteosarcoma surpasses 1,000 novel cases, annually [4]. The 5-year event-free survival rate for patients afflicted with nonmetastatic osteosarcoma is approximately 60 %, whereas for patients presenting with a primary tumor site, multiple pulmonary nodules, or other metastatic lesions, this falls below 20 % [5, 6]. Modifications in TP53 and RB1 tumor suppressor genes, the bone microenvironment, and osteosarcoma predisposition have been posited to be linked to the tumorigenesis and metastasis of osteosarcoma [7]. The optimal, standard-of-care therapeutic intervention for osteosarcoma employs high-dose methotrexate, doxorubicin, and cisplatin (MAP) [3, 7]. Additional treatments, including surgical intervention and radiotherapy, have been shown to provide localized modality control for metastatic sites [6]. Nevertheless, the scarcity of osteosarcoma has resulted in a paucity of dependable markers, thereby impeding the establishment of systematic population-based osteosarcoma screening programs as a component of preventative strategies [6]. Consequently, the discernment of pivotal pathways and genes assumes a crucial role in osteosarcoma diagnosis and the exploration of innovative therapeutics.Disulfides denote organic sulfur compounds encompassing disulfide linkages, predominantly occurring in cystine within living organisms [8]. Elevated intracellular disulfide formation can augment the levels of reactive oxygen species (ROS) within the mitochondria of cancer cells, thereby potentially inducing mitochondrial dysfunction and subsequently activating cellular demise [9]. It has previously been proposed that the overexpression of solute carrier family 7 member 11 (SLC7A11) contributes to cellular apoptosis. SLC7A11 assists cancer cells in the uptake of cysteine, thereby sustaining intracellular glutathione levels and safeguarding cells from oxidative stress-induced cell death [10]. SLC7A11 is commonly overexpressed in cancer and induces cell death by converting cystine to the more soluble cysteine, which necessitates a substantial amount of NADPH and renders cancer cells to be reliant on the pentose phosphate pathway (PPP) [10], [11], [12]. Restricting glucose availability to tumor cells augments cystine accumulation and exacerbates the depletion of NADPH and GSH [10, 13]. Pharmaceutical combinations that simultaneously target GLUT1 and GSH synthesis may induce tumor cell death, thus offering a novel potential therapeutic strategy for osteosarcoma treatment [10, 13].Recently, Liu et al. reported disulfidptosis, which is different from apoptosis and ferroptosis [14]. In this study, genes were procured that encode proteins associated with augmented disulfide bonds under conditions of NADPH depletion. Disulfidptosis-associated genes were applied to facilitate the subtyping of osteosarcoma, yielding three subgroups with a distinct prognosis. Subsequently, differentially-expressed genes were identified between best survival and worst survival groups, reflecting tumor alterations post-disulfidptosis. Subsequently, genes indicative of osteosarcoma patient prognosis were screened and prognostic models were established.MethodsBaseline information of cohortsIn this study, two osteosarcoma assemblages participated. For the Therapeutically Applicable Research to Generate Effective Treatments Osteosarcoma (TARGET-OS) cohort, the gene expression profiles of 84 patients were acquired from the UCSC Xena platform (http://xena.ucsc.edu/), which were accompanied by the patient’s clinical information. Furthermore, the GSE21257 collection, encompassing 53 patients, was obtained from the Gene Expression Omnibus platform (GEO, http://www.ncbi.nlm.nih.gov/geo/). The baseline information for patients in both cohorts is presented in Table 1.Table 1:Population demographics and prognostic outcomes of enrolled patients.GSE21257 (n=53)TARGET-OS (n=84)p-ValueAge, years Mean (SD)18.7 (12.2)15.0 (4.83)0.0379a Median [Min, Max]16.7 [3.08, 79.0]14.4 [3.56, 32.4]Gender Female19 (35.8 %)37 (44.0 %)0.376 Male34 (64.2 %)47 (56.0 %)Tumor site Arm/hand8 (15.1 %)6 (7.1 %)0.145 Leg/foot44 (83.0 %)76 (90.5 %) Pelvis0 (0 %)2 (2.4 %) Unknown1 (1.9 %)0 (0 %)Metastatic at diagnosis No39 (73.6 %)63 (75.0 %)0.844 Yes14 (26.4 %)21 (25.0 %)Alive status Alive30 (56.6 %)57 (67.9 %)0.205 Dead23 (43.4 %)27 (32.1 %)Overall survival time, months Mean (SD)68.5 (59.3)48.9 (34.0)0.0316a Median [Min, Max]45.0 [4.00, 246]47.9 [2.43, 143]ap-Value<0.05.Collection of disulfidptosis-associated genesLiu et al. reported a novel type of disulfidptosis. Cells exhibiting elevated SLC7A11 levels after being subjected to glucose scarcity showed augmented cystine uptake. In line with an inadequate NADPH provision, this culminates in NADPH depletion, aberrant disulfide bonding within actin cytoskeleton proteins, actin lattice disintegration, and consequently cellular demise [14]. Ninety cysteine sites were documented with a 1.5-fold augmentation in disulfide bonds subsequent to glucose deprivation, encoded by 77 genes, the ID of the genes were collected from Supplementary Table 2 of Liu et al.’s study.Identifying distinct disulfidptosis subtypesThe consensus clustering and molecular subtype discernment were conducted with the R package named “ConsensusClusterPlus” [15]. Fifty iterations were implemented to reveal the cluster k-means, of which at least 80 % of samples were enrolled each time. The score of cumulative distribution function (CDF), which means the area that under the CDF curve, was calculated for the selection of the optimal cluster number. The principal component analysis (PCA) was harnessed to substantiate the credibility of the consensus clusters, as well as the t-distributed stochastic neighbor embedding (t-SNE).Gene set enrichment analysisDifferentially-expressed genes (DEGs) are important for identification of different mechanisms among subgroups and were collected by the ‘limma’ package in R. An absolute fold-change exceeding 1.5 and an adjusted p-value below 0.01 were established as threshold values for DEGs filtration. Packages ‘org.Hs.e.g.db’ and ‘msigdbr’ were utilized to conduct Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG), and HALLMARK enrichment analyses [16, 17], while ‘ClusterProfiler’ was employed to explicate representative signaling pathways [18].Constructing the prognostic signatureUnivariate Cox examination was executed to investigate prognostic genes in both TARGET-OS and GSE21257 cohorts. Subsequently, a Venn diagram was amalgamated, depicting concordant prognostic genes for ensuing analysis. The least absolute shrinkage and selection operator (LASSO) regression analysis, executed by the ‘glmnet’ package in R [19], was harnessed to refine and regularize all input factors, ultimately yielding a statistical model with exceptional predictive precision and interpretability. LASSO, a regularization technique, can avert overfitting, conduct feature selection more proficiently, and optimize the structure of the regression model more effectively compared to traditional ordinary least squares subset selection methods. It is particularly suitable for situations involving numerous potential predictors where only a few are presumed pertinent. Genes selected by LASSO analysis were employed to calculate the risk score of each patient, thereby utilizing gene expression values and coefficients to devise the following disulfidptosis-related (DSPR) prognostic signature:DSPR=∑i=1nCi×Ei$$\mathrm{D}\mathrm{S}\mathrm{P}\mathrm{R}=\sum\limits _{i=1}^{n}{C}_{i}\times {E}_{i}$$where Ci${C}_{i}$denotes the coefficient of a gene, and Ei${E}_{i}$represents the expression value of the corresponding gene.Establishing a predictive nomogramIn this study, we embarked on a comprehensive multivariate Cox regression examination as the initial endeavor, with the objective of pinpointing autonomous prognostic determinants. Subsequently, the DSPR signature was integrated into the analytical framework. Next, a nomogram was created showing the autonomous prognostic constituents to amalgamate the prognosticated values. The “rms” and “regplot” R packages were employed for this endeavor. To ascertain the accuracy and reliability of our study, a calibration diagram, decision curve scrutiny (DCA), and clinical impact curve assessment were executed.Personalized chemotherapy predictionBy leveraging the pharmacological susceptibility and phenotypic data from GDSC 2016 (https://www.cancerrxgene.org/), the R package “pRRophetic” was employed to predict the chemotherapeutic responsiveness for each patient [20]. Estimated IC50 values (lower IC50 indicates enhanced therapy susceptibility) for each specimen exposed to a specific chemotherapeutic agent were obtained via ridge regression [21], and the prediction accuracy was ascertained through a tenfold cross-validated process [22].StatisticsThe T-test or Mann-Whitney U test was employed to evaluate disparities between subgroups with continuous variables. Survival curves were created utilizing Kaplan-Meier survival analysis and log-rank analysis, while the hazard ratio (HR) and 95 % confidence interval (CI) for candidate genes were ascertained via a univariate Cox regression model. Multivariate Cox analyses were employed to determine the independent prognostic influence of the risk score, subsequent to adjustments for assorted clinical characteristics. Adopting p<0.05 as the criterion for significance, two-tailed statistical evaluation was conducted. R software (4.1.2) was employed for all analytical analyses.ResultsUnveiling the disulfidptosis-defined subcategoriesAs mentioned in the methods section, 79 disulfidptosis associated genes were collected, among which 40 could be identified in both the TARGET-OS cohort and GSE21257 cohort. Therefore, the gene expression profiles of the 40 genes were collected as the input matrix for clustering. The consensus matrices for k=2, 3, 4, 5, 6 were tested (Figure 1A). With the cumulative distribution function (CDF) curve, it was observed that while the area under the CDF did not change significantly with the number of k after k=3 (Figure 1B), and that the relative change of the CDF between k=3 and k=4 showed the biggest alteration. Therefore, it was recognized that k=3 was the best cluster number for subsequent analyses (Figure 1C).Figure 1:Identification the optimal number of clusters by ConsensusClusterPlus. (A) Heatmaps of the consensus matrices for k=2, 3, 4, 5, 6; (B) the cumulative distribution functions of the consensus matrix for each k; (C) the relative change in area under the CDF curve comparing k and k − 1.In this study, the overall survival (OS) prognosis of patients with distinct subtypes was evaluated. Subtype C1 comprised 29 individuals, among whom 15 patients succumbed, while subtype C2 included 41 individuals with 12 fatalities. Subtype C3 encompassed 15 individuals, all of whom survived until the conclusion of the follow-up period. Consequently, among the three subtypes delineated by disulfidptosis-associated genes, subtype C3 exhibited the most favorable prognosis, while subtype C1 demonstrated the least favorable outcome, with a marked disparity in OS prognosis among the three patient cohorts (overall: p=0.0048; C1 vs. C2: p=0.0550; C1 vs. C3: p=0.0072; C2 vs. C3: p=0.0550 Figure 2A). Patients within the three subtypes were both mutually independent and intersected among other subtypes (Figure 2B). Furthermore, an increased prevalence of tumors situated in the arm/hand region among subtype C2 patients was observed, while all tumors in subtype C3 were located in the leg/foot. Regarding disulfidptosis-associated genes, it was observed that subtype C2 exhibited elevated levels of FLNC and MYH14. High expression of MYH4, MYH6, MYH8, MYH13, CHCHD3, and PRDX1 in subtype C1 could potentially contribute to an unfavorable prognosis. In the context of subtype C3, heightened expression of TLN1, ZHX2, MYH9, and MYH10 was observed (Figure 2C). In addition, immune infiltration among three subtypes was evaluated, and it was observed that the best prognosis C3 subtype presented the highest infiltration of immunocytes, especially for CD56 bright natural killer cells, central memory CD8 T cells, eosinophil,s gamma delta T cells, macrophages, and myeloid-derived suppressor. Cells (MDSC). Regarding the C1 subtype with the worst prognosis, immunocyte infiltration was at an intermediate level (Figure 2D).Figure 2:Revealing three disulfidptosis determined molecular subtypes. (A) Diverse overall survival outcome of three subtypes; (B) separation and intersection of the distribution of three subtypes; (C) expression levels of 40 disulfidptosis genes across different patients; (D) differential immune cell infiltration among three subtypes.Determine the disulfidptosis-influenced pivotal genes and pathwaysBecause of apparent differences in OS prognosis between subtypes C1 and C3, DEGs within these two subcategories were compared (Figure 3A). Employing a screening threshold of an absolute log2 (fold change) greater than 1.5 and a p-value less than 0.05, 56 DEGs were identified; the expression patterns of which in distinct patients are depicted in Figure 3B. Furthermore, to investigate the potential functional mechanisms of these 56 DEGs, pathway enrichment analysis was conducted for these genes, and it was revealed that disulfidptosis might influence the signaling of protein digestion and absorption, the Wnt signaling pathway, ECM-receptor interaction, focal adhesion and the TGF-beta signaling pathway as based on KEGG terms. Moreover, based on the HALLMARK terms, the enrichment of epithelial mesenchymal transition, myogenesis and angiogenesis pathways was also evaluated (Figure 3C).Figure 3:Disulfidptosis impacted key genes and pathways. (A) Different expressed genes between C1 and C3 subtype, which representing the best and worst prognostic; (B) specific expression profiles of differentially expressed genes in C1 and C3 subtypes; (C) activated signaling pathways of the identified genes.Ascertain disulfidptosis-associated genes correlated with patient prognosisThe prognostic significance of all genes within the TARGET-OS and GSE21257 cohorts were evaluated, respectively. A total of 356 risk-associated genes (HR>1, p<0.05) and 704 protective genes (HR<1, p<0.05) were identified in TARGET-OS patients (Figure 4A), and 512 risk-associated genes (HR>1, p<0.05) and 426 protective genes (HR<1, p<0.05) were uncovered in GSE21257 patients (Figure 4B). Following the screening for genes exhibiting similar trends in both two cohorts, nine osteosarcoma prognosis-associated risk genes (Figure 4C) and 16 protective genes were identified (Figure 4D). Ultimately, it was discovered that six of these 25 genes were associated with cellular disulfidptosis (Figure 4E).Figure 4:Identifying prognostic genes associated with disulfidptosis. (A) Isolating risky and protective genes correlated with osteosarcoma patient overall survival in TARGET-OS cohort; (B) isolating risky and protective genes correlated with osteosarcoma patient overall survival in GSE21257 cohort; (C) Venn diagram depicting the coexisting risky genes within TARGET-OS and GSE21257 cohorts; (D) Venn diagram depicting the coexisting protective genes within TARGET-OS and GSE21257 cohorts; (E) illustrating risky genes, protective genes, disulfidptosis-related genes, and their intersections.Construction of a prognostic signature based disulfidptosis-associated genesThe six identified genes were incorporated into a LASSO regression to further screen osteosarcoma prognosis-associated genes and establish a prognostic prediction signature. In the LASSO path plot (Figure 5A), each line represents a gene in the model. The X-axis represents the value of the penalty parameter (lambda) used in the LASSO algorithm, and the Y-axis represents the value of the standardized regression coefficient for each variable. The plot shows how the magnitude of the regression coefficients changes as the penalty parameter increases. When the minimal lambda value is 0.0332, the partial likelihood deviance meets the minimal value, which reflects the optimal formula and prevents overfitting (Figure 5B). Ultimately, of the six inputted genes, five were selected by the LASSO regression model, culminating in the construction of a prognostic signature; a disulfidptosis-related (DSPR) prognostic signature. DSPR=0.47887788 × EI24 + 0.22775366 × FAIM + 0.41722950 × FAM81A − 0.63659820 × RAB8B − 0.05951149 × GNG12. Subsequently, risk scores were computed based on the signature for each patient in the TARGET-OS cohort and were presented alongside their prognostic survival status and the expression patterns of each included gene (Figure 5C).Figure 5:Developing the disulfidptosis-related (DSPR) prognostic signature. (A) The LASSO regression model, coupled with the cross-validation technique, was employed for gene screening; (B) the gene regression coefficient map within the LASSO model; (C) risk distribution diagram depicting the risk scores, alive status, and gene expression profiles for each patient within the TARGET-OS cohort.Patients were subsequently stratified into low and high score groups, and it was observed that those with high-DSPR scores encountered adverse overall survival prognosis, exhibiting a hazard ratio of 4.37 compared to the low-DSPR score cohort, with a 95 % confidence interval ranging from 1.837 to 10.388 (Figure 6A). A subsequent evaluation of the overall clinical outcome predictive value was conducted using ROC curves, revealing AUC values of 0.816 at 1-year, 0.740 at 2-year, and 0.767 at 3-year follow-up (Figure 6B). Upon excluding the potential influence of other clinical parameters, such as patient age, gender, metastatic presence at diagnosis, and primary tumor location, it was determined that the DSPR signature serves as an independent prognostic indicator (HR=4.407, 95 % CI: 1.797–10.810, p=0.0012), in addition to metastatic status at diagnosis (HR=5.479, 95 % CI: 2.515–12.01, p<0.001) (Figure 6C).Figure 6:Evaluating the prognostic predictive utility of the DSPR signature for osteosarcoma patients. (A). Kaplan–Meier survival plot illustrating the prognostic disparities among patients with high or low risk scores; (B) ROC curve demonstrating the DSPR signature’s predictive capacity for overall survival at 1-year, 3-year, and 5-year intervals; (C) forest plot showcasing the prognostic prediction ability of the DSPR signature for osteosarcoma patients after adjusting for clinical and pathological parameters; (D) developing a prognostic prediction nomogram based on clinical indicators and the DSPR signature.Nomogram reflects OS with clinical features and DSPR signatureSubsequently, a nomogram was constructed to delineate the prognosis of osteosarcoma patients, encompassing factors, such as age, gender, metastatic presence at diagnosis, primary tumor location, and DSPR signature (Figure 6D). For instance, a 17-year-old female osteosarcoma patient diagnosed with a tumor in the leg/foot, presenting with metastasis would receive a cumulative score of 270, which corresponds to a 1-year mortality risk of 0.345, a 3-year mortality risk of 0.975, and a 5-year mortality risk of 0.994 (Figure 6D). After calculating the scores of all patients using the nomogram, patients were divided into high-point and low-point groups based on the median value. It was discovered that patients with higher points exhibited a poorer prognosis (HR=4.370, 95 % CI: 1.837–10.388, p=0.0012, Figure 7A). The predictive accuracy of the nomogram was further enhanced, as evidenced by the AUC values of 0.952 at 1-year, 0.890 at 2-year, and 0.873 at 3-year (Figure 7B). Calibration curves also demonstrated a high degree of concordance between the prognostic predictions of the nomogram and actual outcomes, with Hosmer-Lemeshow p values of 0.623, 0.887, and 0.593 at 1-year, 3-year, and 5-year, respectively, thus corroborating the nomogram’s commendable ability to predict clinical outcomes (Figure 7C). DCA revealed that when the threshold probability exceeded 16 %, employing the DSPR signature for predicting mortality risk groups resulted in equivalent, if not superior, benefits compared to treating all or none of the patients (Figure 7D). The clinical impact curve demonstrated that administering treatment to all patients with a predicted probability of clinical death greater than 70 % had the most significant clinical impact (Figure 7E). Moreover, to further assess the prognostic value the C-index was calculated, and it was observed that the nomogram presented a C-index high to 0.830, while the other features were less than 0.7, and the metastatic status at diagnosis (C-index: 0.696), thereby adding more to the overall index after DSPR signature (C-index: 0.719) (Figure 7F). A total of 50 patients were randomly selected from the TARGET-OS cohort and the predictive performance of the nomogram was evaluated using the C-index. The data showed that the newly constructed nomogram demonstrated a stable prognostic predictive ability (Figure 7G).Figure 7:Assessing the clinical applicability of the nomogram. (A) Kaplan–Meier survival plot illustrating the overall survival disparities between high-point and low-point patient groups determined by the nomogram; (B) ROC curve demonstrating the nomogram’s predictive capacity for overall survival at 1-year, 3-year, and 5-year intervals; (C) calibration curves depict the concordance between the nomogram’s predicted probability of mortality and the actual probability at the end of 1-year, 3-year and 5-year; (D) decision curve analysis for the nomogram; (E) clinical impact curve for the nomogram; (F) C-index showing the prognostic value of nomogram and other features; (G) C-index for the randomly selected 50 patients from the TARGET-OS cohort.The DSPR signature meets a preferable prognostic value in the external GSE21257 cohortUtilizing the DSPR signature formula derived from LASSO analysis in the TARGET-OS cohort, the risk scores were computed for each patient in the GSE21257 cohort, the results are presented in Figure 8A, and the expression profiles of the five genes are depicted in Figure 8B. Patients were divided into low- and high-score groups, with the high-score group encountering less favorable overall survival (OS) outcomes compared to the low-score group (HR=4.000, 95 % CI: 1.517–10.183, p=0.004, Figure 8C). Additionally, the AUC values demonstrated satisfactory results, registering 0.824 at 1-year, 0.818 at 2-year, and 0.871 at 3-year (Figure 8D). Upon conducting multifactorial Cox regression analysis to eliminate the influence of other factors on patient prognosis, it was observed that, after adjusting for gender, age, HUVOS grade, metastatic presence at diagnosis, and primary tumor location, the DSPR signature served as an independent prognostic indicator in the GSE21257 cohort as well (HR=4.512, 95 % CI: 1.308–15.56, p=0.0171, Figure 8E).Figure 8:Validating the prognostic prediction performance of the DSPR signature in an external cohort. (A and B) Risk distribution diagram depicting the risk scores and gene expression profiles for each patient within the GSE21257 cohort; (C) Kaplan–Meier survival plot illustrating the prognostic disparities among patients with high or low risk scores; (D) ROC curve demonstrating the DSPR signature’s predictive capacity for overall survival at 1-year, 3-year, and 5-year intervals; (E) Forest plot showcasing the prognostic prediction ability of the DSPR signature for osteosarcoma patients after adjusting for clinical and pathological parameters.Elucidating the potential molecular mechanisms underlying distinct point subgroups of patientsTo comprehend the potential molecular characteristics of patients classified into high-point and low-point subgroups, or in other words, with disparate clinical outcomes, clusterProfiler was utilized to enrich the signaling pathways for each subgroup, respectively. It was observed that the high-point group met the activation of biomineral biomineralization of bone ossification, extracellular encapsulating, structure organization, cartilage neuron fate commitment and the response to copper ion regulation (Figure 9A). Regarding the other patients, the favorable prognosis might be regulated by cell-substrate adhesion system development, adhesion-dependent cellular to nutrient vitamin, immune activation, and complement humoral immunity (Figure 9B).Figure 9:Investigating potential signaling pathways influencing the prognosis of osteosarcoma patients, and predicting the effective chemotherapeutic agents. (A) Enrichment results of highly expressed genes in the high points group within the background of GO terms; (B) enrichment results of highly expressed genes in the low points group within the background of GO terms; (C) predicting analysis of efficacious chemotherapeutic agents conducted in TARGET-OS cohort; (D) predicting analysis of efficacious chemotherapeutic agents conducted in GSE21257 cohort.Identifying chemotherapy drugs suitable for the treatment of patients in distinct subgroupsTo provide guidance for the clinical treatment of osteosarcoma, the ‘pRRophetic’ R package was employed to analyze the therapeutic disparities of 86 potential chemotherapeutic agents across various subgroups. It was discovered that, within the TARGET-OS cohort, patients in the low-point group benefited more from LFM-A13 treatment, while those in the high-point group were better suited for BMS-509744, BMS-536924, Crizotinib, Cyclopamine, GNF-2, Imatinib, KIN001-135, PHA-665752, Salubrinal, and Sorafenib therapies (all p<0.05, Figure 9C). Analogous findings were observed in the GSE21257 cohort (all p<0.05, Figure 9D).Patients in the low-risk group exhibit greater immune cell infiltration and may derive more substantial benefits from immunotherapyIn further studies, we evaluated the differences in immune cell infiltration levels between high- and low-risk groups. Our data showed that in the TARGET-OS cohort, low-risk patients exhibited greater immune cell infiltration, including CD56 bright natural killer cells, central memory CD8 T cells, immature B cells, macrophages, and natural killer cells (Figure 10A). Similarly, in the GSE21257 cohort, it was observed that low-risk patients showed an increase in macrophage infiltration (Figure 10B). Consequently, we attempted to assess the response of patients in different risk groups to immunotherapy. Based on submap analysis, it was found that low-risk patients could benefit more from anti-PD-1 immunotherapy, a result that was validated in both cohorts (Figure 10C, D).Figure 10:Differential immune cell infiltration and response to immunotherapy across subgroups. (A) More immunocytes infiltration presenting in low risk group of TARGET-OS cohort; (B) more immunocytes infiltration presenting in low risk group of GSE21257 cohort; (C) patients in low-risk group of TARGET-OS cohort benefiting more from anti-PD-1 therapy; (D) patients in low-risk group of GSE21257 cohort benefiting more from anti-PD-1 therapy.DiscussionProgrammed cell death pathways, including apoptosis, autophagy, necrosis, and ferroptosis, have received extensive attention for their critical roles in regulating cancer cell growth [23]. It is well-documented that the Bcl-2 family is a major participant in the intrinsic apoptotic pathway, and pro-apoptotic proteins are released from the mitochondrial intermembrane space into the cytoplasm due to mitochondrial membrane disruption [24]. Various molecules, such as miR-9, Lnc-MAP6-1:3, and Siglec-15, have been reported to affect cell apoptosis and promote osteosarcoma progression by regulating Bcl-2 expression [25], [26], [27]. Furthermore, ferroptosis, a non-apoptotic type of cell death, has gained attention, the biosynthesis of glutathione or glutathione-dependent antioxidant enzyme glutathione peroxidase 4 is inhibited upon activation of ferroptosis [28]. FANCD2 has been reported to inhibit ferroptosis in osteosarcoma by regulating the JAK2/STAT3 pathway [29]. On the other hand, miRNA-1287-5p has been shown to directly bind to the 3′-untranslated region of GPX4 to suppress its protein level and activity, leading to ferroptosis induction and tumor suppression [30]. Several death-associated signaling cascades may be activated within neoplastic cells, and the interplay between programmed cellular demise pathways can dictate the destiny of tumoral cells and the clinical consequence of therapeutic endeavors. Recently, Liu et al. reported a novel form of cell death called disulfidptosis. Partifularly, this form of cell death is initiated when cells exhibiting elevated concentrations of SLC7A11 protein undergo glucose deprivation. In preclinical models, administration of a glucose inhibitor can provoke disulfidptosis in neoplastic cells overexpressing SLC7A11, thereby efficaciously impeding tumor proliferation without eliciting substantial toxicity in healthy tissues.Tumor subtyping and prognostic signature are wildly identified for the understanding of underlying tumor heterogeneity and the prediction of diverse clinical outcome for several types of tumors [31, 32]. This study aimed to identify disulfidptosis-associated genes and their potential impact on osteosarcoma prognosis. Using gene expression profiles, patients were clustered into subtypes C1, C2, and C3 based on 40 disulfidptosis-associated genes. The data demonstrated subtype C3 had the most favorable prognosis, while subtype C1 had the least favorable outcome. DEGs between subtypes C1 and C3 were compared, identifying 56 genes with potential functional mechanisms of action. It was observed that the higher expression of disulfidptosis genes PCBP1, DHX9, and IPO7 significantly associated with the poor prognosis of C3 subtype. Huang et al. [33] documented that PCBP1 extensively modulates the alternative splicing and expression of genes abundant in oncogenic pathways, encompassing the extracellular matrix, cellular adhesion, diminutive molecule metabolic processes, and apoptosis. DHX9 is a type of RNA-binding protein, which has been reported to regulate the growth of glioma and impact the infiltration of tumor-associated macrophages in the microenvironment [34]. Overexpression of IPO7 facilitated the malignant phenotypes of pancreatic cancer cells, which might act through the axis of the IPO7/p53/LncRNA MALAT1/MiR-129-5p pathway [35]. Pathway enrichment analysis revealed the influence of disulfidptosis on various signaling pathways and cellular processes. Further analysis identified nine osteosarcoma prognosis-associated risk genes and 16 protective genes; six of these 25 genes being associated with disulfidptosis.Subsequently, a prognostic signature, termed a DSPR prognostic signature, was constructed using LASSO regression, and its predictive value was assessed through risk scores and survival outcomes. The DSPR signature is composed by five genes, including EI24, FAIM, FAM81A, RAB8B, and GNG12. High levels of FAIM, EI24 and FAM81A were observed in high-risk patients, while a low level of GNG12 and RAB8B links to poor prognosis. Yuan et al. [36] likewise conveyed that GNG12 serves as a potential biomarker for osteosarcoma prognosis that is suppressed in metastatic osteosarcoma relative to non-metastatic osteosarcoma. FAIM facilitated the tetrameric assembly of GAC by augmenting PRKCE/PKCε-mediated phosphorylation and concurrently stabilize GAC by sequestering it from ClpXP protease-mediated degradation. These impacts enhanced the generation of α-ketoglutarate, leading to the activation of MTOR [37]. Methylation of the 3′ UTR in FAM81A has been reported to have a close association with the stemness of pancreatic cancer and contributes to tumoral immune infiltration [38]. The DSPR signature was found to be an independent prognostic indicator in TARGET-OS cohort, as well as in the external cohort GSE21257.In recent days, along with the thorough study of tumor mechanisms, more and more precision therapy strategies were applied for the clinical treatment, such as some molecular drug chemotherapy and targeted immunotherapy. To understand the potential molecular characteristics of patients with different clinical outcomes, signaling pathways were enriched for each subgroup. The results suggested that the high-point group showed activation of bone ossification and other processes, while the low-point group had a favorable prognosis regulated by cell-substrate adhesion, immune activation, and other factors. Suitable chemotherapy drugs are essential for clinical treatment based on the target, due to the increase rate of drug resistance to typical medicine [39, 40]. In the current study, potential chemotherapy drugs suitable for distinct subgroups were identified, which indicated that patients in the low-point group benefited more from LFM-A13 treatment, while those in the high-point group were better suited for a range of other therapies. The infiltration of immunocytes, especially the activation status of the immune microenvironment, reflects the response to anti-PD-1 immunotherapy [41], while the immunosuppressive microenvironment has shown to lead to poor prognosis of tumor [42, 43]. In the current study, it was observed that DSPR determined low-risk patients are more suitable for anti-PD-1 immunotherapy.While the study described above provides valuable insights into the potential impact of disulfidptosis-associated genes on osteosarcoma prognosis, there are several limitations to consider. First, the sample size may not be large enough to generalize the findings to the broader population. Although an external cohort was used to validate the results, further studies with larger sample sizes are needed to confirm the findings. Second, this study relies heavily on computational and statistical methods, which may introduce potential bias. For instance, the identification of disulfidptosis-associated genes could be influenced by the choice of the clustering method, gene selection criteria, and the statistical threshold. Third, in this study, potential chemotherapy drugs suitable for distinct subgroups were identified, but these findings are based on gene expression profiling and pathway analysis. Therefore, further studies are needed to validate these findings in preclinical and clinical settings.ConclusionsIn conclusion, in this study, we identified disulfidptosis-associated genes and studied their roles in osteosarcoma prognosis, constructed a prognostic signature, and provided guidance for personalized treatment strategies based on the identified subgroups.

Journal

ONCOLOGIEde Gruyter

Published: Jul 1, 2023

Keywords: disulfidptosis; nomogram; osteosarcoma; overall survival; prognosis

There are no references for this article.