Osteosarcoma, also known as osteogenic sarcoma, is considered to be the most common malignancy in bone among children and adolescents [1, 2]. Although the metaphyseal region is the most common location of osteosarcoma, this disease usually progresses rapidly and is prone to metastasis [3, 4]. The prognosis of osteosarcoma is poor and it seriously threatens the life and health of adolescents [5]. Therefore, it is urgent to explore novel targets or signatures for improving the clinical practice of osteosarcoma patients in the future.

Autophagy is a biological process mediated by certain genes, in which the aberrant organelles as well as macromolecules are digested by lysosomes, and it is involved in multiple processes such as cell metabolism, renewal of organelles and intracellular homeostasis maintenance [4, 6]. In recent years, several studies have indicated the relationship of autophagy with the initiation and progression of various diseases, such as cancer and diseases associated with neurodegeneration and immunization [7]. In bladder cancer, autophagy could inhibit the tumorigenesis via limiting tissue damage and oncogenic signaling [8]. Furthermore, autophagy could eliminate the accumulation of damaged proteins and organelles, indicating its role in the prevention of tumorigenesis. However, other research has reported that autophagy was essential for enhancing the survival ability of tumor cells and suppressing the necrosis in some cancers, including melanoma and breast cancer [9, 10]. Recycling by autophagy is necessary for the maintenance of energy balance and mitochondrial metabolism for tumor growth and proliferation. In addition, suppression of autophagy is considered as a potential modality for tumor treatment [11]. Autophagy also plays a significant role in osteosarcoma. It has been proved that autophagy is induced in osteosarcoma, and several intermediates are implicated in this process [12]. Autophagy could be promoted by the increased expression of high mobility group box 1, contributing to the drug resistance during the treatment of osteosarcoma [13]. Sun et al. reported that the silence of autophagy-related gene 5 reduced the malignancy of osteosarcoma with anti-oncogenic effects [14]. Moreover, Liu et al. found that the expression of autophagy-related 4B was obviously elevated, which accelerated osteosarcoma development and suppressed the apoptosis of osteosarcoma cells [15]. Consequently, analysis and identification of autophagy-related genes are helpful to improve our knowledge on the association of autophagy with osteosarcoma.

Herein, we analyzed 210 genes associated with autophagy and identified 14 optimized autophagy-associated genes related to the survival outcome of osteosarcoma patients. A death risk model based on those 14 autophagy-related genes could effectively predict the prognosis of osteosarcoma patients. Finally, we established the nomogram model by including independent factors of prognosis (risk score and metastasis) and demonstrated its better performance in predicting the long-term prognosis of osteosarcoma patients.

Osteosarcoma is a frequent malignancy in bone among children and adolescents, with an estimated incidence rate of 3 per million per year all over the world. Osteosarcoma is often accompanied with early metastasis, thus being considered as an invasive tumor [26, 27]. Although radiotherapy and neo-adjuvant chemotherapy have been developing over the past years, the OS at 5 years in metastatic cases remains only about 20%. The cure rate for osteosarcoma patients that have focal tumor increases up to 70% [28, 29]. Therefore, it is urgent to explore and identify novel prognostic biomarkers for proper clinical decision making, which may provide a searchable idea to ameliorate treatment status and survival outcome of osteosarcoma patients.

Autophagy is self-degradation targeting defective proteins and organelles, and maintains the function of mitochondria under the condition of stress. Extensive research has reveal that autophagy is associated with tumor initiation and progression. Indeed, autophagy could contribute to suppression of tumorigenesis of liver tumor through cell-intrinsic p62 accumulation [30]. Another study found that autophagy suppressed the pancreatic tumor formation by p53 loss [31]. Nevertheless, autophagy was actually required for growth, survival, and tumorigenesis of pancreatic cancers [32]. Autophagy-related genes were upregulated in breast cancer cells transformed with RAS and promoted the invasion of cancer cells [33]. Hence, autophagy provides a potential alternative for osteosarcoma treatment.

