Construction of an Immune-Related lncRNA Signature That Predicts Prognosis and Immune Microenvironment in Osteosarcoma Patients

Osteosarcoma is one of the most common bone tumors in teenagers. We hope to provide a reliable method to predict the prognosis of osteosarcoma and find potential targets for early diagnosis and precise treatment. To address this issue, we performed a detailed bioinformatics analysis based on the Cancer Genome Atlas (TCGA). A total of 85 osteosarcoma patients with gene expression data and clinicopathological features were included in this study, which was considered the entire set. They were randomly divided into a train set and a test set. We identified six lncRNAs (ELFN1-AS1, LINC00837, OLMALINC, AL669970.3, AC005332.4 and AC023157.3), and constructed a signature that exhibited good predictive ability of patient survival and metastasis. What’s more, we found that risk score calculated by the signature was positively correlated to tumor purity, CD4⁺ naive T cells, and negatively correlated to CD8⁺ T cells. Furthermore, we investigated each lncRNA in the signature and found that these six lncRNAs were associated with tumorigenesis and immune cells in the tumor microenvironment. In conclusion, we constructed and validated a signature, which had good performance in the prediction of survival, metastasis and immune microenvironment. Our study indicated possible mechanisms of these lncRNAs in the development of osteosarcoma, which may provide new insights into the precise treatment of osteosarcoma.

Introduction

Osteosarcoma is one of the most common bone tumors, most commonly occurring in young children and adolescents (1). This tumor is most likely to happen in the metaphyses of the distal femur, proximal tibia, and proximal humerus (2–4). Osteosarcoma is prone to pulmonary metastases, and 20% of patients are found to have pulmonary metastases at the time of initial diagnosis (5). Medical advances have significantly reduced the mortality rate of osteosarcoma patients. However, the lack of specific markers makes early screening for osteosarcoma still difficult. Treatment of osteosarcoma is mainly based on local excision and chemotherapy, but chemotherapy for osteosarcoma is prone to drug resistance and has great toxic side effects (6–8). As a result, overall survival is still not satisfactory. It is urgent to figure out the mechanism of tumorigenesis, metastasis and drug resistance in osteosarcoma, and to discover potential target for earlier diagnosis and gene therapy.

The immune system is an important part of the human body. It can help us fight against pathogen infections and participate in the monitoring and prevention of cancer, playing an essential role in maintaining the integrity of the body. However, some tumor cells can evade the surveillance by immune system, or suppress the immune response, making cancer progression. Meanwhile, more and more evidence has shown that the imbalance of immune state in tumor microenvironment plays a decisive role in tumor development (9). Therefore, the role of immune-related factors in tumor development deserves to be studied.

LncRNAs are highly heterogeneous RNA characterized by their length of more than 200 nucleotides and do not encode proteins. With advances in gene chip technology, lncRNAs are being rapidly identified. Numerous studies have found that lncRNAs are involved in various physiological processes, including cellular differentiation, immune response, and tumor progression (10, 11). LncRNAs appear to play an important role in tumor progression, exhibiting tumor-inhabiting and tumor-promoting functions (12). LncRNA-ATB was reported to functions as a tumor promoter in papillary thyroid cancer (13). LOC285194 was reported to function as a tumor suppressor in non-small cell lung cancer and osteosarcoma (14, 15). As a result, lncRNAs are expected to be new biomarkers and therapeutic targets for cancer (16).

There was little research studied on the immune-related lncRNAs in osteosarcoma. In this study, we screened immune-related lncRNAs in osteosarcoma and made a bioinformatics analysis for them on the basis of TCGA data. The flowchart of our study is as follows (Figure 1).

FIGURE 1

Figure 1 Flowchart of our study.

Materials and Methods

Sample Datasets

