An Immune-Related Long Non-Coding RNA Signature to Predict the Prognosis of Ewing’s Sarcoma Based on a Machine Learning Iterative Lasso Regression

The aim of this study was to construct a new immune-associated long non-coding RNA (lncRNA) signature to predict the prognosis of Ewing sarcoma (ES) and explore its molecular mechanisms. We downloaded transcriptome and clinical prognosis data from the Gene Expression Omnibus (GSE17679, which included 88 ES samples and 18 matched normal skeletal muscle samples), and used it as a training set to identify immune-related lncRNAs with different expression levels in ES. Univariable Cox regression was used to screen immune-related lncRNAs related to ES prognosis, and an immune-related lncRNA signature was constructed based on machine learning iterative lasso regression. An external verification set was used to confirm the predictive ability of the signature. Clinical feature subgroup analysis was used to explore whether the signature was an independent prognostic factor. In addition, CIBERSORT was used to explore immune cell infiltration in the high- and low-risk groups, and to analyze the correlations between the lncRNA signature and immune cell levels. Gene set enrichment and variation analyses were used to explore the possible regulatory mechanisms of the immune-related lncRNAs in ES. We also analyzed the expression of 17 common immunotherapy targets in the high- and low-risk groups to identify any that may be regulated by immune-related lncRNAs. We screened 35 immune-related lncRNAs by univariate Cox regression. Based on this, an immune-related 11-lncRNA signature was generated by machine learning iterative lasso regression. Analysis of the external validation set confirmed its high predictive ability. DPP10 antisense RNA 3 was negatively correlated with resting dendritic cell, neutrophil, and γδ T cell infiltration, and long intergenic non-protein coding RNA 1398 was positively correlated with resting dendritic cells and M2 macrophages. These lncRNAs may affect ES prognosis by regulating GSE17721_CTRL_VS_PAM3CSK4_12H_BMDC_UP, GSE2770_IL4_ACT_VS_ACT_CD4_TCELL_48H_UP, GSE29615_CTRL_VS_DAY3_ LAIV_IFLU_VACCINE_PBMC_UP, complement signaling, interleukin 2-signal transducer and activator of transcription 5 signaling, and protein secretion. The immune-related 11-lncRNA signature may also have regulatory effects on the immunotherapy targets CD40 molecule, CD70 molecule, and CD276 molecule. In conclusion, we constructed a new immune-related 11-lncRNA signature that can stratify the prognoses of patients with ES.

Introduction

Ewing sarcoma (ES) is one of the most common malignant tumors in children, young adults, and adults (Grünewald et al., 2018). In the past two decades, there has been great progress in ES treatment, through surgery, radiotherapy, and intensive chemotherapy (Pappo and Dirksen, 2018), and patient prognosis has significantly improved. The current 5-year survival rate of patients with local ES is >70%; however, the 5-year survival rate of patients with metastatic or recurrent ES tumors remains at only 20–30% (Burdach and Jürgens, 2002; Iwamoto, 2007). Unfortunately, breakthroughs in the treatment of recurrent and metastatic ES have been difficult to achieve. The precise classification of patients with different prognoses is crucial for precise ES treatment. ES prognosis is closely related to immune factors. For example, CD8+ T cells can kill ES cells by specifically recognizing the ET-derived antigens enhancer of zeste 2 polycomb repressive complex 2 subunit 666 and chondromodulin 319 (Thiel et al., 2011; Blaeschke et al., 2016). Natural killer (NK) cells do not recognize specific tumor antigens to cause an immune response, but exert a direct killing effect on ES cells. Studies have shown that allogeneic transplantation of NK cells has a more pronounced killing effect on tumors than autologous NK cells (Ljunggren and Malmberg, 2007; Verhoeven et al., 2008). In addition, macrophages, mast cells, antigen presenting cells, and dendritic cells are also involved in the molecular mechanisms of ES (Dagher et al., 2002; Lau et al., 2007; Inagaki et al., 2016; Fedorova et al., 2019); however, their specific roles remain unclear. Studies have shown that interleukin (IL)-6, IL-10, and killer cell lectin like receptor K1 regulate the ES tumor microenvironment and are closely related to its prognosis (Hempel et al., 1997; Berghuis et al., 2012; Lissat et al., 2015). Therefore, immune-related prediction signatures may provide accurate guidance for ES treatment.

The discovery of the first long non-coding RNA (lncRNA) (Rinn et al., 2007) sparked an entire field of research regarding their effects and molecular mechanisms in disease. LncRNAs can regulate gene expression through signals, decoys, guides, and scaffolds (Ingolia et al., 2011). Increasing studies have shown that lncRNAs can not only regulate immune responses, but also play important roles in the molecular mechanisms controlling tumors, and are closely related to their prognosis (Peng et al., 2017; Ma et al., 2018; Mowel et al., 2018). According to Marques Howarth et al. (2014), the lncRNA Ewing sarcoma associated transcript 1 (EWSAT1) is a downstream target of EWS RNA binding protein 1 (EWSR1), and the proliferation of ES cells can be inhibited by inhibiting EWSAT1 expression. Immune-associated lncRNAs can be used as prognostic biomarkers for glioblastoma multiforme, breast cancer, and bladder cancer (Zhou et al., 2018; Cao et al., 2020; Zhang Y. et al., 2020). However, due to the lack of research on lncRNAs involved in the molecular mechanisms of ES, an immune-related lncRNA prognosis signature has not been reported. In this study, we have identified lncRNAs strongly related to ES prognosis and used machine learning iterative lasso regression to generate and validate an immune-related 11-lncRNA signature that can predict ES prognosis. We also explored its correlations with immune cell infiltration, to provide accurate and reliable guidance for clinical individualized treatment.

