What this is
- This research investigates the role of -related long non-coding RNAs (FRLncs) in osteosarcoma (OS).
- It identifies seven FRLncs that serve as prognostic markers and influence the immune microenvironment in OS.
- The study employs transcriptome and clinical data from 86 OS cases to develop a risk prognostic model.
Essence
- Seven FRLncs were identified as potential prognostic markers for osteosarcoma, influencing both patient outcomes and the immune microenvironment. The model developed can stratify patients into high-risk and low-risk categories based on these markers.
Key takeaways
- Three high-risk FRLncs (APTR, AC105914.2, AL139246.5) and four low-risk FRLncs (DSCR8, LOH12CR2, AC027307.2, AC025048.2) were identified. These markers can help predict patient prognosis and guide treatment decisions.
- The high-risk group showed down-regulation in immune functions, including antigen-presenting cell co-stimulation and cytolytic activity. This indicates a potential immunosuppressive environment in patients with worse prognoses.
- Seven immune checkpoint-associated genes (CD200R1, HAVCR2, LGALS9, CD27, LAIR1, LAG3, TNFSF4) were linked to OS prognosis. Their differential expression may provide insights into therapeutic targets.
Caveats
- The study's sample size of 86 cases may limit the generalizability of the findings. Larger studies are needed to confirm these results.
- The identified FRLncs have not undergone experimental validation, which is crucial for establishing their clinical utility as prognostic markers.
Definitions
- Ferroptosis: A regulated form of cell death driven by iron accumulation and lipid peroxidation.
- Long non-coding RNA (lncRNA): A type of RNA that does not encode proteins but plays a role in regulating gene expression.
AI simplified
Introduction
Osteosarcoma (OS) stands as the prevailing malignant neoplasm affecting bone tissue, with a predilection for the adolescent demographic [1]. This condition is typified by frequent vascular infiltration, adjacent soft tissue infiltration, a notable proclivity for local recurrence and premature distant metastasis [2]. Approximately one-fifth of OS patients experience the emergence of metastatic lesions, while the remainder often develop subclinical micrometastases. Standard therapeutic modalities encompass the deployment of chemotherapy and surgical resection [3]. Notwithstanding the comprehensive implementation of a multidisciplinary regimen, which encompasses chemotherapeutic intervention and extensive surgical excision, discernible enhancements in clinical outcomes have been documented in patients with OS. However, in instances of advanced disease with remote metastasis and local recurrence, even with the rigorous administration of chemotherapy, clinical outcomes and 5-year overall survival rates remain suboptimal [4]. A burgeoning corpus of scientific inquiry posits that a multifaceted interplay of cellular and molecular events may underlie the pathogenesis of OS [1]. Long non-coding RNA (lncRNA) constitutes a pivotal regulator in a myriad of biological processes, including cell proliferation, apoptosis, invasion and migration. Anomalous expression of lncRNA assumes a pivotal role in the orchestration of tumoral metastasis [5, 6]. The nexus between lncRNA and OS is undeniably close, as investigations have unearthed compelling evidence. For instance, lncRNA LOC100129620 has been implicated in the promotion of OS progression through the modulation of cyclin dependent kinase 6 (CDK6) expression, tumor angiogenesis and macrophage polarization [5]. Conversely, LncRNA FTX exerts an inhibitory effect on OS proliferation and migration via its regulatory influence on the miR-320a/TXNRD1 axis [6]. The exploration of the interrelationship between lncRNA and OS bears paramount significance in the context of disease prognosis and therapeutic intervention. Consequently, the quest for novel prognostic markers in the realm of OS and the pursuit of strategies to enhance clinical efficacy remain imperatives in the management of this disease.
Ferroptosis constitutes a tightly regulated mode of cellular demise, triggered by perturbations in the intracellular milieu, primarily governed by the activity of glutathione peroxidase 4 (GPX4) [7]. The pivotal drivers of ferroptosis initiation encompass the accrual of ferrous iron (Fe2+) and the subsequent peroxidation of lipids [7]. Notably, ferroptosis is amenable to inhibition via iron-chelating agents and lipophilic antioxidants [7]. This multifaceted cellular event intricately interweaves processes inclusive of iron homeostasis, lipid metabolism, oxidative stress, as well as the synthesis of essential molecules such as nicotinamide adenine dinucleotide phosphate (NADPH), glutathione (GSH) and coenzyme Q10 (CoQ10) [8]. Accumulating investigations have underscored the pivotal role of ferroptosis in various malignancies. In the context of breast cancer cells, both a lysosome-disrupting agent, silamesine and a tyrosine kinase inhibitor, lapatinib, have been demonstrated to induce ferroptosis [9]. Furthermore, in the realm of pancreatic ductal adenocarcinoma (PDAC), the agent artesunate (ART) exerts its ferroptosis effects through the generation of reactive oxygen species (ROS) in an iron-dependent manner [10]. These findings collectively emphasize the burgeoning potential of ferroptosis as a therapeutic modality in the domain of oncology, consequently intensifying research efforts toward the design and development of anticancer agents capable of inducing ferroptosis [11]. Moreover, an emerging intersection has been observed between ferroptosis and OS. For instance, tirapazamine has been shown to exert partial inhibition of OS cells through the mediation of solute carrier family 7 member 11 (SLC7A11)-associated ferroptosis [12]. Similarly, EF24 has been identified as an inducer of ferroptosis in OS cells, operating via heme oxygenase 1 (HMOX1)-dependent mechanisms [13]. Furthermore, the identification and characterization of ferroptosis-related genes (FRGs) capable of prognosticating outcomes in OS patients have constituted pivotal advancements in this realm of study [14].

