The Value of Immune-Related Genes Signature in Osteosarcoma Based on Weighted Gene Co-expression Network Analysis

Wang, Xin; Gan, Li; Ye, Ju; Tang, Mengjie; Liu, Wei

doi:https://doi.org/10.1155/2021/9989321

Journal of Immunology Research

On this page

Abstract Introduction Materials and Methods Results Discussion Conclusions Data Availability Ethical Approval Conflicts of Interest Authors’ Contributions Acknowledgments Supplementary Materials References Copyright Related Articles

Special Issue

The Crosstalk Between Circadian Clocks, Immunity, and Tumour Genesis

View this Special Issue

Research Article | Open Access

Volume 2021 | Article ID 9989321 | https://doi.org/10.1155/2021/9989321

The Value of Immune-Related Genes Signature in Osteosarcoma Based on Weighted Gene Co-expression Network Analysis

Xin Wang,¹Li Gan,²Ju Ye,³Mengjie Tang,⁴and Wei Liu⁵

Academic Editor: Jialiang Liang

Received16 Mar 2021

Accepted25 Apr 2021

Published15 May 2021

Abstract

Background. Osteosarcoma (OS) is a serious malignant tumor that is more common in adolescents or children under 20 years of age. This study is aimed at obtaining immune-related genes (IRGs) associated with the progression and prognosis of OS. Method. Expression profiling data and clinical data for OS were downloaded from the Therapeutically Applicable Research to Generate Effective Treatments (TARGET) database. ESTIMATE calculates immune scores and stromal scores of samples and performs the prognostic analysis. Weighted gene coexpression network analysis (WGCNA) was used to find modules correlated with immune and stromal scores. Cox regression analysis and least absolute shrinkage and selection operator (LASSO) analysis were used to explore IRGs associated with OS prognosis and construct and validate a hazard score model. Finally, we verified the expression and function of EVI2B in OS. Results. WGCNA selected twenty-eight IRGs, 10 of which were associated with OS prognosis, and LASSO further obtained three key prognostic genes. A prognostic model of EVI2B was constructed, and according to the risk score model, patients in the high-risk group had a worse prognosis than those in the low-risk group, and the prognosis was statistically significant in the high- and low-risk groups. Receiver operating characteristic (ROC) curves were used to assess the prognostic model’s accuracy and externally validate the independent GSE21257 cohort. The results of immunohistochemical staining and qPCR showed that EVI2B was a tumor suppressor gene. The differential genes in the high- and low-risk groups were analyzed by enrichment analysis of GO and KEGG, indicating that the EVI2B model is associated with immune response. Conclusion. In this study, IRG EVI2B is closely related to OS’s prognosis and can be used as a potential biomarker for prognosis and treatment of OS.

1. Introduction

Osteosarcoma is a malignant tumor that occurs more commonly in adolescents or children under 20 years of age [1]. Conventional therapies for osteosarcoma include surgery, adjuvant chemotherapy, and neoadjuvant chemotherapy [2–4]. For patients with osteosarcoma, the use of standard multiagent chemotherapy in combination with surgical resection produces a long-term survival of approximately 70% for the localized disease at diagnosis [5, 6] and 20–30% for metastatic disease at diagnosis or recurrence [7]. Although advances have been made in surgical techniques, targeted therapy and tumor immunity, and complications such as infection and poor survival due to limb salvage surgery, there are still many issues to be addressed in the treatment of OS; therefore, there is an urgent need to develop new prediction methods to improve OS patients’ survival [5, 8].

A large body of evidence suggests that OS has an immune system with multiple therapeutic targets, including receptor T cell therapy [9], HER2-specific immunity [10–12], and adjuvant immune therapy with mifamurtide [13]. Immune checkpoint inhibitors may develop immune tolerance to tumor immunosuppressants, thereby increasing endogenous antitumor activity, and could improve therapy-induced tumor immunogenic chemotherapy, radiotherapy, and targeted therapy [13–17]. To date, the impact of IRGs on OS prognosis is not well clear, so prognostic biomarkers with IRG are clinically essential for OS prognosis identification.

A standard analysis method relates the measurement of these genomic covariates to patients’ survival time, which often censored survival data [18]. A popular strategy is to use these covariates to fit a Cox regression model to the censored survival data and then predict new cancer prognosis based on this fitted model [18–20]. Compared with Cox risk regression analysis, data overfitting can be solved entirely by the LASSO analysis [21]. LASSO variants are a popular strategy to provide variable selection in regression analysis and have extended to the Cox regression model [18, 22, 23]. Regression with LASSO penalty is a commonly used method for selection in a high-dimensional variable. Still, the results depend heavily on the value of the shrinkage setting [23]. Systems biology algorithms for WGCNA have been used to assess the association between gene sets and clinical features by constructing scale-free gene coexpression networks [24–27].

This study focused on revealing IRGs involved in OS prognosis. The WGCNA analysis and PPI network analysis were first used by constructing a Cox proportional hazards (PH) prognostic model for LASSO [18, 19, 22]. Univariate and multivariate Cox regression analyses were performed for IRGs in PPI, and the hub genes identified could be used for prognostic risk assessment in OS patients. Finally, we validated the prediction model’s performance and accuracy by ROC analysis. Furthermore, the prognostic significance of this model was validated on an independent GEO dataset. To reveal the role of differential genes in these training groups, functional and pathway enrichment analysis was performed using differential genes. The results showed that these differential genes were mainly associated with immune response and inflammatory response. We also performed immunohistochemical staining and qPCR experimental validation, and the results were consistent with the bioinformatics results. In this study, IRG is involved in the immune process to affect the prognosis of OS patients, thereby intervening in patients’ immune response to improve the poor prognosis and enhance the quality of life of cancer patients.

2. Materials and Methods

2.1. Data Collection and Processing