In this study, 210 autophagy-related genes were collected which could distinguish the osteosarcoma patients with significantly different prognoses. LASSO Cox regression analysis was carried out and identified 14 optimized autophagy-related genes, including ATG4A, BAK1, BNIP3, CALCOCO2, CCL2, DAPK1, EGFR, FAS, GRID2, ITGA3, MYC, RAB33B, USP10, and WIPI1, for the prognosis of osteosarcoma. ATG4A, a redox-regulated cysteine protease, is a vital autophagy regulator. ATG4A promoted the transition from epithelium to mesenchyme partly by the Notch signaling pathway in osteosarcoma cells [34] and was related to reduced risk for lung cancer [35]. Another report showed that hypomethylation of ATG4A predicted a poor prognosis for ovarian cancer patients [36]. BAK1, a member of the B cell lymphoma family containing BH3 domain, could induce the mitochondria-mediated apoptosis by interacting with other proteins. Studies demonstrated that BAK1 played a role in drug resistance and tumor proliferation in many cancers including breast, lung and cervical cancers [37–39]. An established prognostic signature based on seven genes including BAK1 was able to predict the survival outcome of head and neck squamous cell carcinoma patients [40]. BNIP3, also a member of B cell lymphoma family, could regulate the cell survival, autophagy, and cytoprotection. Moreover, research suggested that the elevated BNIP3 levels were correlated with progression to metastasis and poor prognosis in multiple cancers, such as breast and lung cancers, and uveal melanoma [41–43]. BNIP3 overexpression could induce the apoptosis of osteosarcoma cells, and BNIP3 inhibition plays a suppressive role in osteosarcoma cells apoptosis [44, 45]. CALCOCO2, also known as nuclear domain 10 protein 52, is implicated in autophagy factors recruitment and TANK-binding kinase 1 (TBK1) activation [46]. The model based on 16 autophagy related genes including CALCOCO2 could discriminate the multiple myeloma patients with distinct clinical outcomes, presenting potential prognostic value in multiple myeloma research [47]. CCL2 belongs to the CC chemokine family and is secreted by various cells including endothelial cells, fibroblasts, monocytes and tumor cells [48]. A previous study indicated that CCL2 could promote the invasion of pancreatic ductal adenocarcinoma [49] and the metastasis in cervical cancer [50]. In addition, high levels of CCL2 were related to the inferior survival outcome in gastric cancer [51]. Compared with the low-grade osteosarcoma, CCL2 expression was elevated in the osteosarcoma with high grade, which enhanced the proliferative and invasive abilities of osteosarcoma cells [52]. DAPK1 belongs to the Ser/Thr kinase family and is considered as a key regulator of autophagy and apoptosis [53]. It was found that DAPK1 expression could significantly inhibit the tumor growth and metastasis [54]. Down-regulation of DAPK1 expression may be a prognostic factor in many tumors, such as diffuse large B-cell lymphoma [55] and liver cancer [56]. EGFR is a receptor tyrosine kinase (RTK) for ErbB family, and exhibits over-expression in various tumor cells [57]. Overexpression of EGFR is related to survival, invasion, metastasis, drug resistance, and poor prognosis of tumor [58], for example, EGFR is considered as an indicator of inferior prognosis in node-negative breast cancer [59]. EGFR was reported to be abnormally expressed in osteosarcoma, and the expression as well as amplification of EGFR were observed in the osteosarcoma with high grade (PMID: [60, 61]). FAS is an essential enzyme in the process of lipogenesis, and could effectively maintain the energy homeostasis. FAS expression level is up-regulated in several cancers and exhibits a strong effect on tumor cell proliferation and apoptosis. It was found that inhibition of FAS could obviously inhibited the capacity of growth and migration of bladder cancer cells [62], and Fas was proved to be a significant marker for the prognosis of breast cancer [63]. GRID2, belongs to the ionic glutamate receptor family, and regulates excitatory synaptic transmission [64]. Previous study has confirmed that GRID2 was related to the inferior survival outcome of prostate and gastric cancers [65, 66]. ITGA3 belongs to the integrin family of cell surface receptors and is involved in the survival, proliferation, and migration of cells. In gastric carcinomas, ITGA3 expression could facilitate the invasion [67] and was considered as a key signature for colon cancer [68]. Interestingly, ITGA3 was also a component of the prognostic signature for head and neck squamous cell carcinoma, with a similar role like BAK1, as described above [40]. ITGA3 polymorphisms might influence the osteosarcoma in terms of the incidence rate, metastatic status and prognosis, which was considered as a potential signature for osteosarcoma [69]. MYC, an oncogenic transcription factor, regulates cell proliferation, apoptosis, and carcinogenesis [70]. Research revealed that mutations of c-MYC could result in tumorigenesis [71]. The amplification of c-MYC was observed during the development of hepatocellular carcinoma, which was associated with impaired survival [72]. Moreover, increased expression of c-MYC was proved to enhance the invasive ability of osteosarcoma cells by targeting MEK-ERK pathway [73]. RAB33B belongs to the Rab family of small GTP binding proteins, and regulates the fusion of autophagosomes and membrane trafficking [74]. In addition, RAB33B was identified as a biomarker for lung cancer diagnosis [75]. However, its prognostic value or relationship with the survival outcome of cancer is rarely reported. USP10 belongs to the ubiquitin-specific protease family, modulates DNA damage response and autophagy [76, 77]. Moreover, USP10 inhibited cell growth and invasion in lung cancer [78], and was an independent factor for the prognosis of gastric carcinoma [79]. WIPI1, a member of WD-repeat protein which interacts with phosphoinositides (WIPI) family, participates in the formation of autophagosome [80]. Furthermore, WIPI1 was a relevant novel melanoma marker [81], and the increased expression of WIPI1 indicated poor clinical outcome in uveal melanoma [82] In osteosarcoma, WIPI1 expression was obviously elevated, which promoted the proliferation of osteosarcoma cells through regulating CDKN1A expression [83]. This research further confirmed the potential prognostic value of the identified genes in osteosarcoma.