The analysis flowchart of this study
Materials and methods
Data download and arrangement
The transcriptome data pertaining to 86 cases of OS and accompanying clinical information were obtained from The Cancer Genome Atlas (TCGA) database, accessible at https://portal.gdc.cancer.gov/↗. The OS transcriptome dataset encompassed both messenger RNA (mRNA) and lncRNA data. The clinical dataset includes futime, fustat, gender, age at diagnosis in days, metastatic status (metastatic/non-metastatic), primary tumor site and specific tumor site. Additionally, the GSE19276↗, GSE16088↗ and GSE33382↗ datasets were procured from the Gene Expression Omnibus (GEO) database, available at https://www.ncbi.nlm.nih.gov/geo/↗. The GSE19276↗ dataset encompasses five normal tissue samples and 23 OS tissue samples. The GSE16088↗ dataset comprised six normal tissue samples and 14 OS tissue samples, while the GSE33382↗ dataset included three normal tissue samples and 84 OS tissue samples. Furthermore, a comprehensive list of FRGs was sourced from the FerrDb database, accessible at http://www.zhounan.org/ferrdb/↗.
OS-related differentially FRGs
The list of FRGs was intersected with the gene set from the OS transcriptome data, yielding a set of FRGs may be related to OS. The original study of GSE19276↗ dataset obtained 205 differentially expressed genes (DEGs) through performed differential analysis between five normal tissue samples and 23 OS tissue samples. The screening criteria for difference analysis were as follows: P < 0.05 and |logFC|≥ 1 [17]. Differential analysis of the GSE16088↗ and GSE33382↗ datasets was conducted utilizing the “limma” package in the R programming environment. DEGs in OS were identified based on stringent criteria, with a significance threshold of P < 0.05 and |logFC|≥ 1 [15]. The FRGs may be related to OS, DEGs from GSE19276↗ dataset, and DEGs of GSE16088↗ and GSE33382↗ datasets were subsequently intersected, resulting in the identification of OS-related differentially expressed FRGs.
OS-related FRLncs
The co-expression analysis of OS-related differentially expressed FRGs and lncRNAs within the OS transcriptome dataset was conducted utilizing the “limma” package in the R. The objective of this analysis was to identify co-expression lncRNAs of OS-related differentially expressed FRGs, herein referred to as OS-related FRLncs. The criteria for screening were established as follows: |Person correlation coefficient|> 0.4, P < 0.001, as indicated in a previous study [18]. Cytoscape [19] visualizes the one-to-one correspondence between OS-related differentially expressed FRGs and co-expression lncRNAs.
Construction of risk prognostic model
Perl integrated expression matrices of OS-related FRLncs with OS survival data. The “survival” package in the R was employed to identify statistically significant FRLncs that are associated with OS prognosis, utilizing Univariate Cox regression analysis. To mitigate the potential issue of overfitting, the “glmnet” package in R was utilized to conduct Lasso Cox regression analysis. Construct a risk prognosis model for OS prognosis FRLncs. The main body of the risk prognosis model is the optimal number of FRLncs obtained by Lasso Cox regression analysis. A riskscore for each sample was calculated based on these FRLncs. The risk prognosis model divides the sample into two groups: high-risk and low-risk groups. The essence of the risk prognosis model is to compare the survival differences of patients with OS between high-risk and low-risk groups. Calculate the riskscore: Riskscore = \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\sum }_{i=1}^{n}\left(lncrnaex{p}_{i}\times coe{f}_{i}\right)$$\end{document}∑i=1nlncrnaexpi×coefi, n represents the number of OS prognostic FRLncs, i is the ith FRLncs, and coef is the regression coefficient. Each OS prognostic FRLncs's expression level is multiplied by its corresponding regression coefficient, and then accumulated to obtain the sample riskscore. Subsequently, the samples are dichotomized into high-risk and low-risk groups based on the median value of the sample riskscore [20].
Risk curves and survival analysis
R further generates survival status maps and risk heatmaps to assess differences in survival time and OS prognosis FRLncs among high- and low-risk groups defined by the risk prognostic model [21]. The creation of survival curves, aimed at evaluating potential survival disparities between these high- and low-risk groups within the risk prognostic model, was facilitated using the “survival” and “survminer” packages in R.
Receiver operating characteristic curve (ROC) analysis and independent prognostic analysis
The ROC curves were generated employing the “survival,” “survminer” and “timeROC” packages within the R. Utilizing the “survival” package in R, we conducted independent prognostic assessments via both univariate and multivariate COX regression analyses. The aim was to investigate the viability of the riskscore within the prognostic model as an autonomous prognostic determinant [22].
Model validation for clinical grouping
To evaluate the applicability of the risk prognostic model across diverse clinical patient cohorts, we conducted a comprehensive validation procedure involving the alignment of clinical characteristics with the risk prognostic model. The clinical attributes were stratified as follows: gender dichotomized into male and female categories; age categorized into two groups, namely, individuals aged ≤ 5245(14 years old) and > 5245(14 years old); metastasis status categorized as either metastatic or non-metastatic; primary tumor site classified as upper limb or lower limb + pelvis; and specific tumor site categorized as upper limb or lower limb + pelvis. Subsequently, Perl scripting was employed to merge the categorized clinical data with the corresponding riskscores. To assess the model’s suitability for each specific clinical trait, we leveraged the “survival” and “survminer” packages within the R, conducting a model validation process for each of the clinical subgroups.
Single sample gene set enrichment analysis (ssGSEA)
The enrichment scores for immune cells and immune functions in the OS transcriptome data were computed utilizing the R packages “GSVA,” “limma” and “GSEABase.” Subsequently, an analysis of differences between high- and low-risk groups in terms of immune cell composition and immune function was performed using the R packages “limma,” “reshape2” and “ggpubr”[20].
Differential analysis of immune checkpoints
The analysis of correlation among immune checkpoint-related genes and risk prognostic model was conducted using the R packages “limma,” “reshape2,” “ggplot2” and “ggpubr.” The objective was to assess the difference in immune checkpoint-related genes between high- and low-risk groups [22].
Validation of the risk prognosis model
Importantly, we sought to validate the constructed risk prognosis model. To do so, we evenly partitioned the samples into two distinct groups: a train group and a test group. This partitioning allowed us to assess the model's accuracy effectively. Subsequently, we calculated the riskscore for each sample within the train group, employing the previously described formula. The samples in the train group were then classified into high-risk and low-risk groups based on the median value of riskscore. For the categorization of samples within the test group, we employed a similar approach, dividing them into high-risk and low-risk groups based on the median value of riskscores calculated within the train group. To perform survival analysis, we utilized the “survival” and “survminer” packages of R. Subsequently, we employed the R to generate a risk heatmap.
Results
OS-related differentially FRGs
A total of 382 FRGs were extracted from FerrDb (Additional file 1: Table 1). The transcriptome data of OS from the TCGA database contained 19,262 genes. The intersection of 382 FRGs and 19,262 genes that yielded 239 FRGs may be related to OS. The original study of the GSE19276↗ dataset provided 205 DEGs (Additional file 1: Table 2).