Materials and Methods

ES Source Data and Identification of Differentially Expressed Immune-Related lncRNAs

Transcriptome data and corresponding clinical data from the GSE17679 dataset were downloaded from the Gene Expression Omnibus¹. The dataset includes 88 ES samples and 18 matched healthy skeletal muscle samples, and was used as the training set. The immune scores of the 88 ES samples were calculated using the ESTIMATE algorithm (Yoshihara et al., 2013, https://bioinformatics.mdanderson.org/public-software/estimate/), and they were divided into high and low immune infiltration groups accordingly. The Stromal scores, ESTIMATE scores, and tumor purity levels of the two groups were evaluated. The limma package was used to compare the transcriptome data between the groups to identify immune-related lncRNAs, and differential expression analysis was performed between the ES samples and healthy skeletal muscle samples to identify differentially expressed lncRNAs. The intersection of immune-related lncRNAs and differentially expressed lncRNAs is regarded as immune-related and differentially expressed lncRNAs. Transcriptome and clinical data from 58 cases of ES were downloaded from the International Cancer Genome Consortium database for use as an external validation set.

Construction of an Optimal Immune-Related lncRNA Predictive Signature

A conditional probability survival graph can describe the survival of patients at different time stages in detail (Latenstein et al., 2020). In this study, 88 patients with ES were used to construct a predictive model, and the conditional survival rate was determined using a conditional probability survival curve. Univariate Cox regression was used to identify lncRNAs associated with ES prognosis. The screening criterion was p < 0.05. Lasso regression (Frost and Amos, 2017) is mainly used for the supervised learning of high-dimensional data. Each iteration of the regression produces a gene combination related to prognosis. We conducted 500 lasso regressions on candidate lncRNAs, and considered the lncRNA combination with the largest area under curve (AUC) of the receiver operating characteristic (ROC) as the optimal lncRNA signature (Sveen et al., 2012). We also evaluated the optimal prognosis ability of the lncRNA signature, in terms of overall survival time and lncRNA expression according to risk score.

Verification of the Optimal Immune-Related lncRNA Predictive Signature

To evaluate the reliability of the immune-related lncRNA signature, we evaluated its prognostic value using the external validation set, and calculated the AUC after 3, 5, and 8 years. We also compared the prognostic value of the immune-related lncRNA signature with established ES prognostic biomarkers (BIK, EGFR, CD44, and LGR5) using the external validation set.

Age, gender, and metastasis are common factors affecting ES. To detect whether the optimal immune-associated lncRNA signature was an independent prognostic factor, Kaplan-Meier survival analysis was used to compare the survival of the high and low-risk groups based on these clinical characteristics. Time-dependent AUC analysis can evaluate the consistency indices of different models based on that of the survival model, and was used to verify the accuracy of the lncRNA signature in predicting ES prognosis (compared with individual clinical characteristics and the lncRNA signature and clinical characteristics combined).

Analysis of Correlations Between lncRNAs in the Optimal Signature and Immune Cells

CIBERSORT (Newman et al., 2019) is an online tool for immune cell subtype deconvolution based on the principle of linear support vector regression. It can use transcriptome data to evaluate the infiltration of 22 types of immune cells. We used CIBERSORT to analyze immune cell infiltration in the high- and low-risk groups. PCA clustering was performed on the filtered data to detect differences between the groups, and the ggplot2 package was used for visualization. The corrplot, ggplot2, and igraph packages were used in R to visualize the correlations, infiltration differences, and interactions, respectively, between 22 kinds of immune cells. Kaplan-Meier analysis was used to explore relationships between immune cell infiltration and ES prognosis. To explore correlations between the lncRNA signature and prognosis-related immune cell infiltration, Pearson correlation analysis was performed and the ggplot2 package was used for visualization.

Exploration of Immune Checkpoints and Related Pathways

CD27 molecule, CD40, CD70, TNF receptor superfamily member 14, CD276, V-set domain containing T cell activation inhibitor 1, indoleamine 2,3-dioxygenase 1, programmed cell death 1, CD274 molecule, programmed cell death 1 ligand 2, hepatitis A virus cellular receptor 2, T cell immunoreceptor with Ig and ITIM domains, cytotoxic T-lymphocyte associated protein 4, CD86 molecule, inducible T cell co-stimulator, lymphocyte activating 3, and CD58 molecule are the most common immune checkpoint markers used in tumor research. We explored the expression of these common immune checkpoints in the high- and low-risk groups.

To explore the enrichment of important pathways in the high-risk group, we conducted gene set enrichment analysis (GSEA; Subramanian et al., 2005) and gene set variation analysis (GSVA; Hänzelmann et al., 2013). GSEA was performed in GSEA 4.0.3 using “h.all.v7.1.symbols.gmt” and “c7.all.v7.1.symbols.gmt” as reference gene sets. Nominal p-values < 0.05 and false discovery rates < 0.05 were considered significant. GSVA was performed on the “h.all.v7.1 symbols.gmt” gene set using the cluster profiler and gsva packages.

Results

Identification of Differentially Expressed Immune-Related lncRNAs

The tumor purity of the samples were evaluated using the time-of-life method. According to their immune scores, the 88 ES samples were divided into high and low immune cell infiltration groups (n = 44 each; Figure 1A). The ESTIMATE scores (p < 0.001) and the Stromal scores (p < 0.001) was higher in high immune cell infiltration group, while tumor purity (p < 0.001) was lower (Figures 1A–D). Principal component analysis (PCA; Supplementary Figures 1A–D) revealed dramatic differences between the two groups. We obtained 262 immune-related lncRNAs by differential expression analysis of the high and low immune infiltration groups, and 884 differentially expressed lncRNAs in the 88 ES samples compared to the 18 matched healthy skeletal muscle samples. The intersection of these groups of lncRNAs yielded 171 immune-related and differentially expressed lncRNAs (Supplementary Figure 1E).

FIGURE 1

Figure 1. Classification of 88 ES samples based on immune score. (A) Box plot of the immune scores of the high- and low-risk groups. (B) Box plot of stromal scores. (C) Box plot of ESTIMATE scores. (D) Box diagram of tumor purity.

Construction of an Immune-Related lncRNA Prognostic Signature

The annual conditional survival probability increased with the overall survival time (Figure 2). From a 49% chance of survival immediately post-resection, the probability of 5-year survival 1, 2, 3, and 4 years after resection increased by 58, 72, 84, and 100%, respectively. The probability of surviving the next year decreased from 84% to 81% after 1 year, and then increased to 84 and 89% at 3 and 5 years, respectively. The univariate Cox regression model identified 35 correlations between immune-related differentially expressed lncRNAs and patient prognosis (Figure 3). Prognosis-related lncRNAs were cross-validated via 500 lasso regressions to reveal an optimal immune-related lncRNA prognostic model consisting of 11 lncRNAs (Figure 4A). Receiver operating characteristic (ROC) analysis was used to further evaluate the predictive performance of the immune-related lncRNA prognostic signature. The results show that it has good performance in predicting ES prognosis (AUC = 0.819; Figure 4B). When the patients were divided into high- and low-risk groups according to their risk scores, Kaplan-Meier survival analysis showed that ES prognosis was significantly worse in the high-risk group (log-rank p < 0.001, Figure 4C). Both the risk scores and the number of deaths in the high-risk group were significantly higher than those of the low-risk group (Figure 4D). Of the 11 immune-related lncRNAs, ARHGAP26 antisense RNA 1 (ARHGAP26-AS1), FUT8 antisense RNA 1 (FUT8-AS1), FOXC1 upstream transcript (FOXCUT), and chromosome 5 putative open reading frame 64 (C5orf64) were highly expressed in the high-risk group, while NAV2 antisense RNA 2 (NAV2-AS2), long intergenic non-protein coding RNA (LINC)00408, SEC24B antisense RNA 1 (SEC24B-AS1), LINC01343, LINC01398, LINC01197, and DPP10 antisense RNA 3 (DPP10-AS3) were lowly expressed in the high-risk group (Figure 4D).

FIGURE 2

Figure 2. Conditional survival estimates over time. Each column represents a survival period, and each row represents the percentage to reach a certain survival time from that point (in years).

FIGURE 3

Figure 3. Univariate Cox regression analysis of immune-related differentially expressed lncRNAs.

FIGURE 4

Figure 4. Construction and evaluation of the optimal immune-related lncRNA signature. (A) Line chart of the AUC of different immune-related lncRNA signature models. (B) ROC curve of the immune-related 11-lncRNA signature. (C) Kaplan-Meier analysis of the signature. Survival was compared using the log-rank test. (D) Evaluation of the 11-lncRNA signature based on risk factors in the high- and low-risk groups, the RFS, and gene expression in the signature.

Verification of the Optimal Immune-Related lncRNA Signature

To verify the reliability of the optimal immune lncRNA signature, we evaluated its predictive value in 58 ES samples in the external validation set through ROC analysis. The immune-related lncRNA signature had obvious prognostic value after 3 (AUC = 0.71), 5 (AUC = 0.68), and 8 years (AUC = 0.75; Figure 5A). Compared with prognostic biomarkers such as BCL2 interacting killer (BIK), epidermal growth factor receptor (EGFR), CD44 molecule (Indian blood group) (CD44), and leucine rich repeat containing G protein-coupled receptor 5 (LGR5), the lncRNA signature had better prognostic value (Figure 5B).

FIGURE 5

Figure 5. Verification of the optimal immune-related 11-lncRNA signature. (A) 3-, 5-, and 8-year ROC curves of the 11-lncRNA signature in the external verification set. (B) Comparison of the 11-lncRNA signature and common prognostic biomarkers of ES.

Evaluation of the Immune-Related lncRNA Signature as an Independent ES Prognostic Factor and Its Prediction Accuracy

To assess whether the immune-related lncRNA signature acts is a prognostic factor independent of clinical characteristics (age, sex, and metastasis), Kaplan-Meier survival analysis was performed in different subgroups of the high- and low-risk groups. Prognosis was poor in all subgroups of the high-risk group (p < 0.05, Figure 6), suggesting that the signature is independent of age, sex, and metastasis. The lncRNA signature was more accurate in predicting ES prognosis than age, sex, age + sex, and age + sex + lncRNA signature models (Figure 7).

FIGURE 6

Figure 6. Survival analysis of the immune-related 11-lncRNA signature with different clinical characteristics.

FIGURE 7

Figure 7. Time-dependent ROC curve. Concordance index (cindex) shows measure of concordance of predictor with survival of patients. The minimum threshold for cindex to be considered effective is 0.6.

Correlations Between the Immune-Related lncRNA Prognostic Signature and Immune Cell Subtype Infiltration

Principal component analysis of the high- and low-risk groups revealed differences in immune cell infiltration (Figure 8A). Correlation analysis showed that plasma cells were positively correlated with M1 macrophages and resting mast cells, but negatively correlated with M2 macrophages. activated T cells CD4 memory were positively correlated with γδ T cells and negatively correlated with M0 macrophages (Figure 8B). M2 macrophages, resting NK cells, and activated NK cells had the strongest interactions with other immune cells, while monocytes, naïve B cells, and resting dendritic cells had the weakest interactions (Figure 8C). Memory B cells and activated NK cells showed higher infiltration in the high-risk group compared to the low-risk group (Figure 8D). Infiltration of regulatory T cells (Tregs; p = 0.001) and activated CD4 memory T cells (p = 0.001) indicated good ES prognosis, while infiltration by activated dendritic cells (p = 0.009), M2 macrophages (p = 0.001), monocytes (p < 0.001), resting mast cells (p < 0.001), and γδ T cells (p < 0.001) indicated poor prognosis (Figure 9A). DPP10-AS3 was positively correlated with resting dendritic cell, neutrophil, and γδ T cell infiltration, while LINC01398 was negatively correlated with resting dendritic cell and M2 macrophage infiltration (Figure 9B).

FIGURE 8

Figure 8. Immune cell infiltration in ES. (A) PCA analysis of immune cell infiltration in ES samples and healthy skeletal muscle samples. (B) Correlation heat maps of 22 immune cell types. The size of each colored circle represents the related p-value; the color represents the strength of the correlation (blue and red indicate positive and negative correlations, respectively, with darker colors indicating stronger correlations). (C) Interaction network between 22 immune cell types. The circle size represents the strength of the interaction. (D) Violin chart showing differences in the infiltration of 22 immune cell types in the ES and healthy skeletal muscle groups.

FIGURE 9

Figure 9. Correlation analysis between the immune-related 11-lncRNA signature and prognosis-related immune cells. (A) Kaplan-Meier survival analysis of immune cell infiltration and ES prognosis. (B) Heat map of correlations between the 11-lncRNA signature and prognostic immune cells. The ordinate is the gene name, the abscissa is the immune cell type, and the color represents the correlation coefficient; ^∗p < 0.05, ^∗∗p < 0.01.

Signature-Related Pathways and Immune Checkpoint Markers

The results of gene set enrichment analysis (GSEA) are shown in Figures 10A,B. Pathways enriched in the high-risk group included GSE17721_CTRL_VS_PAM3CSK4_12H_BMDC_UP, GSE2770_IL4_ACT_VS_ACT_CD4_TCELL_48H_UP, GSE29615_CTRL_VS_DAY3_LAIV_STAT_VACCEMENT, REMARK_MARKINHALL_COM_PLTION, REMARK_COM_PL_COM, REMARK_COMP_UP, and REMARK_COMP_PROG. Gene set variation analysis (GSVA) revealed activation of IL2-signal transducer and activator of transcription 5, protein secretion, complement, and phosphoinositide 3-kinase (PI3K)-Akt-mammalian target of rapamycin signaling in the high-risk group (Figure 10C).

FIGURE 10

Figure 10. GSEA and GSVA analyses. (A) GSEA analysis based on h.all.v7.1.symbols.gmt. (B) GSEA analysis based on c7.all.v7.1.symbols.gmt. (C) GSVA analysis.

Among common immune checkpoint markers, the levels of CD40 molecule (CD40; p = 0.01) and CD70 molecule (CD70; p = 0.019) were higher in the high-risk group, while CD276 molecule (CD276; p = 0.019) was higher in the low-risk group (Figure 11).

FIGURE 11

Figure 11. Immune checkpoint expression in the high and low risk groups. (A–C) Expression of (A) CD40, (B) CD70, and (C) CD276.

Discussion

Increasing studies have shown that lncRNAs play important roles in the occurrence and development of various tumors. LncRNAs are not only involved in tumor regulatory mechanisms, but their levels are also closely related to patient prognosis. For example, the lncRNA CBR3 antisense RNA 1 can not only promote the occurrence of osteosarcoma by regulating the proliferation, migration, invasion, and apoptosis of osteosarcoma cells, but is also an independent prognostic factor of the disease (Zhang et al., 2018). This study aimed to identify an optimal immune-related lncRNA signature to predict the prognosis of ES. After screening for prognosis-related lncRNAs, the 11-lncRNA signature was constructed using a machine learning-iterative lasso regression model. Compared with the traditional stepwise regression method for constructing prognostic signatures, this method is based on the penalized lasso regression method, and combines lncRNAs with strong prognostic correlations to obtain optimal lncRNA signatures (Friedman et al., 2010; Goeman, 2010). This method not only considers the prognostic information of each lncRNA, but also removes redundant prognostic information, maximizing the prognostic value of the lncRNA signature. We also used bioinformatic methods to explore relationships between the lncRNA signature and prognosis-related immune cells, and explored the potential regulatory mechanisms involved, providing new research avenues in the study of immune-related lncRNAs in ES.

We identified 11 differentially expressed immune-related lncRNAs: ARHGAP26-AS1, FUT8-AS1, FOXCUT, C5orf64, NAV2-AS2, LINC00408, SEC24B-AS1, LINC01343, LINC01398, LINC01197, and DPP10-AS3. NAV2-AS2 and SEC24B-AS1 are prognostic biomarkers for lung adenocarcinoma (He and Zuo, 2019) and non-small cell lung cancer (Yang et al., 2020), respectively. Zhang X. et al. (2020) showed that FOXCUT promotes the metastasis and proliferation of colorectal cancer by activating the forkhead box C1/PI3K/Akt pathway. In addition, FOXCUT plays important roles in the molecular mechanisms of breast cancer, nasopharyngeal carcinoma, gastric adenocarcinoma, and esophageal squamous cell carcinoma, and can be used as a prognostic biomarker of esophageal squamous cell carcinoma (Pan et al., 2014; Liu et al., 2015; Xu et al., 2017; Zhao and Shen, 2019). The relationships between the 11 lncRNAs and ES are currently unclear, and biological studies will be required to explore their roles in its molecular mechanisms. Our results demonstrate that the immune-related 11-lncRNA signature has a higher prognostic value than other known prognostic biomarkers and is not affected by clinical characteristics. Therefore, the signature has strong prognostic evaluation value. However, large-scale experimental verification will be required before it can be accurately and reliably used in the clinic.

To evaluate immune cell infiltration in ES, we applied the deconvolution method to analyze ES expression data and found that memory B cells and activated NK cells had higher infiltration in the high-risk group than in the low-risk group, and that Tregs, activated CD4 memory T cells, activated dendritic cells, M2 macrophages, monocytes, resting mast cells, and γδ T cells were significantly related to ES prognosis. The killing effects of NK cells on ES cells have been experimentally verified (Verhoeven et al., 2008). Osteosarcoma displays mast cell infiltration (Inagaki et al., 2016); however, their infiltration of ES is currently unclear. The infiltration of tumor-associated macrophages and mast cells indicates poor ES prognosis (Pollard, 2004; Fujiwara et al., 2011), consistent with our results. In immunodeficient mice, γδ T cells can mediate the cytotoxicity of antibody-dependent ES cells with high expression of GD2 (Fisher et al., 2016). Guo et al. (2008) reported that dendritic cells have strong inhibitory effects on the proliferation of subcutaneous ES cells in mice. The roles of memory B cells, activated CD4 memory T cells, and Tregs in the occurrence and development of ES have not been reported. We also found that DPP10-AS3 was negatively correlated with resting dendritic cell, neutrophil, and γδ T cell infiltration, and LINC01398 was positively correlated with resting dendritic cell and M2 macrophage infiltration. However, the relationships between these lncRNAs and these prognostic-related immune cell types remain unclear, and will require further biological experimentation.

In GSEA and GSVA, GSE17721_CTRL_VS_PAM3CSK4_12H _BMDC_UP, GSE2770_IL4_ACT_VS_ACT_CD4_TCELL_48H_ UP, GSE29615_CTRL_VS_DAY3_LAIV_IFLU_VACCINE_ PBMC_UPLING, HSTATALLMARK_COMPLEMENT2_HALL MARK_MARK_UP_HALL_MARK_COMPLTION, and HSTAT ALL_COMPLTION_HALLMARK_COMPLTION_HALLMARK _COMPLTION were significantly enriched in the high risk group compared to the low-risk group. The complement system plays an important role in ES pathogenesis. For example, complement C5 is activated in ES and is positively correlated with better prognosis (Savola et al., 2011). In addition, a decrease in extracellular matrix protein secretion is related to loss of ES cell invasion ability (Javelaud et al., 2002). Therefore, the lncRNAs included in the signature may affect ES prognosis in part by regulating the complement system and the secretion of extracellular matrix proteins. The expression of the EWSR1-WT1, a fusion protein containing sections of EWSR1 and WT1 transcription factor, in proliferative small round cell tumors can induce the expression of IL-2 and IL-15, which are related to the proliferation of these tumor cells (Wong et al., 2002). However, the role of the IL-2-STAT5 signaling pathway in the pathogenesis of ES has not yet been reported. The specific role of protein secretion in ES pathogenesis also remains unclear, and further research is needed. Further study of these pathways will elucidate the immune regulatory mechanisms of ES.

Among common immune checkpoint markers, the levels of CD40, CD70, and CD276 differed between the high- and low-risk groups. Lollini et al. (1998) analyzed the expression of CD40 in 12 human osteosarcoma cell lines, six ES lines, and five rhabdomyosarcoma lines by flow cytometry. CD40 was highly expressed in osteosarcoma and ES, and was closely related to ES prognosis. CD70 is a therapeutic target for osteosarcoma (Pahl et al., 2015); however, its role in ES pathogenesis is currently unclear. CD276 is an immunotherapy target for peritoneal cancer, glioma, and central nervous tumors (Picarda et al., 2016), and its role in ES is also unclear. Our results show that the 11-lncRNA signature is closely related to these therapeutic targets. However, the specific regulatory relationships between the signature and CD40, CD70, and CD276 will require biological experimentation.

In summary, this study reports the first immune-related lncRNA signature related to ES prognosis. The gene signature is closely related to the infiltration of a variety of immune cell types, and reveals pathways and immune checkpoints that may be regulated by the 11 lncRNAs comprising it.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

Author Contributions

E-HR, Y-JD, and W-HY: conception and design, data analysis, and manuscript writing and revision. E-HR, Z-LL, and G-ZZ: data sorting and data checking. Q-QX and C-YL: data acquisition, research design, and picture drawing. 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.

Acknowledgments

We would like to thank editage (www.editage.cn) for English language editing.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fcell.2021.651593/full#supplementary-material

Abbreviations

EWSAT1, Ewing sarcoma associated transcript 1; EWSR1, EWS RNA binding protein 1; PCA, principal component analysis; ROC, receiver operating characteristic; ARHGAP26-AS1, ARHGAP26 antisense RNA 1; FUT8-AS1, FUT8 antisense RNA 1; FOXCUT, FOXC1 upstream transcript; C5orf64, chromosome 5 putative open reading frame 64; NAV2-AS2, NAV2 antisense RNA 2; LINC, long intergenic non-protein coding RNA; SEC24B-AS1, SEC24B antisense RNA 1; DPP10-AS3, DPP10 antisense RNA 3.

Footnotes

1. ^ https://www.ncbi.nlm.nih.gov/geo/

References

Berghuis, D., Schilham, M. W., Vos, H. I., Santos, S. J., Kloess, S., Buddingh, E. P., et al. (2012). Histone deacetylase inhibitors enhance expression of NKG2D ligands in Ewing sarcoma and sensitize for natural killer cell-mediated cytolysis. Clin. Sarcoma Res. 2:8. doi: 10.1186/2045-3329-2-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Blaeschke, F., Thiel, U., Kirschner, A., Thiede, M., Rubio, R. A., Schirmer, D., et al. (2016). Human HLA-A^∗02:01/CHM1+ allo-restricted T cell receptor transgenic CD8+ T cells specifically inhibit Ewing sarcoma growth in vitro and in vivo. Oncotarget 7, 43267–43280. doi: 10.18632/oncotarget.9218

PubMed Abstract | CrossRef Full Text | Google Scholar

Burdach, S., and Jürgens, H. (2002). High-dose chemoradiotherapy (HDC) in the Ewing family of tumors (EFT). Crit. Rev. Oncol. Hematol. 41, 169–189. doi: 10.1016/s1040-8428(01)00154-8

CrossRef Full Text | Google Scholar

Cao, R., Yuan, L., Ma, B., Wang, G., and Tian, Y. (2020). Immune-related long non-coding RNA signature identified prognosis and immunotherapeutic efficiency in bladder cancer (BLCA). Cancer Cell Int. 20:276. doi: 10.1186/s12935-020-01362-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Dagher, R., Long, L. M., Read, E. J., Leitman, S. F., Carter, C. S., Tsokos, M., et al. (2002). Pilot trial of tumor-specific peptide vaccination and continuous infusion interleukin-2 in patients with recurrent Ewing sarcoma and alveolar rhabdomyosarcoma: an inter-institute NIH study. Med. Pediatr. Oncol. 38, 158–164. doi: 10.1002/mpo.1303

PubMed Abstract | CrossRef Full Text | Google Scholar

Fedorova, L., Mudry, P., Pilatova, K., Selingerova, I., Merhautova, J., Rehak, Z., et al. (2019). Assessment of immune response following dendritic cell-based immunotherapy in pediatric patients with relapsing sarcoma. Front. Oncol. 9:1169. doi: 10.3389/fonc.2019.01169

PubMed Abstract | CrossRef Full Text | Google Scholar

Fisher, J. P., Flutter, B., Wesemann, F., Frosch, J., Rossig, C., Gustafsson, K., et al. (2016). Effective combination treatment of GD2-expressing neuroblastoma and Ewing’s sarcoma using anti-GD2 ch14.18/CHO antibody with Vγ9Vδ2+ γδT cells. Oncoimmunology 5:e1025194. doi: 10.1080/2162402x.2015.1025194

PubMed Abstract | CrossRef Full Text | Google Scholar

Friedman, J., Hastie, T., and Tibshirani, R. (2010). Regularization paths for generalized linear models via coordinate descent. J. Stat. Softw. 33, 1–22.

Google Scholar

Frost, H. R., and Amos, C. I. (2017). Gene set selection via LASSO penalized regression (SLPR). Nucleic Acids Res. 45:e114. doi: 10.1093/nar/gkx291

PubMed Abstract | CrossRef Full Text | Google Scholar

Fujiwara, T., Fukushi, J., Yamamoto, S., Matsumoto, Y., Setsu, N., Oda, Y., et al. (2011). Macrophage infiltration predicts a poor prognosis for human ewing sarcoma. Am. J. Pathol. 179, 1157–1170. doi: 10.1016/j.ajpath.2011.05.034

PubMed Abstract | CrossRef Full Text | Google Scholar

Goeman, J. J. (2010). L1 penalized estimation in the Cox proportional hazards model. Biom. J. 52, 70–84. doi: 10.1002/bimj.200900028

PubMed Abstract | CrossRef Full Text | Google Scholar

Grünewald, T. G. P., Cidre-Aranaz, F., Surdez, D., Tomazou, E. M., de Álava, E., Kovar, H., et al. (2018). Ewing sarcoma. Nat. Rev. Dis. Primers 4:5. doi: 10.1038/s41572-018-0003-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Guo, W., Guo, Y., Tang, S., Qu, H., and Zhao, H. (2008). Dendritic cell-Ewing’s sarcoma cell hybrids enhance antitumor immunity. Clin. Orthop. Relat. Res. 466, 2176–2183. doi: 10.1007/s11999-008-0348-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Hänzelmann, S., Castelo, R., and Guinney, J. (2013). GSVA: gene set variation analysis for microarray and RNA-seq data. BMC Bioinformatics 14:7. doi: 10.1186/1471-2105-14-7

PubMed Abstract | CrossRef Full Text | Google Scholar

He, R., and Zuo, S. (2019). A robust 8-gene prognostic signature for early-stage non-small cell lung cancer. Front. Oncol. 9:693. doi: 10.3389/fonc.2019.00693

PubMed Abstract | CrossRef Full Text | Google Scholar

Hempel, L., Körholz, D., Nussbaum, P., Bönig, H., Burdach, S., and Zintl, F. (1997). High interleukin-10 serum levels are associated with fatal outcome in patients after bone marrow transplantation. Bone Marrow Transplant. 20, 365–368. doi: 10.1038/sj.bmt.1700902

PubMed Abstract | CrossRef Full Text | Google Scholar

Inagaki, Y., Hookway, E., Williams, K. A., Hassan, A. B., Oppermann, U., Tanaka, Y., et al. (2016). Dendritic and mast cell involvement in the inflammatory response to primary malignant bone tumours. Clin. Sarcoma Res. 6:13. doi: 10.1186/s13569-016-0053-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Ingolia, N. T., Lareau, L. F., and Weissman, J. S. (2011). Ribosome profiling of mouse embryonic stem cells reveals the complexity and dynamics of mammalian proteomes. Cell 147, 789–802. doi: 10.1016/j.cell.2011.10.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Iwamoto, Y. (2007). Diagnosis and treatment of Ewing’s sarcoma. Jpn. J. Clin. Oncol. 37, 79–89. doi: 10.1093/jjco/hyl142

PubMed Abstract | CrossRef Full Text | Google Scholar

Javelaud, D., Poupon, M. F., Wietzerbin, J., and Besançon, F. (2002). Inhibition of constitutive NF-kappa B activity suppresses tumorigenicity of Ewing sarcoma EW7 cells. Int. J. Cancer 98, 193–198. doi: 10.1002/ijc.10192

PubMed Abstract | CrossRef Full Text | Google Scholar

Latenstein, A. E. J., van Roessel, S., van der Geest, L. G. M., Bonsing, B. A., Dejong, C. H. C., Groot Koerkamp, B., et al. (2020). Conditional survival after resection for pancreatic cancer: a population-based study and prediction model. Ann. Surg. Oncol. 27, 2516–2524. doi: 10.1245/s10434-020-08235-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Lau, Y. S., Adamopoulos, I. E., Sabokbar, A., Giele, H., Gibbons, C. L., and Athanasou, N. A. (2007). Cellular and humoral mechanisms of osteoclast formation in Ewing’s sarcoma. Br. J. Cancer 96, 1716–1722. doi: 10.1038/sj.bjc.6603774

PubMed Abstract | CrossRef Full Text | Google Scholar

Lissat, A., Joerschke, M., Shinde, D. A., Braunschweig, T., Meier, A., Makowska, A., et al. (2015). IL6 secreted by Ewing sarcoma tumor microenvironment confers anti-apoptotic and cell-disseminating paracrine responses in Ewing sarcoma cells. BMC Cancer 15:552. doi: 10.1186/s12885-015-1564-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, J., Shen, L., Yao, J., Li, Y., Wang, Y., Chen, H., et al. (2015). Forkhead box C1 promoter upstream transcript, a novel long non-coding RNA, regulates proliferation and migration in basal-like breast cancer. Mol. Med. Rep. 11, 3155–3159. doi: 10.3892/mmr.2014.3089

PubMed Abstract | CrossRef Full Text | Google Scholar

Ljunggren, H. G., and Malmberg, K. J. (2007). Prospects for the use of NK cells in immunotherapy of human cancer. Nat. Rev. Immunol. 7, 329–339. doi: 10.1038/nri2073

PubMed Abstract | CrossRef Full Text | Google Scholar

Lollini, P. L., Landuzzi, L., Frabetti, F., Rossi, I., Nicoletti, G., Scotlandi, K., et al. (1998). Expression of functional CD40 on human osteosarcoma and Ewing’s sarcoma cells. Clin. Cancer Res. 4, 1843–1849.

Google Scholar

Ma, Y., Zhang, J., Wen, L., and Lin, A. (2018). Membrane-lipid associated lncRNA: a new regulator in cancer signaling. Cancer Lett. 419, 27–29. doi: 10.1016/j.canlet.2018.01.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Marques Howarth, M., Simpson, D., Ngok, S. P., Nieves, B., Chen, R., Siprashvili, Z., et al. (2014). Long noncoding RNA EWSAT1-mediated gene repression facilitates Ewing sarcoma oncogenesis. J. Clin. Invest. 124, 5275–5290. doi: 10.1172/jci72124

PubMed Abstract | CrossRef Full Text | Google Scholar

Mowel, W. K., Kotzin, J. J., McCright, S. J., Neal, V. D., and Henao-Mejia, J. (2018). Control of immune cell homeostasis and function by lncRNAs. Trends Immunol. 39, 55–69. doi: 10.1016/j.it.2017.08.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Newman, A. M., Steen, C. B., Liu, C. L., Gentles, A. J., Chaudhuri, A. A., Scherer, F., et al. (2019). Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat. Biotechnol. 37, 773–782. doi: 10.1038/s41587-019-0114-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Pahl, J. H., Santos, S. J., Kuijjer, M. L., Boerman, G. H., Sand, L. G., Szuhai, K., et al. (2015). Expression of the immune regulation antigen CD70 in osteosarcoma. Cancer Cell Int. 15:31. doi: 10.1186/s12935-015-0181-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Pan, F., Yao, J., Chen, Y., Zhou, C., Geng, P., Mao, H., et al. (2014). A novel long non-coding RNA FOXCUT and mRNA FOXC1 pair promote progression and predict poor prognosis in esophageal squamous cell carcinoma. Int. J. Clin. Exp. Pathol. 7, 2838–2849.

Google Scholar

Pappo, A. S., and Dirksen, U. (2018). Rhabdomyosarcoma, Ewing Sarcoma, and other round cell sarcomas. J. Clin. Oncol. 36, 168–179. doi: 10.1200/jco.2017.74.7402

PubMed Abstract | CrossRef Full Text | Google Scholar

Peng, W. X., Koirala, P., and Mo, Y. Y. (2017). LncRNA-mediated regulation of cell signaling in cancer. Oncogene 36, 5661–5667. doi: 10.1038/onc.2017.184

PubMed Abstract | CrossRef Full Text | Google Scholar

Picarda, E., Ohaegbulam, K. C., and Zang, X. (2016). Molecular pathways: targeting B7-H3 (CD276) for human cancer immunotherapy. Clin. Cancer Res. 22, 3425–3431. doi: 10.1158/1078-0432.Ccr-15-2428

PubMed Abstract | CrossRef Full Text | Google Scholar

Pollard, J. W. (2004). Tumour-educated macrophages promote tumour progression and metastasis. Nat. Rev. Cancer 4, 71–78. doi: 10.1038/nrc1256

PubMed Abstract | CrossRef Full Text | Google Scholar

Rinn, J. L., Kertesz, M., Wang, J. K., Squazzo, S. L., Xu, X., Brugmann, S. A., et al. (2007). Functional demarcation of active and silent chromatin domains in human HOX loci by noncoding RNAs. Cell 129, 1311–1323. doi: 10.1016/j.cell.2007.05.022

PubMed Abstract | CrossRef Full Text | Google Scholar

Savola, S., Klami, A., Myllykangas, S., Manara, C., Scotlandi, K., Picci, P., et al. (2011). High expression of complement component 5 (C5) at tumor site associates with superior survival in Ewing’s sarcoma family of tumour patients. ISRN Oncol. 2011:168712. doi: 10.5402/2011/168712

PubMed Abstract | CrossRef Full Text | Google Scholar

Subramanian, A., Tamayo, P., Mootha, V. K., Mukherjee, S., Ebert, B. L., Gillette, M. A., et al. (2005). Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl. Acad. Sci. USA 102, 15545–15550. doi: 10.1073/pnas.0506580102

PubMed Abstract | CrossRef Full Text | Google Scholar

Sveen, A., Ågesen, T. H., Nesbakken, A., Meling, G. I., Rognum, T. O., Liestøl, K., et al. (2012). ColoGuidePro: a prognostic 7-gene expression signature for stage III colorectal cancer patients. Clin. Cancer Res. 18, 6001–6010. doi: 10.1158/1078-0432.Ccr-11-3302

PubMed Abstract | CrossRef Full Text | Google Scholar

Thiel, U., Pirson, S., Müller-Spahn, C., Conrad, H., Busch, D. H., Bernhard, H., et al. (2011). Specific recognition and inhibition of Ewing tumour growth by antigen-specific allo-restricted cytotoxic T cells. Br. J. Cancer 104, 948–956. doi: 10.1038/bjc.2011.54

PubMed Abstract | CrossRef Full Text | Google Scholar

Verhoeven, D. H., de Hooge, A. S., Mooiman, E. C., Santos, S. J., ten Dam, M. M., Gelderblom, H., et al. (2008). NK cells recognize and lyse Ewing sarcoma cells through NKG2D and DNAM-1 receptor dependent pathways. Mol. Immunol. 45, 3917–3925. doi: 10.1016/j.molimm.2008.06.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Wong, J. C., Lee, S. B., Bell, M. D., Reynolds, P. A., Fiore, E., Stamenkovic, I., et al. (2002). Induction of the interleukin-2/15 receptor beta-chain by the EWS-WT1 translocation product. Oncogene 21, 2009–2019. doi: 10.1038/sj.onc.1205262

PubMed Abstract | CrossRef Full Text | Google Scholar

Xu, Y. Z., Chen, F. F., Zhang, Y., Zhao, Q. F., Guan, X. L., Wang, H. Y., et al. (2017). The long noncoding RNA FOXCUT promotes proliferation and migration by targeting FOXC1 in nasopharyngeal carcinoma. Tumour Biol. 39:1010428317706054. doi: 10.1177/1010428317706054

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, D., He, Y., Wu, B., Liu, R., Wang, N., Wang, T., et al. (2020). Predictions of the dysregulated competing endogenous RNA signature involved in the progression of human lung adenocarcinoma. Cancer Biomark. 29, 399–416. doi: 10.3233/cbm-200133

PubMed Abstract | CrossRef Full Text | Google Scholar

Yoshihara, K., Shahmoradgoli, M., Martínez, E., Vegesna, R., Kim, H., Torres-Garcia, W., et al. (2013). Inferring tumour purity and stromal and immune cell admixture from expression data. Nat. Commun. 4:2612. doi: 10.1038/ncomms3612

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, X., Yi, S., Xing, G., Wu, H., Zhu, Y., Guo, X., et al. (2020). FOXCUT promotes the proliferation and invasion by activating FOXC1/PI3K/AKT pathway in colorectal cancer. Cancer Manag. Res. 12, 6269–6278. doi: 10.2147/cmar.S259801

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Y., Li, Z., Chen, M., Chen, H., Zhong, Q., Liang, L., et al. (2020). lncRNA TCL6 correlates with immune cell infiltration and indicates worse survival in breast cancer. Breast Cancer 27, 573–585. doi: 10.1007/s12282-020-01048-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Y., Meng, W., and Cui, H. (2018). LncRNA CBR3-AS1 predicts unfavorable prognosis and promotes tumorigenesis in osteosarcoma. Biomed. Pharmacother. 102, 169–174. doi: 10.1016/j.biopha.2018.02.081

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhao, D. L., and Shen, G. (2019). Verification of expressions of lncRNA FOXCUT in gastric adenocarcinoma patients and its effects on cell biological function based on TCGA database. Eur. Rev. Med. Pharmacol. Sci. 23, 6139–6147. doi: 10.26355/eurrev_201907_18427

CrossRef Full Text | Google Scholar

Zhou, M., Zhang, Z., Zhao, H., Bao, S., Cheng, L., and Sun, J. (2018). An immune-related six-lncRNA signature to improve prognosis prediction of glioblastoma multiforme. Mol. Neurobiol. 55, 3684–3697. doi: 10.1007/s12035-017-0572-9

PubMed Abstract | CrossRef Full Text | Google Scholar