Several prognostic models of osteosarcoma have been established in previous research. For example, Qu et al. constructed a 5-gene-signature for the prognosis prediction of osteosarcoma based on the super-enhancer-associated genes [84] The model basing on the biomarkers including RBC, PNI, CRE, Ca²⁺and LSR in blood presented good performance in predicting the overall survival of osteosarcoma patients [85]. Lin et al. established a predictive model with five differentially expressed genes related to metastasis between the metastatic and non-metastatic samples for the prognosis of osteosarcoma patients [86]. To our knowledge, we are the first to construct a prognostic model of osteosarcoma with autophagy related genes via integrated methods of bioinformatics and machine learning.

Except for risk score, metastasis was also an independent factor for osteosarcoma prognosis. We applied risk score and metastasis as independent prognostic factors to construct a nomogram model for OS prediction in osteosarcoma patients. The result indicated that the nomogram model with two independent factors showed better performance in OS prediction than that with one factor. Furthermore, three immune cell types presented obvious distinction in infiltration proportion between samples with high and low risk scores. The infiltration proportions of activated dendritic cells, M2 macrophages and CD8 T cells were significantly higher in the low risk group than those in the high risk group. Zhang et al. indicated that the osteosarcoma patients with superior survival outcome had higher levels of M2 macrophages, compared with those with inferior survival outcome [87], which was consistent with our finding that the low-risk osteosarcoma patients with improved prognosis had higher proportion of M2 macrophage. Gomez-Brouchet et al. found that CD8 T cells were related to the non-metastatic osteosarcoma [88], and higher infiltration rate of CD8 T cells indicated improved survival outcome [89], showing consistency with our result that the low-risk osteosarcoma patients with superior prognosis had higher infiltration proportion of CD8 T cells. In terms of activated dendritic cells, Wang et al. demonstrated that activated dendritic cell was an independent predictor of osteosarcoma and the model basing on several immune cell types including activated dendritic cell could reliably predict the prognosis of osteosarcoma patients [90], suggesting the potentially important role of activated dendritic cells in osteosarcoma. However, the specific functions of these significantly distinct infiltrating immune cells between the osteosarcoma patients with high and low risk scores and the underlying mechanism still need further research. In addition, analysis of five immune checkpoints in osteosarcoma patients showed that the expressions of CTLA4, TIM3, LAG3, and TIGIT in the samples with low risk scores were markedly elevated compared with those with high risk scores, suggesting that patients with low risk scores might be sensitive to the treatment targeting immune checkpoints.

However, there are some limitations in our study. First, it lacks experimental work, which would be helpful for further validation of the results. Second, the sample size is relatively small. Further research with more samples is needed to better evaluate the performance of the model and elucidate the underlying mechanism in the future.

In brief, we identified a 14-autophagy-gene-based prognostic signature in osteosarcoma. Based on these 14 genes associated with autophagy, a model was established and risk score was able to predict the prognosis of osteosarcoma patients independently. Importantly, a nomogram model based on risk score and metastasis was established and exhibited better performance to predict the OS at 1, 3, and 5 years for osteosarcoma patients.