OS-related differentially expressed FRGs and OS-related FRLncs.The intersection of DEGs fromdataset and FRGs may be related to OS that obtained OS-related differentially expressed FRGs.Volcanomap ofand GSE3338 datasets difference analysis, five OS-related differentially expressed FRGs overlap with DEGs ofand GSE3338 datasets.By co-expression analysis, OS-related differentially expressed FRGs yielded a total of 48 OS-related FRLncs. The screening criteria are: |Person correlation coefficient|> 0.4,< 0.001 a b c GSE19276 GSE16088 GSE16088 P
OS-related FRLncs
Subsequently, we analyzed the co-expression of five OS-related differentially expressed FRGs and lncRNAs contained in the OS transcriptome data. Through co-expression analysis, a subset consisting of three out of the five OS-related differentially expressed FRGs produced a combined total of 48 OS-related FRLncs (Fig. 2c).
Construction of risk prognostic model

Construction of risk prognostic model.Univariate Cox regression analysis obtained eight candidate prognostic FRLncs for OS, including four high-risk FRLncs and four low-risk FRLncs.LASSO Cox regression analysis for eight candidate prognostic FRLncs.Selection of the optimal penalty parameter for LASSO regression a b c
Risk curves and survival analysis

Comprehensive analysis of risk prognosis model.Survival status plot, from the low-risk group to the high-risk group, the survival time of patients with OS decreased.Risk heatmap. APTR,andare high-risk FRLncs. DSCR8, LOH12CR2,andare low-risk FRLncsSurvival analysis, there was a significant survival difference of patients with OS between the high- and low-risk groups.ROC curves, the risk prognosis model could well predict the survival rate of patients with OS at one, three and 5 years a b c d AC105914.2 AL139246.5 AC027307.2 AC025048.2
ROC curve analysis and independent prognostic analysis