Gene expression data (FPKM, fragments per kilobase million) and clinical data for OS obtained from the TARGET database (https://ocg.cancer.gov/programs/target). Log₂ (FPKM+1) conversion for the expression value of the TARGET database. Clinicopathological data of the corresponding patients retrieved from the database, including gender, ethnicity, age, tumor location, time to recurrence, and survival information. In this study, data from the GEO database were used (http://www.ncbi.nlm.nih.gov/geo/) of an independent dataset for external validation, downloaded, and collected a high-throughput gene expression dataset, numbered GSE21257 (series matrix file). The Ensemble IDs of the mRNAs in the TARGET database were extracted from the GENCODE project. The ID conversion for the GSE21257 data set is provided by GPL10295 (Illumina human-6 v2.0 expression beadchip).

2.2. Correlation between Prognosis and Stromal/Immune Score

The ESTIMATE algorithm was applied to the normalized expression stromal [28] to estimate each osteosarcoma sample’s stromal and immune score. The overall survival deeds the primary prognostic endpoint, the immune score, and the stromal score of each sample calculated with the R package ESTIMATE; it analyzed overall prognosis survival using the R package survival.

2.3. Construction of Coexpression Network

TARGET’s data set used as the training set, the WGCNA program package [29, 30] in R software used to remove outlier samples and find the modules related to immune score and stromal score toolkit of heat map drawn to analyze the strength of interacting. Module-trait correlations were estimated using correlations between module signature genes and traits (immune score and stromal score). The soft-threshold power of calculates using a scale-free topology criterion. A weighted adjacency stromal generates, and the soft-threshold power can emphasize weak and robust correlations between penalized genes [24, 31]. Modules of RNA were obtained using the dynamic tree cutting method. It is classifying genes with similar expression profiles by modules; the average linkage was performed to perform hierarchical clustering by TOM-based dissimilarity [24]. Finally, the difference of module eigengenes (MES) of the module dendrogram was calculated, and some modules were merged. The pink module with the highest correlation coefficient was selected as the next research object. The conditions for the screening of IRGs in the pink module were that the correlation between genes and the pink module was more significant than 0.8, the correlation coefficient between genes and the immune score was greater than 0.5, and a total of 28 IRGs were selected.

2.4. PPI Network Analysis

Using search tools for reciprocal gene/protein retrieval (https://string-db.org/) database (Version 10.0) [32], interactions with a were considered statistically significant. PPI network results were visualized using Cytoscape (3.4.0); the most important modules were identified by Cytoscape’s [33] Plugin Molecular Complex Detection (MCODE) (version 1.4.2).

2.5. Construction of Risk Assessment Model

Then, we performed univariate Cox regression analysis of 28 IRGs followed by screening for IRGs with significant prognostic differences (); IRGs were further analyzed by LASSO to reduce genes. Finally, a multivariate Cox prognostic model was built by selecting the IRGs closely related to survival. The sample was divided into two groups based on the median risk score. Besides, we validated the model in the training set. External validation was performed in an independent GEO dataset, GSE21257.

2.6. Establishment and Validation of the Prognostic Risk Scoring System

To generate a risk scoring system for genes, we performed multivariate Cox proportional hazards regression. First, we used the survival package of R to obtain the regression coefficients for each gene. The coefficient of each selected gene (parameter coefficient R) represents the estimated logarithm of the hazard ratio (HR, parameter exp (coefficient R)) in R. Then, a risk score formula was established for all patients. ROC curve (AUC) predictions over time calculated for 1-, 3-, and 5-year survival using the survival ROC package in R. The optimal cut-off point was selected as the maximum sensitivity and specificity. According to the optimal cut-off point, patients were divided into high- and low-risk groups; the survival difference between the two groups were assessed with the R package survival analysis.

2.7. GO and KEGG Enrichment Analysis

Differential genes between TARGET high- and low-risk groups (R for differential analysis wrapped as limma package, , ) were subjected to enrichment analysis. GO analysis is a standard method to define genes and their RNA or protein products to identify high-throughput transcriptome or genomic data [34]. KEGG is a collection of databases dealing with genomes, diseases, biological pathways, drugs, and chemicals [35]. We used the R package ClusterProfiler [36] to enrich differential genes in the training set’s high- and low-risk groups.

2.8. Sample Collection

Osteosarcoma samples were collected from Hunan Cancer Hospital in 2020. OS tissues and adjacent normal tissues from 20 patients were collected, immediately placed in liquid nitrogen, and stored at −80°C. None of the OS patients received antineoplastic therapy. Both the patients and their families in this study were fully informed, and informed consent was obtained from the participants. The Hunan Cancer Hospital ethics committee approved the study.

2.9. Immunohistochemical Staining and Analysis

Immunohistochemistry was performed in 20 OS tissues and adjacent normal tissues. Paraffin-embedded sections were stained to determine the expression level of EVI2B protein. Sections were incubated overnight at 4°C with anti-EVI2B 1 antibody (Invitrogen) at a 1 : 100 dilution. After washing with phosphate-buffered saline (PBS), the slides were incubated with a goat anti-mouse IgG secondary antibody conjugated to fluorescein isothiocyanate (ZSDB-BIO, Beijing, China) for 30 mins. They were washed with PBS and then incubated with an antifade reagent (Invitrogen, Carlsbad, USA). Finally, staining was observed using an Olympus CX41 fluorescence microscope (Olympus, Tokyo, Japan). Positive staining intensity in OS and adjacent normal tissues were analyzed by integrated optical density (IOD) using Image-Pro Plus software (version 6.0; Media Cybernetics, United States). All images were taken using the same microscope and camera set. Image-Pro Plus software was used to calculate the mean IOD (μm²) per stained area (IOD/area) for positive staining. -test was used to analyze the results. A value < 0.05 was considered statistically significant.

2.10. Real-Time Quantitative PCR

Total RNA was isolated from tissues from the above patients using Trizol reagent (TaKaRa, Japan) according to the manufacturer’s instructions. One microgram of RNA was reverse transcribed into cDNA using the Revert Aid First Strand cDNA Synthesis Kit (Thermo, USA). Quantitative RT-PCR was then performed with Pro Taq HS Premix Probe qPCR Kit (Accurate, Hunan, China). The amplification program consisted of one cycle of predenaturation at 95°C for 5 mins, 37 cycles of denaturation at 95°C for 30 s, annealing at 60°C for 30 s, and extension at 72°C for 10 mins. The GAPDH gene was used as an endogenous control gene for normalizing the expression of target genes. Each sample was analyzed in triplicate. The sequences of primers were used for RT-qPCR and annealing temperature (Table S1).

2.11. Statistical Analyses

All statistical analyses were performed using the R software 3.5.0. Statistical significance set at a probability value of . Univariate, LASSO and multivariate Cox regression analyses were applied to predict the overall survival of OS patients. The effect of ROC and calibration curve on the prediction accuracy of the prognostic model was compared. Enrichment analysis of differential genes in the high- and low-risk groups in the training set was performed with the R package cluster profile.

3. Results

3.1. Relationship between Immune Score and Stromal Score for Prognosis in OS

The flowchart of this study design is shown in Figure 1. Overall survival analysis did with R package survival analysis of the immune score’s prognosis and stromal score. The immune score and stromal score were calculated based on the ESTIMATE algorithm can promote the quantification of immune and stromal components in tumors; in this algorithm, the immune and stromal score is calculated by analyzing the specific gene expression characteristics of the immune and stromal cells to predict the infiltration of nontumor cells. Patients with high immune scores had better overall survival than those with low immune scores (Figure 2(a)), log-rank value = 0.001. Patients with a high stromal score had a more favorable prognosis than those with a low stromal score (Figure 2(b)), log-rank value = 0.008. The above results indicated that the immune score and stromal score significantly correlated with OS patients’ prognosis. Patients with high immune scores and the stromal score had a good prognosis.

(a)

(b)

3.2. IRGs of OS Screened with WGCNA and PPI Network Diagram

The sample cluster tree had an outlying sample, and the red line was the cut-off value to filter data (Figure 3(a)). All samples were in clusters after removing one outlying sample based on immune score and stromal score (Figure 3(b)). Sample dendrogram and trait heatmap were drawn according to the immune score and stromal score (Figure 3(b)). The above results were obtained by further analysis of the modules associated with immune score and stromal score by WGCNA. In this study, and connectivity were highest when the power set at 15 (Figures 3(c) and 3(d)).

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

Figure 3

The samples clustering dendrogram, immune score, stromal score, and the determination of soft-thresholding power analyzed by WGCNA. (a) Sample clustering was performed to detect outlying samples based on immune score and stromal score. The red line indicates the cut-off point for data filtering in the data preprocessing step. (b) Clustering based on incorporating immune scores and stromal scores of expression data from osteosarcomas with clinical data and color intensity is proportional to immune and stromal scores. (c) Analysis of the scale-free fit index for various soft-thresholding powers (). (d) Analysis of the mean connectivity for various soft-thresholding powers. (e) Clustering dendrogram of genes based on a dissimilarity measure (1-TOM). (f) Pink coexpression modules based on immune scores yielded scatter plots of gene significance (GS) versus module membership (MM). The screening conditions were the correlation between genes, and the pink module was more significant than 0.8, and the correlation coefficient between genes and the immune score was greater than 0.5. (g) The DEG PPI network module analysis via Cytoscape software. (h) Heatmap of the association with MES and immune score and stromal score of OS. OS: osteosarcoma; TOM: topological overlap stromal; MES: module eigengenes; PPI: protein-protein interaction; DEGs: differentially expressed genes.

Thus, identifies distinct gene coexpression modules in OSs. The cluster dendrogram of the selected genes was clustered with an adjacency stromal, and we obtained a module visualization of the genes (Figure 3(e)). We also received the correlation between modules and traits (Figure 3(h)). The pink module has the highest correlation coefficient with an immune score, with a correlation coefficient of 0.58, , which was statistically significant, indicating that the immune score significantly correlated with clinical traits—the above results were based on the WGCNA analysis performed by the immune system stromal scoring. Therefore, we chose this module as the next research object. We screened the pink module’s IRGs with a correlation between genes and the pink module greater than 0.8 and a correlation coefficient between genes and immune scores greater than 0.5 (Figure 3(f)). With conditions, we screened a total of 28 IRGs PTPRJ, SCARF1, CD300C, ANPEP, TBXAS1, RNASE6, HSD3B7, EVI2B, NCKAP1L, PTK2B, SPI1, GRN, SQOR, GLRX, APOBR, GNA15, LILRB4, LAPTM5, CCR1, FERMT3, CYTH4, LCP1, CD4, CSF1R, PSTPIP1, LAT2, SELPLG, and RAC2 shown in Table 1.

A PPI regulatory network of 28 IRGs was obtained based on the STRING online database (https://string-db.org/). We found that the direct interaction between hub genes was through the PPI network analysis (Figure 3(g)). Nineteen genes (PTPRJ, CD300C, RNASE6, EVI2B, NCKAP1L, PTK2B, SPI1, GNA15, LILRB4, LAPTM5, CCR1, FERMT3, CYTH4, LCP1, CD4, CSF1R, PSTPIP1, SELPLG, and RAC2) were correlated well with immune score.

3.3. Recognize IRG Prognostic Gene and Survival Analysis and ROC Analysis

The 28 IRGs were subjected to univariate Cox analysis, resulting in 10 IRGs with prognostic significance (Table 2) ( value < 0.05). These ten genes, EVI2B, GRN, NCKAP1L, SELPLG, LILRB4, FERMT3, SQOR, CYTH4, GNA15, and RNASE6, were subjected to Lasso regression analysis further to screen prognostic genes (Figures 4(a) and 4(b)) after 1,000 stimuli run through the crossvalidation possibility, the optimal lambda determined, and three genes, EVI2B, GRN, and NCKAP1L, were screened. Next, multivariate Cox analysis was performed to construct a prognostic model based on EVI2B, GRN, and NCKAP1L, resulting in a model containing one IRG EVI2B. EVI2B is highly expressed in normal bone tissue controls and lowly expressed in osteosarcoma (Table 2), which can be regarded as a tumor suppressor gene consistent with the results of a colorectal cancer study by Yuan et al. [37].

(a)

(b)

(c)

(d)

By dividing the sample into high- and low-risk groups according to the model’s median risk score, we found a significant difference in survival between the high- and low-risk groups (log-rank value = 0.016) (Figure 4(c)). ROC analysis showed AUC values of 0.676, 0.697, and 0.708 for ROC curves at 1, 3, and 5 years, respectively (Figure 4(d)). All results showed that the prediction effect of the IRG EVI2B model was moderately accurate.

3.4. Evaluating the Accuracy of the IRG EVI2B Model

The time-dependent ROC curve analysis evaluated the accuracy of the OS prediction model constructed by EVI2B. The IRG EVI2B model evaluates in an independent GEO dataset GSE21257. By dividing GSE21257 samples into high- and low-risk groups according to the TARGET dataset samples’ median risk score in the model, we found a significant difference in survival between the high- and low-risk groups (log-rank value = 0.003) (Figure 5(a)). ROC analysis showed that the AUC values for 1, 3, and 5 years of the ROC curve were 0.663, 0.690, and 0.658, respectively (Figure 5(b)). The survival rate of the low-risk group of the GSE21257 sample in the validation set based on the IRG EVI2B model was higher than that of the high-risk group, and the prolonged survival time was consistent with the results in the training set. The AUC values for the validation set GSE21257 samples were like the AUC results in the training set. The result fully validates the IRG EVI2B prognostic model’s accuracy and illustrates that the IRG EVI2B is feasible as a prognostic indicator in OS patients.

(a)

(b)

3.5. GO and KEGG Enrichment Analysis

Next, to further elucidate the molecular functions and signaling pathways in which differential genes between the TARGET high- and low-risk groups are involved, we performed functional enrichment analysis of all DEGs. GO results: these genes are involved not only in T cell activation, regulation of leukocyte differentiation, regulation of lymphocyte activation, and several other biological processes but also in cytokine binding, cytokine receptor complex activity, MHC protein binding, and molecular functions in other tumors. Some of the encoded proteins are important components of the plasma membrane’s external side, secretory granule membrane, and vacuolar membrane (Figure 6(a)). We further analyzed these genes’ signaling pathways, and KEGG results showed that these genes play a regulatory role in important immune-related pathways (such as the B cell receptor signaling pathway and chemokine signaling pathway) (Figure 6(b)). These results suggest that the OS prognostic model is significantly associated with immune response and inflammatory response, affecting cancer patients’ prognosis by participating in immune processes, indicating the correctness of our analysis results.

(a)

(b)

3.6. IRG EVI2B Expression Is Downregulated in OS Tissues

Finally, immunostaining analysis of IRG EVI2B showed that EVI2B expression was higher in adjacent normal tissues (Figure 7(a)). EVI2B-positive staining was significantly different between the control and OS, and positive staining was low and statistically significant in OS (Figure 7(a)). We examined IRG EVI2B expression in clinical samples from OS patients by qPCR. The qPCR results showed that IRG EVI2B was highly expressed at adjacent normal bone tissue levels (Figure 7(b)). The above results showed that IRG EVI2B was highly expressed in adjacent normal bone tissue (control) and lowly expressed in osteosarcoma, and Yuan et al. included colorectal cancer study [37].

(a)

(b)

4. Discussion

In the past decades, traditional cancer treatments, including surgical treatment and chemotherapy/radiotherapy, have incredibly prolonged osteosarcoma patients’ survival time [38, 39]. Although the 5-year overall survival rate is high [9, 40, 41], the survival rate of patients with metastatic disease (most commonly in the lung parenchyma and distal bone) is meager, 19-30% [39, 42]. Also, chemotherapy/radiotherapy acquired resistance to new antitumor drugs and severe side effects of OS treatment’s drug resistance characteristics after extensive surgical resection of OS treatment reach the platform [39, 43].

The use of targeted therapy [44–46] also benefits OS patients; for example, therapies targeting the unfolded protein response (UPR) pathway may be a feasible approach for treating osteosarcoma. The components targeting the UPR may prove beneficial for patients refractory to conventional chemotherapy [47]. In addition to being a marker of metastatic disease, therapeutic targeting of cathepsin D (CTSD) may be a promising novel approach for osteosarcoma treatment; it may produce a good response in metastatic disease [47]. However, targeted therapy also has some limitations and cannot be suitable for all types of osteosarcoma. In recent years, immunotherapy [48] has dramatically developed in cancer. Immune checkpoint-based therapy has been shown to encouraging the antitumor effect by restoring immune response in the tumor microenvironment [49, 50]. At the basement of immune checkpoint inhibitors, ipilimumab (monoclonal antibody) anticytotoxic T lymphocyte antigen four antibodies (mAb) (CTLA4) and monoclonal antibody death protein 1 (PD1) or PD1 ligand (PDL1) against programmed cells show that immune checkpoints may be immune-tolerant to the tumor [13–17, 51]. A study has shown that for osteosarcoma, the immune checkpoint inhibitor PD-L1 is negatively correlated with prognosis, and PD-1 is negatively correlated with overall survival [13–17, 51]. To date, although immunotherapy [10, 11, 13, 51–53] has used in OS, the impact on OS prognosis is not well understood, so there is an urgent need for a new biomarker that can use to predict the prognosis of OS patients, and our study is to affect the prognosis of cancer patients by participating in the immune process.

At present, the use of EVI2B [54] is still relatively rare. The EVI2B is located in the intron of the neurofibromatosis type 1 (NF1) gene and transcribed in the opposite direction to the NF1 gene [55–58]. Still, their expression is not related, indicating that they are independently regulated [59]. EVI2B is a transmembrane protein [60], while NF1 is a tumor suppressor gene [61]. EVI2B is expressed in many different cell types, including myeloid cells [62]. The EVI2B is also located within EVI2, a common viral integration site found in retroviral induced myeloid tumors [54]. It has postulated that viral integration on EVI2 alters EVI2B expression and that this altered expression predisposes mice to the development of bone marrow tumors [54].

The successful construction of the IRG EVI2B prognostic model and the evaluation and validation of the accuracy of the IRG EVI2B prognostic model with the ROC curve method and in the GEO dataset GSE21257 give us hope to use the IRG EVI2B gene as a new biomarker and therapeutic target for predicting the prognosis of osteosarcoma patients. In the same study, EVI2B was found to correlate with the prognosis of patients significantly negatively with colon cancer (), and EVI2B may be involved in the prognosis of colon cancer patients [37, 63]. EVI2B can similarly be used as a prognosis for osteosarcoma compared to the study by Yang et al. [64].

EVI2B is expressed in various tumors [62] and breast cancer mutations [65]. EVI2B is involved in the differentiation of melanocytes and keratinocytes [66, 67]. In fibroblast-like cells derived from neurofibromas, increased EVI2B mRNA levels were found [66]. EVI2B gene, their related studies indicate that the EVI2B gene is a direct target of CCAAT/enhancer-binding protein (C/EBP). The product of this gene, the transmembrane glycoprotein EVI2B (CD361) [68], shows to be abundantly expressed on the surface of primary hematopoietic cells, reaching the highest expression level in mature granulocytes [69]. They use shRNA-mediated downregulation of EVI2B in human and murine cell lines and primary hematopoietic stem and progenitor cells (HSPC), demonstrating impaired myeloid lineage development and altered progenitor function in EVI2B-depleted cells. Its study suggests that EVI2B is an essential modulator of HSPC and bone marrow differentiation and that low levels of EVI2B may contribute to the differentiation inhibitory profile of acute myelocytic leukemia (AML) [69]. It has shown that EVI2B appears highly expressed in patients with 11q deletion in chronic lymphocytic leukemia (CLL). It speculated that the EVI2B gene might be associated with poor prognosis in patients with CLL 11q deletion [70]. It also studied that EVI2B was transiently expressed on cellular adjuvants and used to directly immunize BALB/c mice to obtain a higher specific anti-EVI2B Ab response in immunized mice than stimulated by other cellular adjuvants, which used to accelerate the development of Ab drugs [71]. Our study also showed that EVI2B is associated with the prognosis of osteosarcoma. These results suggest that it is feasible for us to use EVI2B as a marker and potential target for the diagnosis and treatment of osteosarcoma, but further studies should be on its possible action mechanism.

5. Conclusions

In conclusion, we selected IRG EVI2B by WGCNA analysis and constructed an IRG-based prognostic model, which was validated and experimentally analyzed to obtain that IRG EVI2B can significantly predict the prognosis of osteosarcoma patients.

Data Availability

The publicly available datasets were analyzed in this study; these can be found in the TARGET and the NCBI Gene Expression Omnibus (GSE21257).

Ethical Approval

The Hunan Cancer Hospital ethics committee approved the study.

Conflicts 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.

Authors’ Contributions

Li Gan and Xin Wang contributed equally and are co-first authors to this work.

Acknowledgments

This study was supported by the Scientific Research Project of Hunan Health Commission (Nos. 202102041763 and 20200985), Changsha Municipal Natural Science Foundation (No. kq2014267), and Hunan Cancer Hospital Climb Plan (No. 2020QH001).

Supplementary Materials

Table S1: the sequences of primers used for RT-qPCR and annealing temperature. (Supplementary Materials)

References

S. J. Cotterill, C. M. Wright, M. S. Pearce, and A. W. Craft, “Stature of young people with malignant bone tumors,” Pediatric Blood & Cancer, vol. 42, no. 1, pp. 59–63, 2004.
View at: Publisher Site | Google Scholar
M. S. Isakoff, S. S. Bielack, P. Meltzer, and R. Gorlick, “Osteosarcoma: current treatment and a collaborative pathway to success,” Journal of Clinical Oncology, vol. 33, no. 27, pp. 3029–3035, 2015.
View at: Publisher Site | Google Scholar
A. Misaghi, A. Goldin, M. Awad, and A. A. Kulidjian, “Osteosarcoma: a comprehensive review,” SICOT J, vol. 4, p. 12, 2018.
View at: Publisher Site | Google Scholar
S. S. Bielack, S. Hecker-Nolting, C. Blattmann, and L. Kager, “Advances in the management of osteosarcoma,” F1000Research, vol. 5, 2016.
View at: Publisher Site | Google Scholar
F. Liu, L. Xing, X. Zhang, and X. Zhang, “A four-pseudogene classifier identified by machine learning serves as a novel prognostic marker for survival of osteosarcoma,” Genes (Basel), vol. 10, no. 6, p. 414, 2019.
View at: Publisher Site | Google Scholar
N. M. Bernthal, N. Federman, F. R. Eilber et al., “Long-term results (>25 years) of a randomized, prospective clinical trial evaluating chemotherapy in patients with high-grade, operable osteosarcoma,” Cancer, vol. 118, no. 23, pp. 5888–5893, 2012.
View at: Publisher Site | Google Scholar
C. C. Wu, H. C. Beird, J. Andrew Livingston et al., “Immuno-genomic landscape of osteosarcoma,” Nature Communications, vol. 11, no. 1, 2020.
View at: Google Scholar
A. Luetke, P. A. Meyers, I. Lewis, and H. Juergens, “Osteosarcoma treatment - where do we stand? A state of the art review,” Cancer Treatment Reviews, vol. 40, no. 4, pp. 523–532, 2014.
View at: Publisher Site | Google Scholar
Y. Wang, W. Yu, J. Zhu et al., “Anti-CD166/4-1BB chimeric antigen receptor T cell therapy for the treatment of osteosarcoma,” Journal of Experimental & Clinical Cancer Research, vol. 38, no. 1, 2019.
View at: Publisher Site | Google Scholar
N. Ahmed, V. S. Brawley, M. Hegde et al., “Human epidermal growth factor receptor 2 (HER2) -specific chimeric antigen receptor-modified T cells for the immunotherapy of HER2-positive sarcoma,” Journal of Clinical Oncology, vol. 33, no. 15, pp. 1688–1696, 2015.
View at: Publisher Site | Google Scholar
N. J. Mason, J. S. Gnanandarajah, J. B. Engiles et al., “Immunotherapy with a HER2-targeting listeria induces HER2-specific immunity and demonstrates potential therapeutic effects in a phase I trial in canine osteosarcoma,” Clinical Cancer Research, vol. 22, no. 17, pp. 4380–4390, 2016.
View at: Publisher Site | Google Scholar
N. Rainusso, V. S. Brawley, A. Ghazi et al., “Immunotherapy targeting HER2 with genetically modified T cells eliminates tumor-initiating cells in osteosarcoma,” Cancer Gene Therapy, vol. 19, no. 3, pp. 212–217, 2012.
View at: Publisher Site | Google Scholar
M. Kansara, M. W. Teng, M. J. Smyth, and D. M. Thomas, “Translational biology of osteosarcoma,” Nature Reviews. Cancer, vol. 14, no. 11, pp. 722–735, 2014.
View at: Publisher Site | Google Scholar
P. W. Kantoff, C. S. Higano, N. D. Shore et al., “Sipuleucel-T immunotherapy for castration-resistant prostate cancer,” The New England Journal of Medicine, vol. 363, no. 5, pp. 411–422, 2010.
View at: Publisher Site | Google Scholar
J. R. Brahmer, S. S. Tykodi, L. Q. Chow et al., “Safety and activity of anti-PD-L1 antibody in patients with advanced cancer,” The New England Journal of Medicine, vol. 366, no. 26, pp. 2455–2465, 2012.
View at: Publisher Site | Google Scholar
A. D. Garg and P. Agostinis, “Cell death and immunity in cancer: from danger signals to mimicry of pathogen defense responses,” Immunological Reviews, vol. 280, no. 1, pp. 126–148, 2017.
View at: Publisher Site | Google Scholar
G. Kroemer, L. Galluzzi, O. Kepp, and L. Zitvogel, “Immunogenic cell death in cancer therapy,” Annual Review of Immunology, vol. 31, no. 1, pp. 51–72, 2013.
View at: Publisher Site | Google Scholar
C. Liu, X. Wang, G. Z. Genchev, and H. Lu, “Multi-omics facilitated variable selection in Cox-regression model for cancer prognosis prediction,” Methods, vol. 124, pp. 100–107, 2017.
View at: Publisher Site | Google Scholar
Y. Zhang, H. Li, W. Zhang, Y. Che, W. Bai, and G. Huang, “LASSO‑based Cox‑PH model identifies an 11‑lncRNA signature for prognosis prediction in gastric cancer,” Molecular Medicine Reports, vol. 18, no. 6, pp. 5579–5593, 2018.
View at: Publisher Site | Google Scholar
M. Bagnoli, S. Canevari, D. Califano et al., “Development and validation of a microRNA-based signature (MiROvaR) to predict early relapse or progression of epithelial ovarian cancer: a cohort study,” Lancet Oncology, vol. 17, no. 8, pp. 1137–1146, 2016.
View at: Publisher Site | Google Scholar
Y. Xiong, L. Yuan, J. Xiong et al., “An outcome model for human bladder cancer: a comprehensive study based on weighted gene co-expression network analysis,” Journal of Cellular and Molecular Medicine, vol. 24, no. 3, pp. 2342–2355, 2020.
View at: Publisher Site | Google Scholar
R. Tibshirani, “The lasso method for variable selection in the Cox model,” Statistics in Medicine, vol. 16, no. 4, pp. 385–395, 1997.
View at: Publisher Site | Google Scholar
N. Ternes, F. Rotolo, and S. Michiels, “Empirical extensions of the lasso penalty to reduce the false discovery rate in high-dimensional Cox regression models,” Statistics in Medicine, vol. 35, no. 15, pp. 2561–2573, 2016.
View at: Publisher Site | Google Scholar
P. F. Chen, F. Wang, J. Y. Nie et al., “Co-expression network analysis identified CDH11 in association with progression and prognosis in gastric cancer,” OncoTargets and Therapy, vol. Volume 11, pp. 6425–6436, 2018.
View at: Publisher Site | Google Scholar
P. Langfelder and S. Horvath, “WGCNA: an R package for weighted correlation network analysis,” BMC Bioinformatics, vol. 9, no. 1, 2008.
View at: Publisher Site | Google Scholar
F. Wang, Y. Chang, J. Li et al., “Strong correlation between _ASPM_ gene expression and HCV cirrhosis progression identified by co-expression analysis,” Digestive and Liver Disease, vol. 49, no. 1, pp. 70–76, 2017.
View at: Publisher Site | Google Scholar
Z. Zhou, S. Liu, M. Zhang et al., “Overexpression of topoisomerase 2-alpha confers a poor prognosis in pancreatic adenocarcinoma identified by co-expression analysis,” Digestive Diseases and Sciences, vol. 62, no. 10, pp. 2790–2800, 2017.
View at: Publisher Site | Google Scholar
H. Wang, X. Wu, and Y. Chen, “Stromal-immune score-based gene signature: a prognosis stratification tool in gastric cancer,” Frontiers in Oncology, vol. 9, 2019.
View at: Publisher Site | Google Scholar
M. J. Mason, G. Fan, K. Plath, Q. Zhou, and S. Horvath, “Signed weighted gene co-expression network analysis of transcriptional regulation in murine embryonic stem cells,” BMC Genomics, vol. 10, no. 1, p. 327, 2009.
View at: Publisher Site | Google Scholar
X. Zhai, Q. Xue, Q. Liu, Y. Guo, and Z. Chen, “Colon cancer recurrence associated genes revealed by WGCNA coexpression network analysis,” Molecular Medicine Reports, vol. 16, no. 5, pp. 6499–6505, 2017.
View at: Publisher Site | Google Scholar
Z. Liu, M. Li, Q. Hua, Y. Li, and G. Wang, “Identification of an eight-lncRNA prognostic model for breast cancer using WGCNA network analysis and a Cox‑proportional hazards model based on L1-penalized estimation,” International Journal of Molecular Medicine, vol. 44, no. 4, pp. 1333–1343, 2019.
View at: Publisher Site | Google Scholar
D. Szklarczyk, A. Franceschini, S. Wyder et al., “STRING v10: protein-protein interaction networks, integrated over the tree of life,” Nucleic Acids Research, vol. 43, no. D1, pp. D447–D452, 2015.
View at: Publisher Site | Google Scholar
P. Shannon, A. Markiel, O. Ozier et al., “Cytoscape: a software environment for integrated models of biomolecular interaction networks,” Genome Research, vol. 13, no. 11, pp. 2498–2504, 2003.
View at: Publisher Site | Google Scholar
M. Ashburner, C. A. Ball, J. A. Blake et al., “Gene Ontology: tool for the unification of biology,” Nature Genetics, vol. 25, no. 1, pp. 25–29, 2000.
View at: Publisher Site | Google Scholar
D. W. Huang, B. T. Sherman, and R. A. Lempicki, “Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources,” Nature Protocols, vol. 4, no. 1, pp. 44–57, 2009.
View at: Publisher Site | Google Scholar
G. Yu, L. G. Wang, Y. Han, and Q. Y. He, “clusterProfiler: an R package for comparing biological themes among gene clusters,” OMICS, vol. 16, no. 5, pp. 284–287, 2012.
View at: Publisher Site | Google Scholar
Y. Yuan, J. Chen, J. Wang et al., “Identification hub genes in colorectal cancer by integrating weighted gene co-expression network analysis and clinical validation in vivo and vitro,” Frontiers in Oncology, vol. 10, 2020.
View at: Publisher Site | Google Scholar
P. A. Meyers, C. L. Schwartz, M. Krailo et al., “Osteosarcoma: a randomized, prospective trial of the addition of ifosfamide and/or muramyl tripeptide to cisplatin, doxorubicin, and high-dose methotrexate,” Journal of Clinical Oncology, vol. 23, no. 9, pp. 2004–2011, 2005.
View at: Publisher Site | Google Scholar
M. W. Bishop, K. A. Janeway, and R. Gorlick, “Future directions in the treatment of osteosarcoma,” Current Opinion in Pediatrics, vol. 28, no. 1, pp. 26–33, 2016.
View at: Publisher Site | Google Scholar
N. Ahmed, V. S. Salsman, E. Yvon et al., “Immunotherapy for osteosarcoma: genetic modification of T cells overcomes low levels of tumor antigen expression,” Molecular Therapy, vol. 17, no. 10, pp. 1779–1787, 2009.
View at: Publisher Site | Google Scholar
D. C. Allison, S. C. Carney, E. R. Ahlmann et al., “A meta-analysis of osteosarcoma outcomes in the modern medical era,” Sarcoma, vol. 2012, Article ID 704872, 10 pages, 2012.
View at: Publisher Site | Google Scholar
F. Fagioli, E. Biasin, O. M. Mereuta et al., “Poor prognosis osteosarcoma: new therapeutic approach,” Bone Marrow Transplantation, vol. 41, Supplement 2, pp. S131–S134, 2008.
View at: Publisher Site | Google Scholar
J. Luo, Y. Xia, Y. Yin et al., “ATF4 destabilizes RET through nonclassical GRP78 inhibition to enhance chemosensitivity to bortezomib in human osteosarcoma,” Theranostics., vol. 9, no. 21, pp. 6334–6353, 2019.
View at: Publisher Site | Google Scholar
P. Angulo, G. Kaushik, D. Subramaniam et al., “Natural compounds targeting major cell signaling pathways: a novel paradigm for osteosarcoma therapy,” Journal of Hematology & Oncology, vol. 10, no. 1, p. 10, 2017.
View at: Publisher Site | Google Scholar
L. C. Sayles, M. R. Breese, A. L. Koehne et al., “Genome-informed targeted therapy for osteosarcoma,” Cancer Discovery, vol. 9, no. 1, pp. 46–63, 2019.
View at: Publisher Site | Google Scholar
A. B. Shaikh, F. Li, M. Li et al., “Present advances and future perspectives of molecular targeted therapy for osteosarcoma,” International Journal of Molecular Sciences, vol. 17, no. 4, p. 506, 2016.
View at: Publisher Site | Google Scholar
J. Burns, C. P. Wilding, R. L. Jones, and P. H. Huang, “Proteomic research in sarcomas - current status and future opportunities,” Seminars in Cancer Biology, vol. 61, pp. 56–70, 2020.
View at: Publisher Site | Google Scholar
M. F. Wedekind, L. M. Wagner, and T. P. Cripe, “Immunotherapy for osteosarcoma: where do we go from here?” Pediatric Blood & Cancer, vol. 65, no. 9, p. e27227, 2018.
View at: Publisher Site | Google Scholar
K. M. Mahoney, P. D. Rennert, and G. J. Freeman, “Combination cancer immunotherapy and new immunomodulatory targets,” Nature Reviews. Drug Discovery, vol. 14, no. 8, pp. 561–584, 2015.
View at: Publisher Site | Google Scholar
W. Wu, D. Jing, Z. Meng et al., “FGD1 promotes tumor progression and regulates tumor immune response in osteosarcoma via inhibiting PTEN activity,” Theranostics., vol. 10, no. 6, pp. 2859–2871, 2020.
View at: Publisher Site | Google Scholar
B. Zheng, T. Ren, Y. Huang et al., “PD-1 axis expression in musculoskeletal tumors and antitumor effect of nivolumab in osteosarcoma model of humanized mouse,” Journal of Hematology & Oncology, vol. 11, no. 1, 2018.
View at: Publisher Site | Google Scholar
M. F. Heymann, K. Schiavone, and D. Heymann, “Bone sarcomas in the immunotherapy era,” British Journal of Pharmacology, vol. 178, no. 9, 2021.
View at: Publisher Site | Google Scholar
M. Burgess and H. Tawbi, “Immunotherapeutic approaches to sarcoma,” Current Treatment Options in Oncology, vol. 16, no. 6, 2015.
View at: Publisher Site | Google Scholar
A. M. Buchberg, H. G. Bedigian, N. A. Jenkins, and N. G. Copeland, “Evi-2, a common integration site involved in murine myeloid leukemogenesis,” Molecular and Cellular Biology, vol. 10, no. 9, pp. 4658–4666, 1990.
View at: Publisher Site | Google Scholar
D. E. Jenne, S. Tinschert, H. Reimann et al., “Molecular characterization and gene content of breakpoint boundaries in patients with neurofibromatosis type 1 with 17q11.2 microdeletions,” The American Journal of Human Genetics, vol. 69, no. 3, pp. 516–527, 2001.
View at: Publisher Site | Google Scholar
L. M. Kayes, W. Burke, V. M. Riccardi et al., “Deletions spanning the neurofibromatosis 1 gene: identification and phenotype of five patients,” American Journal of Human Genetics, vol. 54, no. 3, pp. 424–436, 1994.
View at: Google Scholar
H. Kehrer-Sawatzki, C. Maier, E. Moschgath, G. Elgar, and W. Krone, “Genomic characterization of the neurofibromatosis type 1 gene of _Fugu rubripes_,” Gene, vol. 222, no. 1, pp. 145–153, 1998.
View at: Publisher Site | Google Scholar
D. Viskochil, R. Cawthon, P. O'Connell et al., “The gene encoding the oligodendrocyte-myelin glycoprotein is embedded within the neurofibromatosis type 1 gene,” Molecular and Cellular Biology, vol. 11, no. 2, pp. 906–912, 1991.
View at: Publisher Site | Google Scholar
G. Serra, V. Antona, G. Corsello, F. Zara, E. Piro, and R. Falsaperla, “NF1 microdeletion syndrome: case report of two new patients,” Italian Journal of Pediatrics, vol. 45, no. 1, p. 138, 2019.
View at: Publisher Site | Google Scholar
R. M. Cawthon, L. B. Andersen, A. M. Buchberg et al., “cDNA sequence and genomic structure of _EVI2B_ , a gene lying within an intron of the neurofibromatosis type 1 gene,” Genomics, vol. 9, no. 3, pp. 446–460, 1991.
View at: Publisher Site | Google Scholar
A. Pemov, NISC Comparative Sequencing Program, H. Li et al., “The primacy of _NF1_ loss as the driver of tumorigenesis in neurofibromatosis type 1-associated plexiform neurofibromas,” Oncogene, vol. 36, no. 22, pp. 3168–3177, 2017.
View at: Publisher Site | Google Scholar
B. Andreopoulos and D. Anastassiou, “Integrated analysis reveals hsa-miR-142 as a representative of a lymphocyte-specific gene expression and methylation signature,” Cancer Information, vol. 11, pp. 61–75, 2012.
View at: Publisher Site | Google Scholar
M. Y. Huang, H. M. Wang, T. S. Tok et al., “EVI2B, ATP2A2, S100B, TM4SF3, and OLFM4 as potential prognostic markers for postoperative Taiwanese colorectal cancer patients,” DNA and Cell Biology, vol. 31, no. 4, pp. 625–635, 2012.
View at: Publisher Site | Google Scholar
B. Yang, Z. Su, G. Chen et al., “Identification of prognostic biomarkers associated with metastasis and immune infiltration in osteosarcoma,” Oncology Letters, vol. 21, no. 3, p. 180, 2021.
View at: Publisher Site | Google Scholar
E. Y. Leung, M. E. Askarian-Amiri, D. C. Singleton et al., “Derivation of breast cancer cell lines under physiological (5%) oxygen concentrations,” Frontiers in Oncology, vol. 8, 2018.
View at: Publisher Site | Google Scholar
D. Kaufmann, S. Gruener, F. Braun et al., “EVI2B, a gene lying in an intron of the neurofibromatosis type 1 (NF1) gene, is as the NF1 gene involved in differentiation of melanocytes and keratinocytes and is overexpressed in cells derived from NF1 neurofibromas,” DNA and Cell Biology, vol. 18, no. 5, pp. 345–356, 1999.
View at: Publisher Site | Google Scholar
J. Voisey, G. Kelly, and A. Van Daal, “Agouti signal protein regulation in human melanoma cells,” Pigment Cell Research, vol. 16, no. 1, pp. 65–71, 2003.
View at: Publisher Site | Google Scholar
J. Matesanz-Isabel, J. Sintes, L. Llinas, J. de Salort, A. Lazaro, and P. Engel, “New B-cell CD molecules,” Immunology Letters, vol. 134, no. 2, pp. 104–112, 2011.
View at: Publisher Site | Google Scholar
P. Zjablovskaja, M. Kardosova, P. Danek et al., “Correction to: EVI2B is a C/EBP _α_ target gene required for granulocytic differentiation and functionality of hematopoietic progenitors,” Cell Death & Differentiation, vol. 26, no. 1, p. 198, 2019.
View at: Publisher Site | Google Scholar
Y. Aalto, W. El-Rifa, L. Vilpo et al., “Distinct gene expression profiling in chronic lymphocytic leukemia with 11q23 deletion,” Leukemia, vol. 15, no. 11, pp. 1721–1728, 2001.
View at: Publisher Site | Google Scholar
C. C. Huang, K. W. Cheng, Y. C. Hsieh et al., “Use of syngeneic cells expressing membrane-bound GM-CSF as an adjuvant to induce antibodies against native multi-pass transmembrane protein,” Scientific Reports, vol. 9, no. 1, p. 9931, 2019.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2021 Xin Wang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

456

Downloads

729

Citations