The RNAseq (level 3) data and corresponding clinical data of osteosarcoma were come from the TCGA database (https://portal.gdc.cancer.gov/), and all data was normalized by TMM method. We used the Ensembl Genome Browser website to annotate the mRNAs and lncRNAs in the expression matrix and subsequently extracted the mRNA and lncRNA expression matrix separately. To make the analysis accurate, we removed the samples with incomplete clinical information and overall survival time lower than 30 days (17). In total, 85 cases were contained for the following study. All TCGA data is available to the public, so there is no further approval needed from the Ethics Committee.

Immune-related lncRNA Acquisition

We acquired immune-related genes from the Molecular Signature Database v 7.1 (MSigDB) (http://www.broad.mit.edu/gsea/msigdb/) (17, 18). IMMUNE_RESPONSE.gmt and IMMUNE_SYSTEM_PROCESS.gmt were chosen as the annotated gene sets. We removed overlapped genes and ultimately obtained 331 immune-related genes. The expression information of these 331 immune-related genes was extracted from the mRNA expression matrix. Then Pearson correlation analysis was performed to calculate the correlation coefficient of each lncRNA with immune-related mRNAs. LncRNAs with coefficient > 0.6 and P value < 0.001 were defined as immune-related lncRNAs. Finally, the expression matrix of these immune-related lncRNAs was extracted.

Signature Construction

These 85 cases were separated into a train set and a test set randomly by “caret” package of R software (v 3.6.2). Meanwhile, all 85 cases were taken as the entire set. The train set was used to construct the prognostic signature. The test set and the entire set were used to verify the accuracy of the signature.

In the train set, we performed univariate Cox regression analysis to identify prognostic-related lncRNAs with p<0.01. Then Least absolute shrinkage and selection operator (LASSO) regression was performed to remove lncRNAs that could lead to phenomenon of overfitting. Finally, the lncRNAs was further screened by using multivariate Cox regression analysis. The coefficient, HR value and P value for each lncRNA were calculated. The risk score of each sample was calculated by the signature:

$$Risk score = \beta \left(gene1\right) \ast expr\left(gene1\right) + \beta \left(gene2\right)\ast expr\left(gene2\right) +···+ \beta \left(genen\right) \ast expr(genen)$$

β is the coefficient of lncRNA, and expr(genen) stands for the expression value of lncRNA.

Signature Application and Validation

Firstly, we did the analysis in the train set. We used the formula above to calculate the risk score of each patient and divided patients into high or low risk groups according to the risk score by using “survminer” package. We ranked the patients by the risk score and plotted the dot-plot for the survival status of each patient. We mapped the lncRNA expression heatmap to observe the expression of lncRNAs in the high and low risk groups. We applied Kaplan-Meier method to explore the survival differences between high-risk group and low-risk group. Finally, receiver operating characteristic (ROC) curve of 5-years was applied for identifying the diagnostic value of the risk scoring signature.

Subsequently, the same analysis was performed in the test set and the entire set to verify the accuracy of the signature.

Comprehensive Analysis in the Entire Set

The entire set was used for the subsequent analysis. Firstly, principal component analysis (PCA) was used to test whether the signature can better distinguish the risk status. Secondly, survival analysis for each lncRNA in the signature was performed to estimate the effect of individual lncRNA on survival. Thirdly, metastasis correlation analysis was carried out to investigate the potential correlation between risk score and metastasis. Fourthly, Gene set enrichment analysis (GSEA; https://www.gsea-msigdb.org/gsea/index.jsp) (17) was performed based on two gene sets (IMMUNE_RESPONSE.gmt and IMMUNE_SYSTEM_PROCESS.gmt) to analyze immune response enrichment in the two groups. Fifthly, ESTIMATE algorithm (19) was used to assess the proportion of immune and stromal components, and to infer the tumor purity of the samples. We subsequently analyzed the relationship between tumor purity and risk scores. Finally, CIBERSORT method (20) was used to analyze the infiltration of 22 immune cells between high and low risk groups. The correlation between the infiltrating proportion of immune cells and risk scores and the lncRNAs were also analyzed.

All analyses were performed on the R software (version 3.6.2, https://www.r-project.org/). And P < 0.05 was considered statistically significant.

Results

Construction of Signature

By comparing mRNA expression data with two gene sets from MSigDB, we matched 331 immune-related genes. We calculated the expression correlation between lncRNAs and immune-related genes to obtain the immune-related lncRNAs. Subsequently, 423 immune-related lncRNAs was identified. Eight of them were screened out by univariate Cox regression analysis (Table 1). The result of the LASSO regression analysis showed that the partial likelihood deviation was smallest when -3 < lambda < -2 (Figure 2A), at which point AL137002.1 was excluded (Figure 2B). Finally, six lncRNAs were screened from the remaining seven lncRNAs by using multivariate Cox regression analysis, and their coefficients, HR values and p-values were calculated respectively (Table 2). The signature was constructed with the coefficients of lncRNAs as below.

$$Riskscore=(0.4707\ast ELFN1-AS1)+(0.1634\ast LINC00837)+(0.4341\ast AL669970.3)+(0.2696\ast OLMALINC)+(-1.1261\ast AC005332.4)+(-2.2230\ast AC023157.3).$$

TABLE 1

Table 1 Eight lncRNAs obtained after univariable Cox regression analysis.

FIGURE 2

Figure 2 The LASSO regression analysis. When -3 < lambda < -2, the partial likelihood deviation was smallest (A), at which point AL137002.1 was excluded (B).

TABLE 2

Table 2 Six lncRNAs obtained after multivariable Cox regression analysis.

Among these lncRNAs, five immune-related lncRNAs were independent prognostic factors. One immune-related lncRNAs acted as a complement to others.

Validation of the Signature

We calculated the risk score of each sample by the signature and divided patients into the low-risk and the high-risk groups by the median value of risk score. We ranked patient by the risk score (Figure 3A) and plotted the dot-plot for the survival status of each patient (Figure 3B). The results showed that the higher risk score patients had, the shorter survival time patients might have. According to the heatmap, the expression levels of ELFN1-AS1, LINC00837, OLMALINC and AL669970.3 were higher in the high-risk group, while AC005332.4 and AC023157.3 were higher in the low-risk group (Figure 3C). The survival curve shows that the patients in the high-risk group had a lower survival time (Figure 3D). The area under curve (AUC) of 5 years in the train set was 0.937 (Figure 3E).

FIGURE 3

Figure 3 Construction of immune-related lncRNA signature for osteosarcoma. (A) The risk curve of each patient reordered by risk score in train set. (B) The scatter plot of all patient’s survival state in train set. (C) The heatmap showed the expression levels of six lncRNAs between the low-risk group and high-risk group in train set. (D) Patients in the high-risk group indicated worse overall survival than those in the low-risk group in train set. (E) AUC for risk score of 5-year survival according to the ROC curves in train set.

The same signature was used in the test set and the entire set. And we got similar results in the two sets as we expected. Patients in the low-risk group have better overall survival. The results in the test set (Figures 4A, C, E, G) and the entire set (Figures 4B, D, F, H) were presented in the figure below. The area under curve (AUC) of 5 years in the test set was 0.797, and 0.879 in the entire set (Figures 4I, J).

FIGURE 4

Figure 4 Validation of immune-related lncRNA signature for osteosarcoma. (A, B) The risk curve of each patient reordered by risk score in test set and entire set. (C, D) The scatter plot of all patient’s survival state in test set and entire set. (E, F) Patients in the high-risk group indicated worse overall survival than those in the low-risk group in test set and entire set. (G, H) The heatmap showed the expression levels of six lncRNAs between the low-risk group and high-risk group in test set and entire set. (I, J) AUC for risk score of 5-year survival according to the ROC curves in test set and entire set.

In order to test whether the signature can better distinguish the risk status, PCA analysis was carried out using the signature and genome-wide expression. When using the signature, the risk status of the patients was separated well (Figure 5A). While it did not display a clear separation when using the whole genome expression (Figure 5B).

FIGURE 5

Figure 5 PCA analysis showed that the use of this signature (A) could better distinguish the risk status of patients than the use of the whole genome (B).

Survival Analysis for lncRNAs in the Signature

In order to estimate the effect of expression of individual lncRNAs in the signature, we performed survival analysis for each lncRNA in the signature. Patients with higher expression of AC005332.4 and AC023157.3 indicated better overall survival than those with lower expression (Figures 6A, B). However, the results for other lncRNAs were not statistically significant (Figures 6C–F).

FIGURE 6

Figure 6 The survival curves for each of the six lncRNAs. (A) AC005332.4. (B) AC023157.3. (C) ELFN1-AS1. (D) LINC00837. (E) OLMALINC. (F) AL669970.3.

Metastasis Correlation Analysis

The prognosis of osteosarcoma is strongly related to the presence of metastasis, so it is important to verify whether the signature is predictive of metastasis. We carried out metastasis correlation analysis to investigate the potential correlation between risk score and metastasis in the entire set. The result showed that the risk score was correlated to metastatic status (P=0.0049) (Figure 7A). And one of the lncRNA, AC023157.3, was down-regulated in metastatic patients (P<0.05) (Figure 7B).

FIGURE 7

Figure 7 Metastasis correlation analysis in the entire sets. (A) The risk score was correlated with metastatic (P=0.0049) (Some outliers are not shown in the figure). (B) AC023157.3 was highly correlated with metastasis (P<0.05). *: P<0.05, ns: P>0.05.

Analysis of Tumor Environment and Immune Infiltration

Further analysis by GSEA showed that immune response pathways were enriched in the low-risk groups (Figures 8A, B). The immune score calculated by ESTIMATE algorithm suggested that immune component was significantly higher in the low-risk group (P<0.01) (Figure 9A). The stromal score was also higher in the low-risk group, but it did not achieve statistical significance (P>0.05) (Figure 9B). The ESTIMATE score was higher in the low-risk group, suggesting the lower tumor purity in low-risk group (Figure 9C). We also found that the risk score was negatively correlated with immune score, stromal score and ESTIMATE score by correlation analysis (Figures 9D–F). We calculated the proportion of 22 immune cell in each sample using the CIBERSORT algorithm. The results showed the main components of the immune environment were T lymphocytes and macrophages (Figure 10A). We then group the samples according to the risk scores, and found that the high-risk group was enriched with CD4⁺ naive T cells (P=0.043), while the low-risk group was enriched with CD8⁺ T cells (p=0.041) and CD4⁺ activated memory T cells (p=0.033) (Figure 10B). In addition, we also analyzed the correlation between risk scores and 22 immune cells, and found risk scores was negatively correlated to CD8⁺ T cells and positively correlated to CD4⁺ naive T cells (Figures 10C, D). We then also analyzed the relationship between each gene in the signature and the proportion of immune cells. ELFN1-AS1 was negatively related to plasma cells (Figure 11A). LINC00837 was positively related to resting dendritic cells (Figure 11B). AL669970.3 was positively related to activated CD4⁺ memory T cells (Figure 11C). AC023157.3 was positively related to monocytes and M2 macrophages, and negatively related to memory B cells and M0 macrophages (Figures 11D–G). OLMALINC was positively related to activated mast cells, CD4⁺ naive T cells and negatively related to memory B cells (Figures 10H–J). AC005332.4 did not show association with immune cells.

FIGURE 8

Figure 8 GSEA analysis showed that immune-related responses were enriched in the low-risk groups. (A) GSEA analysis based on IMMUNE_RESPONSE gene set. (B) GSEA analysis based on IMMUNE_SYSTEM_PROCESS gene set. The green line was the running enrichment score (ES) for the gene set as the analysis walks down the ranked list. A positive value of ES indicates that the gene set was enriched in the high risk group, and a negative value of ES indicates that the gene set was enriched in the low risk group. The ranking metric score indicated a gene’s correlation with a phenotype. A positive value refers to high risk and a negative value refers to low risk. NES, normalized enrichment score; FDR, false discovery rates.

FIGURE 9

Figure 9 Tumor microenvironment analysis showed risk score was negatively correlated with tumor purity. (A) Immune score between high-risk group and low-risk group. (B) Stromal score between high-risk group and low-risk group. (C) ESTIMATE score between high-risk group and low-risk group. (D–F) Immune score, stromal score and ESTIMATE score were negatively correlated with risk score.

FIGURE 10

Figure 10 Immune infiltration analysis showed the relationship between risk score and immune cells. (A) The main components of the immune environment were T-lymphocytes and macrophages in osteosarcoma. (B) The composition of immune cells between high-risk group and low-risk group. (C, D) Risk score was negatively correlated to CD8⁺ T cells and positively correlated to CD4⁺ naive T cells.

FIGURE 11

Figure 11 The relationship of each gene in the signature to the immune cells. (A) ELFN1-AS1. (B) LINC00837. (C) AL669970.3. (D–G) AC023157.3. (H–J) OLMALINC.

Discussion

Although the application of new treatment has increased overall survival in osteosarcoma, the prognosis for patients, especially the metastatic and recurrent patients, remains poor due to the lack of specific biomarkers and therapeutic target for early diagnosis and precise treatment. Therefore, it is urgent to find specific biomarkers and therapeutic target for osteosarcoma. Studies have emphasized that the immune response in the microenvironment plays a crucial role in the development of a variety of cancers, and lncRNAs are an important regulator of the immune response (21, 22). Many studies have also reported that lncRNAs were involved in the proliferation, metastasis and drug resistance of osteosarcoma (14, 23–26). There are also emerging evidence indicating that immune-related lncRNAs are valuable in predicting prognosis and also maybe targets for specific treatment (27, 28).

Based on a TCGA dataset, we included 85 samples with complete clinical information in our study. Then we identified a potential prognostic six-lncRNAs signature. This signature included ELFN1-AS1, LINC00837, OLMALINC, AL669970.3, AC005332.4 and AC023157.3. It was of great value in predicting the prognosis of patients in all set as we hope. The 5-years AUC of ROC was 0.879 in the entire set. The prognosis of almost all tumors is highly correlated with metastasis, especially osteosarcoma. We found that the risk score based on the signature was highly correlated with metastasis, suggesting that signature is also a better predictor of osteosarcoma metastasis. It might be one of the reasons that this signature is such a good predictor of patient outcomes. Principal component analysis showed that the patient’s risk status can be well differentiated when using this signature. Moreover, we found that AC005332.4 and AC023157.3 were mainly expressed in the low-risk groups, and they may improve the overall survival of patient. While ELFN1-AS1, LINC00837, OLMALINC and AL669970.3 were mainly expressed in the high-risk groups, and may reduce the overall survival.

Previous research has identifed ELFN1-AS1 was overexpressed in tumors and indicated that it might play an essential role in carcinogenesis (29). Importantly, Down-regulation of ELFN1-AS1 could inhibit the proliferation and migration of tumor cells in esophageal cancer and colorectal cancer. ELFN1-AS1 promoted tumor cell proliferation and metastasis by acting as a sponge of miR-183-3p to upregulate GFPT1 in esophageal cancer (30). ELFN1-AS1 could promote proliferation, metastasis and exert anti-apoptosis effect by up-regulating TRIM44 by sponging miR-4644 in colorectal cancer (31). However, till now the specific function of the lncRNA remains unknown in osteosarcoma. In our study, we found that ELFN1-AS1 was negatively related to plasma cells. Plasma cells have been shown to be positively associated with patient prognosis in colon cancer, non-small cell lung cancer and triple-negative breast cancer (32–34). We speculate that ELFN1-AS1 may influence the development of osteosarcoma by affecting the proliferation and function of plasma cells, and the exact mechanism remains to be further investigated. OLMALINC was reported to overexpress in the white matter of the human brain, which played a role in maintaining the maturation of oligodendrocytes (35). In our study, we found that OLMALINC was positively related to activated mast cells and CD4⁺ naive T cells, and negatively related to memory B cells. Mast cells have been observed to increase in tumor and peritumor tissues (36). However, Mast cells have different roles in different types of tumors. Mast cells were associated with promoting tumorigenesis in bladder cancer, colorectal cancer, lung cancer and hepatocellular cancer (37–40). Mast cells have antitumor activities in breast cancer, diffuse large B-cell lymphoma and ovarian cancer (41–43). Naive T cells were reported to express functional CXCL8 and promote tumorigenesis (44). Blocking naive CD4⁺ T cell recruitment into tumors reversed immunosuppression in breast cancer, which may be an attractive strategy for antitumor treatment (45). Memory B cells were associated with a good prognosis in gastric cancer and non-small cell lung cancer (46, 47). These studies were consistent with our results. OLMALINC may affect osteosarcoma development by inhibiting memory B cell proliferation and function, recruiting mast cells, and naive CD4⁺ T cells.

While, there was little research studied on the role of LINC00837, AL669970.3, AC005332.4 and AC023157.3. In our study, we found that LINC00837 was positively related to resting dendritic cells. Dendritic cells are the most important antigen-presenting cells in the immune system, playing a key role in regulating immunity. However, resting dendritic cells can induce immune tolerance through T cell deletion and induction of regulatory T cells (48). AL669970.3 was positively related to activated T CD4⁺ memory cells. But this is contrary to our common perception that activated T CD4⁺ memory cells can positively modulate immune function and inhibit tumor growth (49). It warrants further investigation. AC023157.3 was positively related to monocytes and M2 macrophages, while negatively related to memory B cells and M0 macrophages. Using a three-dimensional vascularized microfluidic model, Boussommier-Calleja et al. demonstrated that monocytes were able to directly reduce cancer cell extravasation, but that such an effect was lost when monocytes were converted to macrophages (50). Many studies have suggested that M2 polarized tumor-associated-macrophages (TAM) are associated with tumor growth, invasion, and metastasis (51, 52). However, Anne Gomez-Brouchet et al. found that the presence of CD163-positive M2-polarized macrophages was critical to inhibit OS progression (53). This suggests that M2-polarized macrophages may be influenced by other factors to exert distinctly different pro- or anti-tumor effects. Zhang et al. indicated that the level of polarization of M0 to M1 or M2 macrophages may be an important factor (54). It was also reported that the balance between M1 and M2 macrophage would affect the PD-1/PDL-1, which is an important immune regulatory system (55). In addition, we found that among these six RNAs, AC023157.3 was highly expressed in patients without metastasis. We speculate that AC023157.3 may inhibit tumor growth and metastasis by promoting the proliferation and function of monocytes and by regulating the level of polarization of M0 to M1 or M2 macrophages. The specific mechanism deserves further investigation. Combined with the above findings, we found that all the six RNAs in the signature had some relationship with tumorigenesis or the immune system, which could provide a theoretical basis for immunotherapy of osteosarcoma.

We analyzed the enrichment of immune-related pathways in high- and low-risk groups by GSEA analysis, and found that the immune-related pathways were mainly enriched in the low-risk group, suggesting that the development of osteosarcoma may be highly associated with a lack of immune response. ESTIMATE analysis showed that the low-risk group contained more immune cells with a higher immune score and ESTIMATE score, indicating the lower tumor purity. The risk score was also negatively correlated with immune score, stromal score and ESTIMATE score. The results indicated that risk scores were positively correlated with tumor purity. There were studies reported that a lower tumor purity and more immune cell components mean a poor prognosis (19, 54). However, it was different in osteosarcoma that patients with increased microenvironmental immune cell infiltration had a better prognosis (54), which was consistent with our results. It was reported that the main components of the immune environment in osteosarcoma are T-lymphocytes and macrophages (55). CIBERSORT analysis in our study showed similar results. There were more mature lymphocytes, including CD8⁺ T cells (p=0.041) and CD4⁺ activated memory T cells, enriched in the low-risk group. Meanwhile, CD4⁺ naive T cells were enriched in the high-risk group. The risk scores were also negatively correlated to CD8⁺ T cells, while positively correlated to CD4⁺ naive T cells. These results suggested immature immune function in the high-risk group might be an important factor for the development of osteosarcoma. Therefore, targeting the immune system and changing the composition and proportion of immune cells in the tumor microenvironment is a promising treatment option for osteosarcoma.

However, our study has limitations. First, the sample sizes of osteosarcoma are small. More data is needed to verify the accuracy of the signature. Second, because the microarray data for osteosarcoma are relatively scarce and the platforms used are relatively old, which are not sensitive enough for lncRNA, we did not find other external cohorts with survival information to validate the signature. Third, we need further experiments, such as immunohistochemical analysis, PCR, or western blot, to support these findings. Despite these limitations, we constructed a six lncRNAs signature with a good prognostic value in osteosarcoma. These lncRNAs may play essential roles in the progression of osteosarcoma. They deserve further study.

Conclusion

In conclusion, we identified an immune-related lncRNAs for osteosarcoma, which possess good performance in the prediction of survival, metastasis and immune microenvironment. Our study not only has great significance in predicting the prognosis but also has potential to guide future immunotherapy in osteosarcoma.

Data Availability Statement

Publicly available datasets were analyzed in this study. This data can be found here: https://portal.gdc.cancer.gov/.

Author Contributions

HC and YH performed the conception and design of this manuscript. HZ provided useful suggestions in methodology. HZ and HX performed data analysis and prepared the figures. YH drafted the manuscript. HC and HY revised the manuscript. All authors read and approved the final manuscript.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Abbreviations

LncRNA, Long non-coding RNA; TCGA, The Cancer Genome Atlas; ROC, Operating characteristic curve; PCA, Principal component analysis; GSEA, Gene set enrichment analysis; AUC, Area under curve; LASSO, Least absolute shrinkage and selection operator; K-M, Kaplan–Meier.

References

1. Luetke A, Meyers PA, Lewis I, Juergens H. Osteosarcoma Treatment - Where Do We Stand? A State of the Art Review. Cancer Treat Rev (2014) 40(4):523–32. doi: 10.1016/j.ctrv.2013.11.006

PubMed Abstract | CrossRef Full Text | Google Scholar

2. Ritter J, Bielack SS. Osteosarcoma. Ann Oncol (2010) 21 Suppl 7:vii320–5. doi: 10.1093/annonc/mdq276

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Ottaviani G, Jaffe N. The Epidemiology of Osteosarcoma. Cancer Treat Res (2009) 152:3–13. doi: 10.1007/978-1-4419-0284-9_1

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Rickel K, Fang F, Tao J. Molecular Genetics of Osteosarcoma. Bone (2017) 102:69–79. doi: 10.1016/j.bone.2016.10.017

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Evola FR, Costarella L, Pavone V, Caff G, Cannavo L, Sessa A, et al. Biomarkers of Osteosarcoma, Chondrosarcoma, and Ewing Sarcoma. Front Pharmacol (2017) 8:150. doi: 10.3389/fphar.2017.00150

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Chen R, Wang G, Zheng Y, Hua Y, Cai Z. Drug Resistance-Related microRNAs in Osteosarcoma: Translating Basic Evidence Into Therapeutic Strategies. J Cell Mol Med (2019) 23(4):2280–92. doi: 10.1111/jcmm.14064

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Chindamo G, Sapino S, Peira E, Chirio D, Gonzalez MC, Gallarate M. Bone Diseases: Current Approach and Future Perspectives in Drug Delivery Systems for Bone Targeted Therapeutics. Nanomaterial (Basel) (2020) 10(5):875. doi: 10.3390/nano10050875

CrossRef Full Text | Google Scholar

8. Izadpanah S, Shabani P, Aghebati-Maleki A, Baghbanzadeh A, Fotouhi A, Bisadi A, et al. Prospects for the Involvement of Cancer Stem Cells in the Pathogenesis of Osteosarcoma. J Cell Physiol (2020) 235(5):4167–82. doi: 10.1002/jcp.29344

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Yu WD, Wang H, He QF, Xu Y, Wang XC. Long Noncoding RNAs in Cancer-Immunity Cycle. J Cell Physiol (2018) 233(9):6518–23. doi: 10.1002/jcp.26568

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Huarte M. The Emerging Role of lncRNAs in Cancer. Nat Med (2015) 21(11):1253–61. doi: 10.1038/nm.3981

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Chen Y, Li C, Pan Y, Han S, Feng B, Gao Y, et al. The Emerging Role and Promise of Long Noncoding RNAs in Lung Cancer Treatment. Cell Physiol Biochem (2016) 38(6):2194–206. doi: 10.1159/000445575

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Bhan A, Soleimani M, Mandal SS. Long Noncoding RNA and Cancer: A New Paradigm. Cancer Res (2017) 77(15):3965–81. doi: 10.1158/0008-5472.Can-16-2634

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Fu XM, Guo W, Li N, Liu HZ, Liu J, Qiu SQ, et al. The Expression and Function of Long Noncoding RNA lncRNA-ATB in Papillary Thyroid Cancer. Eur Rev Med Pharmacol Sci (2017) 21(14):3239–46.

PubMed Abstract | Google Scholar

14. Pasic I, Shlien A, Durbin AD, Stavropoulos DJ, Baskin B, Ray PN, et al. Recurrent Focal Copy-Number Changes and Loss of Heterozygosity Implicate Two Noncoding RNAs and One Tumor Suppressor Gene at Chromosome 3q13.31 in Osteosarcoma. Cancer Res (2010) 70(1):160–71. doi: 10.1158/0008-5472.Can-09-1902

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Zhou H, Chen A, Shen J, Zhang X, Hou M, Li J, et al. Long Non-Coding RNA LOC285194 Functions as a Tumor Suppressor by Targeting P53 in Non-Small Cell Lung Cancer. Oncol Rep (2019) 41(1):15–26. doi: 10.3892/or.2018.6839

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Renganathan A, Felley-Bosco E. Long Noncoding RNAs in Cancer and Therapeutic Potential. Adv Exp Med Biol (2017) 1008:199–222. doi: 10.1007/978-981-10-5203-3_7

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, et al. Gene Set Enrichment Analysis: A Knowledge-Based Approach for Interpreting Genome-Wide Expression Profiles. Proc Natl Acad Sci USA (2005) 102(43):15545–50. doi: 10.1073/pnas.0506580102

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Liberzon A, Birger C, Thorvaldsdóttir H, Ghandi M, Mesirov JP, Tamayo P. The Molecular Signatures Database (MSigDB) Hallmark Gene Set Collection. Cell Syst (2015) 1(6):417–25. doi: 10.1016/j.cels.2015.12.004

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Yoshihara K, Shahmoradgoli M, Martínez E, Vegesna R, Kim H, Torres-Garcia W, et al. Inferring Tumour Purity and Stromal and Immune Cell Admixture From Expression Data. Nat Commun (2013) 4:2612. doi: 10.1038/ncomms3612

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Chen B, Khodadoust MS, Liu CL, Newman AM, Alizadeh AA. Profiling Tumor Infiltrating Immune Cells With CIBERSORT. Methods Mol Biol (2018) 1711:243–59. doi: 10.1007/978-1-4939-7493-1_12

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Rooney MS, Shukla SA, Wu CJ, Getz G, Hacohen N. Molecular and Genetic Properties of Tumors Associated With Local Immune Cytolytic Activity. Cell (2015) 160(1-2):48–61. doi: 10.1016/j.cell.2014.12.033

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Wang L, Felts SJ, Van Keulen VP, Scheid AD, Block MS, Markovic SN, et al. Integrative Genome-Wide Analysis of Long Noncoding RNAs in Diverse Immune Cell Types of Melanoma Patients. Cancer Res (2018) 78(15):4411–23. doi: 10.1158/0008-5472.Can-18-0529

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Li M, Chen H, Zhao Y, Gao S, Cheng C. H19 Functions as a ceRNA in Promoting Metastasis Through Decreasing miR-200s Activity in Osteosarcoma. DNA Cell Biol (2016) 35(5):235–40. doi: 10.1089/dna.2015.3171

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Sun J, Wang X, Fu C, Wang X, Zou J, Hua H, et al. Long Noncoding RNA FGFR3-AS1 Promotes Osteosarcoma Growth Through Regulating Its Natural Antisense Transcript Fgfr3. Mol Biol Rep (2016) 43(5):427–36. doi: 10.1007/s11033-016-3975-1

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Zhang CL, Zhu KP, Shen GQ, Zhu ZS. A Long non-Coding RNA Contributes to Doxorubicin Resistance of Osteosarcoma. Tumour Biol (2016) 37(2):2737–48. doi: 10.1007/s13277-015-4130-7

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Shen P, Cheng Y. Long Noncoding RNA lncARSR Confers Resistance to Adriamycin and Promotes Osteosarcoma Progression. Cell Death Dis (2020) 11(5):362. doi: 10.1038/s41419-020-2573-2

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Li X, Meng Y. Survival Analysis of Immune-Related lncRNA in Low-Grade Glioma. BMC Cancer (2019) 19(1):813. doi: 10.1186/s12885-019-6032-3

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Khadirnaikar S, Kumar P, Pandi SN, Malik R, Dhanasekaran SM, Shukla SK. Immune Associated LncRNAs Identify Novel Prognostic Subtypes of Renal Clear Cell Carcinoma. Mol Carcinog (2019) 58(4):544–53. doi: 10.1002/mc.22949

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Polev DE, Karnaukhova IK, Krukovskaya LL, Kozlov AP. ELFN1-AS1: A Novel Primate Gene With Possible microRNA Function Expressed Predominantly in Human Tumors. BioMed Res Int (2014) 2014:398097. doi: 10.1155/2014/398097

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Zhang C, Lian H, Xie L, Yin N, Cui Y. LncRNA ELFN1-AS1 Promotes Esophageal Cancer Progression by Up-Regulating GFPT1 via Sponging miR-183-3p. Biol Chem (2020) 401(9):1053–61. doi: 10.1515/hsz-2019-0430

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Lei R, Feng L, Hong D. ELFN1-AS1 Accelerates the Proliferation and Migration of Colorectal Cancer via Regulation of miR-4644/TRIM44 Axis. Cancer Biomark (2020) 27(4):433–43. doi: 10.3233/cbm-190559

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Berntsson J, Nodin B, Eberhard J, Micke P, Jirström K. Prognostic Impact of Tumour-Infiltrating B Cells and Plasma Cells in Colorectal Cancer. Int J Cancer (2016) 139(5):1129–39. doi: 10.1002/ijc.30138

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Yeong J, Lim JCT, Lee B, Li H, Chia N, Ong CCH, et al. High Densities of Tumor-Associated Plasma Cells Predict Improved Prognosis in Triple Negative Breast Cancer. Front Immunol (2018) 9:1209. doi: 10.3389/fimmu.2018.01209

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Gholiha AR, Hollander P, Hedstrom G, Sundstrom C, Molin D, Smedby KE, et al. High Tumour Plasma Cell Infiltration Reflects an Important Microenvironmental Component in Classic Hodgkin Lymphoma Linked to Presence of B-Symptoms. Br J Haematol (2019) 184(2):192–201. doi: 10.1111/bjh.15703

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Mills JD, Kavanagh T, Kim WS, Chen BJ, Waters PD, Halliday GM, et al. High Expression of Long Intervening Non-Coding RNA OLMALINC in the Human Cortical White Matter Is Associated With Regulation of Oligodendrocyte Maturation. Mol Brain (2015) 8:2. doi: 10.1186/s13041-014-0091-9

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Marone G, Varricchi G, Loffredo S, Granata F. Mast Cells and Basophils in Inflammatory and Tumor Angiogenesis and Lymphangiogenesis. Eur J Pharmacol (2016) 778:146–51. doi: 10.1016/j.ejphar.2015.03.088

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Takanami I, Takeuchi K, Naruke M. Mast Cell Density Is Associated With Angiogenesis and Poor Prognosis in Pulmonary Adenocarcinoma. Cancer (2000) 88(12):2686–92. doi: 10.1002/1097-0142(20000615)88:12<2686::AID-CNCR6>3.0.CO;2-6

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Ammendola M, Sacco R, Sammarco G, Donato G, Montemurro S, Ruggieri E, et al. Correlation Between Serum Tryptase, Mast Cells Positive to Tryptase and Microvascular Density in Colo-Rectal Cancer Patients: Possible Biological-Clinical Significance. PloS One (2014) 9(6):e99512. doi: 10.1371/journal.pone.0099512

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Rao Q, Chen Y, Yeh CR, Ding J, Li L, Chang C, et al. Recruited Mast Cells in the Tumor Microenvironment Enhance Bladder Cancer Metastasis via Modulation of Erβ/CCL2/CCR2 EMT/MMP9 Signals. Oncotarget (2016) 7(7):7842–55. doi: 10.18632/oncotarget.5467

PubMed Abstract | CrossRef Full Text | Google Scholar

40. Tu JF, Pan HY, Ying XH, Lou J, Ji JS, Zou H. Mast Cells Comprise the Major of Interleukin 17-Producing Cells and Predict a Poor Prognosis in Hepatocellular Carcinoma. Med (Baltimore) (2016) 95(13):e3220. doi: 10.1097/md.0000000000003220

CrossRef Full Text | Google Scholar

41. Chan JK, Magistris A, Loizzi V, Lin F, Rutgers J, Osann K, et al. Mast Cell Density, Angiogenesis, Blood Clotting, and Prognosis in Women With Advanced Ovarian Cancer. Gynecol Oncol (2005) 99(1):20–5. doi: 10.1016/j.ygyno.2005.05.042

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Hedström G, Berglund M, Molin D, Fischer M, Nilsson G, Thunberg U, et al. Mast Cell Infiltration Is a Favourable Prognostic Factor in Diffuse Large B-Cell Lymphoma. Br J Haematol (2007) 138(1):68–71. doi: 10.1111/j.1365-2141.2007.06612.x

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Rajput AB, Turbin DA, Cheang MC, Voduc DK, Leung S, Gelmon KA, et al. Stromal Mast Cells in Invasive Breast Cancer are a Marker of Favourable Prognosis: A Study of 4,444 Cases. Breast Cancer Res Treat (2008) 107(2):249–57. doi: 10.1007/s10549-007-9546-3

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Crespo J, Wu K, Li W, Kryczek I, Maj T, Vatan L, et al. Human Naive T Cells Express Functional CXCL8 and Promote Tumorigenesis. J Immunol (2018) 201(2):814–20. doi: 10.4049/jimmunol.1700755

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Su S, Liao J, Liu J, Huang D, He C, Chen F, et al. Blocking the Recruitment of Naive CD4(+) T Cells Reverses Immunosuppression in Breast Cancer. Cell Res (2017) 27(4):461–82. doi: 10.1038/cr.2017.34

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Centuori SM, Gomes CJ, Kim SS, Putnam CW, Larsen BT, Garland LL, et al. Double-Negative (CD27(-)IgD(-)) B Cells Are Expanded in NSCLC and Inversely Correlate With Affinity-Matured B Cell Populations. J Transl Med (2018) 16(1):30. doi: 10.1186/s12967-018-1404-z

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Ni Z, Xing D, Zhang T, Ding N, Xiang D, Zhao Z, et al. Tumor-Infiltrating B Cell Is Associated With the Control of Progression of Gastric Cancer. Immunol Res (2021) 69(1):43–52. doi: 10.1007/s12026-020-09167-z

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Adema GJ. Dendritic Cells From Bench to Bedside and Back. Immunol Lett (2009) 122(2):128–30. doi: 10.1016/j.imlet.2008.11.017

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Broderick L, Yokota SJ, Reineke J, Mathiowitz E, Stewart CC, Barcos M, et al. Human CD4+ Effector Memory T Cells Persisting in the Microenvironment of Lung Cancer Xenografts are Activated by Local Delivery of IL-12 to Proliferate, Produce IFN-Gamma, and Eradicate Tumor Cells. J Immunol (2005) 174(2):898–906. doi: 10.4049/jimmunol.174.2.898

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Boussommier-Calleja A, Atiyas Y, Haase K, Headley M, Lewis C, Kamm RD. The Effects of Monocytes on Tumor Cell Extravasation in a 3D Vascularized Microfluidic Model. Biomaterials (2019) 198:180–93. doi: 10.1016/j.biomaterials.2018.03.005

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Yamaguchi T, Fushida S, Yamamoto Y, Tsukada T, Kinoshita J, Oyama K, et al. Tumor-Associated Macrophages of the M2 Phenotype Contribute to Progression in Gastric Cancer With Peritoneal Dissemination. Gastric Cancer (2016) 19(4):1052–65. doi: 10.1007/s10120-015-0579-8

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Yao RR, Li JH, Zhang R, Chen RX, Wang YH. M2-Polarized Tumor-Associated Macrophages Facilitated Migration and Epithelial-Mesenchymal Transition of HCC Cells via the TLR4/STAT3 Signaling Pathway. World J Surg Oncol (2018) 16(1):9. doi: 10.1186/s12957-018-1312-y

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Gomez-Brouchet A, Illac C, Gilhodes J, Bouvier C, Aubert S, Guinebretiere JM, et al. CD163-Positive Tumor-Associated Macrophages and CD8-Positive Cytotoxic Lymphocytes are Powerful Diagnostic Markers for the Therapeutic Stratification of Osteosarcoma Patients: An Immunohistochemical Analysis of the Biopsies Fromthe French OS2006 Phase 3 Trial. Oncoimmunology (2017) 6(9):e1331193. doi: 10.1080/2162402x.2017.1331193

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Zhang C, Zheng JH, Lin ZH, Lv HY, Ye ZM, Chen YP, et al. Profiles of Immune Cell Infiltration and Immune-Related Genes in the Tumor Microenvironment of Osteosarcoma. Aging (Albany NY) (2020) 12(4):3486–501. doi: 10.18632/aging.102824

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Heymann MF, Lézot F, Heymann D. The Contribution of Immune Infiltrates and the Local Microenvironment in the Pathogenesis of Osteosarcoma. Cell Immunol (2019) 343:103711. doi: 10.1016/j.cellimm.2017.10.011

PubMed Abstract | CrossRef Full Text | Google Scholar