Independent prognostic analysis and model validation for clinical grouping.Univariate independent prognostic analysis, riskscore and metastasis were independent prognostic factors.Multivariate independent prognosis analysis, riskscore and metastasis were independent prognostic factors.Model validation for clinical grouping, the risk prognostic model across various stratifications, encompassing age, gender and tumor metastasis a b c
Model validation for clinical grouping
The clinical traits were grouped to observe whether the risk prognostic model was suitable for patients in different clinical groups. Model validation for clinical categorizations demonstrated the appropriateness of the risk prognostic model across various stratifications, encompassing age, gender and tumor metastasis (Fig. 5c).
Single sample gene set enrichment analysis (ssGSEA)

Single sample gene set enrichment analysis and differential analysis of immune checkpoint.Differential analysis of immune cells in high- and low-risk groups in risk prognostic model.Differential analysis of immune function in high- and low-risk groups in risk prognostic model.Differential analysis of immune checkpoint-related genes in high- and low-risk groups in risk prognostic model. * means< 0.05, ** means< 0.01, *** means< 0.001 a b c P P P
Differential analysis of immune checkpoints
Immune checkpoint differential analysis revealed disparities in the expression levels of 22 immune checkpoint-associated genes between high-risk and low-risk groups. Notably, among these genes, hypocretin receptor 2 (CD200R1), hepatitis A virus cellular receptor 2 (HAVCR2), galectin-9 (LGALS9), cluster of differentiation 27 (CD27), leukocyte associated immunoglobulin like receptor 1 (LAIR1), lymphocyte activating 3 (LAG3) and TNF superfamily member 4 (TNFSF4) displayed a strikingly high degree of statistical significance (P < 0.001) (Fig. 6c).
Validation of the risk prognosis model

Validation of the risk prognosis model.Survival analysis of train group.Risk heatmap of train group.Survival analysis of test group.Risk heatmap of test group a b c d
Discussion
We developed a risk prognostic model for assessing the OS-related FRLncs based on a comprehensive analysis of data from both TCGA and GEO databases. This model incorporates information from seven FRLncs classified into two categories: three high-risk FRLnc and four low-risk FRLncs. The risk prognosis model constructed in this study can well predict the survival of patients with OS, and has a certain impact on the immune microenvironment of patients with OS. Compared with the low-risk group, the immune cells and immune function of patients with OS in the high-risk group showed a downward trend. Immune checkpoint-related genes also differed between the two group.
Previous study has shown that LncRNA APTR participates in the progression of OS by repression of miR-132-3p and upregulation of Yes-associated protein 1 (YAP1) [23]. The expression of APTR in OS tumor tissues and four OS cell lines (MG63, 143B, Saos-2 and HOS) was significantly up-regulated compared with that of in neighboring tissues and human osteoblast cell lines hFOB1.19, respectively. MiR-132-3p is the target of APTR, and its expression is inhibited by APTR. Both knockdown of APTR and overexpression of miR-132-3p can significantly inhibit the proliferation, invasion and migration of human OS cells, and induce cell apoptosis. In addition, YAP1 was identified as a target for the inhibition of miR-132-3p [23]. LncRNA DSCR8 promotes the proliferation of liver cancer cells and inhibits apoptosis via the miR-22-3p/ARPC5 axis [24]. LncRNA DSCR8 mediates miR-137/Cdc42 to regulate gastric cancer cell proliferation, invasion and cell cycle as a competitive endogenous RNA [25]. The study found that smoking promoted the development of lung adenocarcinoma and lung squamous cell carcinoma by affecting gender-specific lncRNA changes. In women with lung squamous cell carcinoma, changes in LOH12CR2 were positively associated with smoking index [26]. AC027307.2↗ can not only be used as an enhancer-associated lncRNA to become a specific prognostic marker for breast cancer, but also as an immune-associated lncRNA to predict the survival rate of patients with colorectal adenocarcinoma [27, 28]. Notably, of the seven FRLncs that can guide OS prognosis and immunity, APTR, DSCR8, LOH12CR2 and AC027307.2↗ were all associated with tumor progression and prognosis. Among them, APTR has been found to be involved in the development of OS in previous studies. These findings underscore the robustness of our results. In the existing literature, AC105914.2↗, AL139246.5↗ and AC025048.2↗ remain to be elucidated and documented.
The immune system was one of the major components in the tumor microenvironment and was often suppressed in hypoxia [29]. The study found that the pre-treatment neutrophil-to-lymphocyte ratio was of diagnostic and prognostic value for OS [30]. Clinical data suggest that NK cells may play an important role in the prevention and therapeutic response of OS. In patients with OS, the number of circulating NK cells in peripheral blood was lower than in normal controls, suggesting that NK cells play a preventive role in the development of OS tumors [31]. pDCs exhibit remarkable proficiency in the rapid and robust production of type I interferon (IFN-I/α) [32]. Nevertheless, in the context of cancer, pDCs demonstrate a diminished responsiveness to Toll-like receptor 7 and 9 (TLR7/9) activation, resulting in a marked reduction or outright loss of IFN-α production. This diminished IFN-α contributes to the establishment of an immunosuppressive tumor microenvironment [33]. Notably, pDCs play a pivotal role in the regulation of both innate and adaptive immune systems, rendering them essential actors in the realm of cancer immunity [33]. Investigations focusing on the interplay between immune cell populations and the OS microenvironment have elucidated that cancer progression is most expedited when anti-tumor immune cells, including DCs, helper T cells, cytotoxic cells and IFN-γ, exhibit a decline in abundance, while regulatory T cells (Treg) undergo an increase [34]. Several studies have posited that the secretion of cytokines and the proliferative capacity of T follicular helper (Tfh) cells are markedly impaired in OS patients. This diminishment in functionality may underlie the body's compromised ability to resist OS development [35]. Importantly, miR-138 serves as an inhibitor of the PI3K/Akt/mTOR pathway by specifically targeting and negatively regulating pyruvate dehydrogenase kinase 1 (PDK1), thereby mitigating Tfh dysfunction in OS [36]. The TIL tasked with combating tumor cells within the OS microenvironment experience depletion, thus hastening tumor recurrence [37]. Notably, the incorporation of adjuvant chemotherapy in conjunction with TIL therapy has demonstrated the capacity to extend survival in OS patients with a suboptimal response to neoadjuvant chemotherapy [37]. Within the cancer microenvironment, innate immunity is systematically suppressed in the context of immune tolerance. A pivotal mechanism underlying immune tolerance is the immune checkpoint mechanism, which operates to restrain T cell activity in order to prevent undue immune responses [38].
The high-risk group exhibited significant down-regulation of immune functions, including APC-co-stimulation, immune checkpoint regulation, cytolytic activity and T cell co-inhibition. Study revealed that the administration of the anti-epidermal growth factor receptor (EGFR) antibody, cetuximab, resulted in an augmentation of the cytolytic activity exhibited by NK cells in the context of OS [39]. Immune checkpoint inhibitors (ICIs) have been recognized for their capacity to reshape the course of various malignancies by eliciting a disruption in immune homeostasis, thereby empowering the host's immune system to combat the tumor [40]. In a comprehensive analysis focusing on hypoxic prognostic indicators associated with OS metastasis and immune cell infiltration, a pronounced down-regulation of immune checkpoint mechanisms was identified among high-risk populations [29]. Regrettably, the precise correlation between APC-co-stimulation and T cell co-inhibition with respect to OS remains an area of investigation that warrants further elucidation in the scientific literature.
CD200R1 is not only differentially expressed in non-small cell lung cancer and has a prognostic effect, but also predicts survival in patients with head and neck squamous cell carcinoma [41, 42]. Previous study has found that HAVCR2 can be used as immune signature to accurately predict the prognosis of patients with OS, high expression of HAVCR2 is associated with improved prognosis [43, 44]. A comprehensive investigation has revealed significant alterations in Galectin-9, encoded by the LGALS9, across various cancer types [45]. Remarkably, LGALS9 not only exhibits associations with mRNA expression levels in cervical cancer cells but also emerges as a potential prognostic biomarker in pancreatic cancer [46]. The CD70-CD27 interaction is important for the regulation of adaptive immunity. Phospholipase C epsilon 1(PLCE1) is a marker of poor prognosis and may promote immune escape of OS through the CD70-CD27 signaling pathway [47]. LAIR1 overexpression inhibits the epithelial–mesenchymal transformation of OS through glucose transporter (Glut) 1-related energy metabolism [48]. The construction and validation of an oxidative stress-related prognostic risk model for OS suggests that LAG3 is a potential immunotherapeutic target for patients with OS [49]. In all melanoma patients and in the stage III and IIIc-IV patient cohorts, low expression of TNFSF4 was associated with a poorer prognosis. In the subgroup of patients with low lymphocytic infiltration, low expression of TNFSF4 was also associated with a poorer prognosis [50]. It is suggested that TNFSF4 is associated with tumor prognosis.
Nevertheless, it is imperative to acknowledge the presence of several limitations within the scope of this investigation. Firstly, the study's sample size remained relatively small, and the sampling methodology employed did not successfully mitigate the potential confounding influence of gender and underlying medical conditions. Secondly, the FRLncs identified in this study, which bear the potential to prognostic outcomes in the context of OS and contribute insights into the immune microenvironment, have not yet undergone experimentally validation. However, it is noteworthy that this validation process constitutes a focal point for our forthcoming research endeavors.
Conclusion
In this study, a prognostic risk model for OS was formulated by integrating data from the TCGA and GEO databases. Concurrently, immune analysis was conducted, yielding seven FRLncs identified as potential prognostic markers for OS and influencers of the immune microenvironment. Differential profiles of immune cells and functional characteristics were delineated within the context of the risk prognostic model. Furthermore, we identified seven immune checkpoint-associated genes with notable implications for OS prognosis. The discovery of these FRLncs, their pivotal roles in OS prognosis, the modulation of the immune microenvironment, as well as the identification of immune checkpoint-related genes, collectively furnish a solid theoretical foundation for advancing research in OS survival and facilitating informed clinical decision making.
Supplementary Information
Additional file 1. Table 1: ferroptosis-related genes. Table 2: differentially expressed genes of GSE19276 dataset. Table 3: differentially expressed genes of GSE16088 and GSE33382 datasets.