Abstract
Osteosarcoma is one of the most aggressive malignant bone tumors worldwide. Although great advancements have been made in its treatment owing to the advent of neoadjuvant chemotherapy, the problem of lung metastasis is a major obstacle in the improvement of survival outcomes. Thus, the aim of the present study is to screen novel and key biomarkers, which may act as potential prognostic markers and therapeutic targets in osteosarcoma. We utilized the robust rank aggregation (RRA) method to integrate three osteosarcoma microarray datasets downloaded from the Gene Expression Omnibus (GEO) database, and we identified the robust differentially expressed genes (DEGs) between primary and metastatic osteosarcoma tissues. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were performed to explore the functions of robust DEGs. The results of enrichment analysis showed that the robust DEGs were closely associated with osteosarcoma development and progression. Immune cell infiltration analysis was also conducted by CIBERSORT algorithm, and we found that macrophages are the most principal infiltrating immune cells in osteosarcoma, especially macrophages M0 and M2. Then, the protein–protein interaction network and key modules were constructed by Cytoscape, and 10 hub genes were selected by plugin cytoHubba from the whole network. The survival analysis of hub genes was also carried out based on the Therapeutically Applicable Research to Generate Effective Treatments (TARGET) database. The integrated bioinformatics analysis was utilized to provide new insight into osteosarcoma development and metastasis and identified EGR1, CXCL10, MYC, and CXCR4 as potential biomarkers for prognosis of osteosarcoma.
Introduction
Osteosarcoma is one of the most aggressive malignant tumors in the bone (), which derived from mesenchymal tissue and showed osteoblastic differentiation (). Annually, the incidence rate of osteosarcoma is approximately four to five cases per million (). In addition, osteosarcoma is most prevalent in children and adolescents (), and 15–20% of osteosarcoma patients have lung metastasis at the initial diagnosis (, ). With the assistance of neoadjuvant chemotherapy, the treatment of osteosarcoma has been greatly improved, but the overall survival of patients with lung metastasis or relapse has not improved and remains low at approximately 20% (, ). Therefore, it is extremely necessary to seek novel prognostic factors and therapeutic targets for osteosarcoma.
In recent years, public databases including Gene Expression Omnibus (GEO), The Cancer Genome Atlas (TCGA), and Therapeutically Applicable Research to Generate Effective Treatments (TARGET) are widely used to explore diagnostic and prognostic biomarkers in osteosarcoma. In the previous studies, the limited number of samples and inappropriate analysis methods of multiple datasets led to deviation of the results. Zhang et al. () and Diao et al. () described the role of gene copy number alterations and methylation changes in the malignant progression of osteosarcoma using bioinformatics analyses, respectively. However, these studies were only based on single datasets and had a limited sample size, which may have biased the final results. To analyze more samples and avoid the sample heterogeneity of each independent experiment and the error caused by different technology platforms and different data processing methods, we used the robust rank aggregation (RRA) method to obtain robust differentially expressed genes (DEGs). RRA was used to compare the ranking of multiple gene lists. If a gene ranked the highest in all gene lists, then the smaller its calculated P-value is, the more likely it is to be a robust DEG (). This method has been widely used in integrated analysis of multiple datasets, and it is robust to errors and noise (–). There have been no reports of the use of RRA in osteosarcoma.
In our study, microarray datasets GSE14827, GSE21257, and GSE32981 from GEO database were downloaded and analyzed by RRA method to identify robust DEGs between primary and metastatic osteosarcoma tissues. A total of 524 robust DEGs were determined, including 272 upregulated and 252 downregulated genes. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were performed to explore the functions of robust DEGs. Immune cell infiltration analysis was also conducted by CIBERSORT algorithm. In addition, we also constructed the protein–protein interaction (PPI) network and key modules, and finally 10 hub genes were selected by plugin cytoHubba from the whole network. Survival analysis of hub genes was carried out using R packages. In conclusion, the integrated bioinformatics analysis was utilized to identify the significant robust DEGs and hub genes, which may act as novel and potential prognostic biomarkers in osteosarcoma.
Materials and Methods
Data Collection and Data Processing
We selected three gene chips of osteosarcoma from GEO database1, including GSE14827, GSE21257, and GSE32981. The selection criteria were as follows: (1) inclusion of primary and metastatic osteosarcoma tissue samples; (2) expression profiling by array as the experiment type; (3) Homo sapiens; and (4) 20 samples as the minimal size. Among them, GSE14827 contained 18 primary samples and 9 metastatic osteosarcoma samples. GSE21257 contained 19 non-metastatic and 34 metastatic osteosarcoma samples. The GSE32981 dataset contained 23 samples. One sample was not available and was excluded. Twenty-two samples were selected for further study including four primary and 18 metastatic samples. The matrix files and platform annotation document of three microarray datasets were downloaded. The names of microarray probes were converted to the gene symbols by Perl. The DEGs were identified between primary and metastatic osteosarcoma samples in each dataset by limma package in R () with the cutoff criteria of |log2 fold change (FC)| > 0.585 and P-value < 0.05.
Robust Rank Aggregation Analysis
To integrate the three microarray datasets, RRA method was used to determine the robust DEGs (), which is a standard method to minimize the bias and errors among several datasets. Before RRA analysis, the upregulated and downregulated genes were ranked by their FC in each dataset. Then, the RobustRankAggreg R package was performed to get robust DEGs on the basis of the ranked genes in the three datasets. Genes with FC > 1.3 and P-value < 0.05 were considered as the significant robust DEGs.
Gene Ontology Function and Kyoto Encyclopedia of Genes and Genomes Pathway Enrichment Analyses
To identify the functional roles of the robust DEGs indicated above, GO enrichment results of biological process (BP), cellular component (CC), and molecular function (MF) were obtained using the R package “clusterprofiler.” The KEGG pathway analysis of robust DEGs was also conducted using the R package (). P < 0.05 was considered statistically significant.
Immune Infiltration by CIBERSORT Analysis
The CIBERSORT algorithm is commonly used to predict the infiltration of 22 types of immune cells in each tissue sample (). The 22 kinds of immune cells include seven types of T cells [CD8+ T cells, naĂŻve CD4+ T cells, resting memory CD4+ T cells, activated memory CD4+ T cells, follicular helper T cells, regulatory T cells (Tregs), gamma delta T cells], three types of macrophages (M0, M1, and M2), naĂŻve B cells, memory B cells, plasma cells, resting natural killer (NK) cells, activated NK cells, monocytes, resting dendritic cells, activated dendritic cells, resting mast cells, activated mast cells, eosinophils, and neutrophils. Normalized gene expression matrix was converted to 22 types of immune cell matrix by the CIBERSORT algorithm. The immune cell matrix was filtered according to the criteria of P < 0.05, and then the relative expression of 22 types of immune cells was identified between primary and metastatic osteosarcoma samples by R packages. The principal component analysis (PCA) was also performed to determine the difference between primary and metastatic samples.
Protein–Protein Interaction Network Construction and Module Analysis
We uploaded the robust DEGs to the STRING online database2, and we chose confidence >0.9 as the screening criteria. The visualized PPI network was performed by Cytoscape (version 3.6.1) software3. Cytoscape plugin-MCODE was used to screen the significant modules in the PPI network.
Hub Gene Identification
cytoHubba, a plugin of Cytoscape, provides several topological analysis algorithms, including Degree, Edge Percolated Component (EPC), Maximum Neighborhood Component (MNC), Density of Maximum Neighborhood Component (DMNC), Maximal Clique Centrality (MCC), and six centralities that include BottleNeck, EcCentricity, Closeness, Radiality, Betweenness, and Stress. These algorithms can be used to identify hub genes ().
Survival Analysis
The RNA-seq-FPKM data and prognostic information of osteosarcoma patients were downloaded from TARGET database4. TARGET is a database that only includes children’s tumors. Presently, TARGET database contains six kinds of tumors, including ALL (Acute Lymphoblastic Leukemia), AML (Acute Myeloid Leukemia), KT (Kidney Tumors), MDLS (Model Systems), NBL (Neuroblastoma), and OS (Osteosarcoma). The survival analyses of hub genes were conducted by R package survival and survminer. P < 0.05 was considered to be statistically significant.
Results
Identification of Differentially Expressed Genes in Each Dataset
In the present study, the biological characteristics of DEGs were identified by integrated bioinformatics analysis. The overall workflow of this study is showed in Figure 1. The osteosarcoma microarray data GSE14827, GSE21257, and GSE32981 were selected and analyzed using limma package in R. A total of 102 osteosarcoma samples including 41 primary and 61 metastatic tissues were identified in our study. According to the cutoff criteria of | log2 FC| > 0.585 and P < 0.05, there were 83 DEGs in GSE14827 including 23 upregulated and 60 downregulated genes. A total of 464 DEGs were screened from the GSE21257 dataset, including 247 upregulated and 217 downregulated genes. Additionally, a total of 650 DEGs were selected in GSE32981, including 299 upregulated and 351 downregulated genes. The distribution of DEGs is shown in volcano plots (Figures 2A–C), where the red and green dots represent the upregulated and downregulated genes, respectively.
FIGURE 1
FIGURE 2
Identification of Robust Differentially Expressed Genes by Robust Rank Aggregation Method
To integrate the three datasets with minimal bias, the RRA method was used. A total of 524 robust DEGs were determined, including 272 upregulated and 252 downregulated genes (Supplementary Table S1). According to the P-value of robust DEGs, we assigned the top 20 upregulated and downregulated robust DEGs in the visualized heatmap (Figure 2D).
Functional Enrichment Analyses of Robust Differentially Expressed Genes
To explore the functions of robust DEGs, the GO and KEGG enrichment analyses were conducted by R packages. The results of GO analysis included three categories: BP, CC, and MF. For BP, the upregulated robust DEGs were mainly enriched in embryonic organ development, multicellular organismal homeostasis, and transmembrane receptor protein serine/threonine kinase signaling pathway. In the CC part, the upregulated genes were particularly enriched in lamellar body, cell–cell junction, and cell–cell adherens junction. The top three significantly enriched terms were DNA-binding transcription activator activity, growth factor receptor binding, and cell adhesion molecule binding in the MF group (Figure 3A and Supplementary Table S2). Moreover, the most significantly enriched GO BP terms of downregulated genes were extracellular matrix (ECM) organization, extracellular structure organization, and regulation of cellular response to growth factor stimulus. For CC category, the downregulated genes were enriched in collagen-containing ECM. In addition, the downregulated DEGs were mainly enriched in ECM structural constituent conferring compression resistance, ECM structural constituent, and glycosaminoglycan binding in the MF group (Figure 3B and Supplementary Table S3). The above results indicated that the robust DEGs were mostly associated with cancer-related functions.
FIGURE 3
The result of KEGG pathway enrichment analysis is also shown in Figure 3C. TGF-beta signaling pathway, wnt signaling pathway, and IL-17 signaling pathway were highly associated with tumor progression.
Immune Cell Infiltration Analysis
With the use of CIBERSORT algorithm, the infiltration of 22 kinds of immune cells in 102 osteosarcoma tissues is shown in Figure 4A. There was no significant difference in the infiltration of immune cells between primary and metastatic osteosarcoma tissues. However, compared with other immune cells, such as T cells and B cells, macrophage infiltration dominated, whether in primary or metastatic osteosarcoma tissues (Figure 4B). The above results demonstrated that macrophages may play an important role in the development and progression of osteosarcoma. The visualized violin plot was also constructed to prove the above findings (Figure 4C). The PCA of Figure 4D showed nothing individual difference between primary and metastasis samples.
FIGURE 4
Protein–Protein Interaction Network Construction and Module Analysis
To further study the interaction of robust DEGs, we constructed the PPI network by STRING database. With the confidence >0.9 and hiding the disconnected nodes, a visualized PPI network was created by Cytoscape (Figure 5A). In the final network, there were 148 nodes and 302 edges, including 84 upregulated and 64 downregulated genes. By using MCODE plugin, three key modules were screened from the whole network (Figures 5B–D). The robust DEGs in module 1 were mainly enriched in type I interferon signaling pathway. BP of genes in module 2 was particularly enriched in chemokine-mediated signaling pathway. In addition, genes in module 3 were mainly enriched in cell–cell adhesion mediated by cadherin, cell–cell junction assembly, and adherens junction organization (Supplementary Table S4).
FIGURE 5
Hub Gene Identification
cytoHubba is a Cytoscape plugin that allows the use several topological analysis algorithms, including MCC, DMNC, MNC, Degree, EPC, BottleNeck, EcCentricity, Closeness, Radiality, Betweenness, and Stress. These approaches can be used to predict and explore important nodes in PPI networks. Scores from topological algorithms are assigned to each node in a PPI network. According to the gene score, the top ranked genes can be considered as the hub genes. In the present study, we used 10 topological analysis algorithms (MCC, DMNC, MNC, Degree, EPC, BottleNeck, EcCentricity, Closeness, Radiality, and Betweenness) to rank the top 50 genes of the whole network. The intersection of these 50 genes from the 10 algorithms revealed the 10 hub genes: POMC, EGR1, CXCL10, SERPINA1, OAS1, MYC, CXCR4, CXCL2, CHRDL1, and GNAI1 (Figure 5E). The description of the 10 hub genes is shown in Table 1, including full names, synonyms and primary functions.
TABLE 1
| Gene | Full name | Synonyms | Function |
| POMC | Pro-opiomelanocortin | Regulation of cytokine-mediated pathway and signal transduction | |
| EGR1 | Early growth response protein 1 | KROX24, ZNF225 | Transcriptional regulator. Regulation of cell survival, proliferation and cell death |
| CXCL10 | C-X-C motif chemokine 10 | INP10, SCYB10 | Pro-inflammatory cytokine that is involved in a wide variety of processes such as chemotaxis, differentiation, and activation of peripheral immune cells |
| SERPINA1 | Alpha-1 antitrypsin | AAT, PI | Inhibitor of serine proteases |
| OAS1 | 2′-5′-oligoadenylate synthase 1 | O1AS | Regulation of interferon-gamma-mediated pathway |
| MYC | Myc proto-oncogene protein | BHLHE39 | Activating the transcription of growth-related genes |
| CXCR4 | C-X-C chemokine receptor type 4 | Enhancing MAPK1/MAPK3 activation and involving in the AKT signaling cascade | |
| CXCL2 | C-X-C motif chemokine 2 | GRO2, GROB, MIP2A, SCYB2 | Chemokine-mediated signaling pathway |
| CHRDL1 | Chordin-like protein 1 | NRLN1 | Cell differentiation and negative regulation of BMP signaling pathway |
| GNAI1 | Guanine nucleotide-binding protein G(i) subunit alpha-1 | GTPase activity and regulation of cell cycle and cell division |
Description of the 10 hub genes.
Survival Analysis
Association between 10 hub genes and the overall survival of osteosarcoma patient were analyzed using R package. Based on each hub gene’s best-separation cutoff value, osteosarcoma patients’ samples within the TARGET-osteosarcoma dataset were divided into two groups to get the Kaplan–Meier (K-M) survival curves. The results demonstrated that gene changes of CXCL10 (P = 0.044), GNAI1 (P = 0.048), MYC (P = 0.011), and OAS1 (P = 0.0091) were significantly correlated with the overall survival of osteosarcoma patients (Figure 6).
FIGURE 6
Discussion
More and more research based on public database such as GEO and TARGET database has been done to determine biomarkers in osteosarcoma. For example, Wang et al. () used the microarray data of 42 different age groups (<20- and >20-year-old) osteosarcoma samples from GSE39058 to find 2113 DEGs, including 1476 upregulated and 637 downregulated genes. Similarly, they also identified 15 differentially expressed miRNAs (DEMs) in GSE39040, and functional enrichment analysis showed that upregulated DEMs were mainly enriched in cell growth and response to growth factor, and downregulated DEMs were involved in cytokine receptor activity. Moreover, using GEO database, Dai et al. () screened candidate genes for predicting the response to chemoresistance in osteosarcoma by miRNA–mRNA network. However, the differentially selected genes in these studies are all based on a single dataset, and the small size of samples will cause the instability of results. We integrated three datasets using RRA method, which is standard and robust, compared with other studies on osteosarcoma.
In our study, a total of 524 robust DEGs were determined by RRA method, including 272 upregulated and 252 downregulated genes. The results of GO and KEGG pathway enrichment analyses indicated that the robust DEGs were significantly correlated with ECM organization, cell adhesion molecule binding, cell–cell adherens junction, collagen-containing ECM, and TGF-beta signaling pathway, which were associated with tumorigenesis and metastasis. Through immune cell infiltration analysis, we compared the infiltration of immune cells in primary and metastatic osteosarcoma specimens. We also constructed the PPI network by STRING database and module analysis; finally, we screened 10 hub genes by cytoHubba including POMC, EGR1, CXCL10, SERPINA1, OAS1, MYC, CXCR4, CXCL2, CHRDL1, and GNAI1. Survival analysis of hub genes based on the TARGET database was also performed in our study.
Based on the results of enrichment analyses, the GO terms and KEGG pathways were explored in osteosarcoma. A substantial body of studies indicated that epithelial-to-mesenchymal transition (EMT) is a process needed for metastasis, during which the loss of cell–cell junction such as adherens junctions, allows tumor cells dissociate from the primary site and acquire the motility to invade stroma (–). EMT is also reported to confer resistance to anoikis, which is necessary to survive in the circulation (). In addition to the involvement of EMT process, the tumor-related ECM is also a key factor in tumor progression. In fact, ECM re-organization such as collagen deposition mediated by collagen-binding integrins may be a general cue for prognosis of tumors (). Consistent with the above conclusions, the results of GO enrichment analysis, such as ECM organization, cell adhesion molecule binding, cell–cell adherens junction, and collagen-containing ECM, indicate their involvement in the progression of osteosarcoma. Additionally, enrichment of robust DEGs in some KEGG pathways, including TGF-beta pathway and wnt pathway, also demonstrates their relationship with osteosarcoma development. The overexpression of TGF-βs is related with the presence of lung metastasis () and is associated with high-grade osteosarcoma (). Inhibition of wnt pathway can reduce osteosarcoma invasiveness by reversing the EMT (). On the basis of the above results, we showed that the robust DEGs were highly associated with pathogenesis and progression of osteosarcoma. Furthermore, on the basis of the analysis of modules, we found that three key modules are mainly related to type I interferon signaling pathway, chemokine-mediated signaling pathway, and cell–cell adhesion functions. In the type I interferon pathway, interferon-α has been widely studied. Interferon-α reportedly enhanced the apoptosis of osteosarcoma cells mediated by etoposide () and doxorubicin (). Whether the weakening of interferon signaling pathway plays an important role in the development and metastasis of osteosarcoma deserves further study. It has been reported that chemokine-mediated pathways induce the metastasis of primary tumors to distant target organs. In osteosarcoma, the interaction between chemokine CXCL12 and its receptor CXCR4 drives the metastasis of osteosarcoma cells to the lung (). CXCR3 and its ligands have the same role in promoting lung metastasis of osteosarcoma (). In addition, the CXCR7 receptor promotes osteosarcoma lung metastasis () and has been recognized as a second receptor with high affinity to CXCL12 (). The interaction of CXCR7/CXCL12 in the progression of osteosarcoma needs to be further investigated.
In the past few decades, accumulating evidence has indicated that cancer initiation and progression are related with not only cancer itself but also tumor microenvironment (TME) (, ). TME is a complex including ECM, exosomes, and stromal cells (). Among stromal cells, tumor-associated macrophages (TAMs), namely, the M2 type macrophages, have been reported to promote angiogenesis, matrix remodeling () and are closely associated with osteosarcoma progression and prognosis (). In our study, we found that macrophages are the most principal infiltrating immune cells in osteosarcoma including undifferentiated macrophage M0 and macrophage M2; thus, the role of macrophages, especially M2 type macrophages, in the microenvironment of osteosarcoma needs to be further clarified.
Based on the PPI network construction, 10 hub genes were identified. Among these hub genes, six key genes were screened to explore their roles. Serpin peptidase inhibitor clade A member 1 (SERPINA1), a protease inhibitor, was reported to be a predictor in breast cancer () and colorectal cancer (). However, the diagnostic and prognostic roles of SERPINA1 in osteosarcoma were still obscure. Chemokine CXCL2 was reported to be related with prognosis of bladder cancer (), but its role in the progression of osteosarcoma was still unclear. Early growth response protein 1 (EGR1), a zinc-finger transcription factor, was reported to be involved in cell proliferation and migration (, ). And an increasing number of studies have shown that EGR1 is highly associated with cancer development and progression. Liu et al. reported that EGR1 was essential for HNF1A-AS1-mediated cell growth and invasion of gastric cells (). EGR1 was also reported to promote prostate cancer bone and brain metastasis, as demonstrated by the reduction of blood vessel density in brain and bone caused by decreased EGR1 expression (). In our study, the expression of EGR1 was higher in metastatic osteosarcoma tissues than primary tissues, but its role in osteosarcoma remained unclear. C-X-C motif ligand 10 (CXCL10) is a member of the CXC subfamily of chemokines and acts through CXC receptor 3 (CXCR3) (, ). The prognostic role of serum CXCL10 was proved by Yu et al. in colorectal cancer, and the authors also indicated that the high levels of serum CXCL10 were highly related with liver metastasis (). Similarly, a study was also reported that high circulating levels of CXCL10 are a biomarker for worse survival in osteosarcoma (). According to our result, high expression of CXCL10 in osteosarcoma tissues predicted a better survival. The differences in the above conclusions may be due to the different sources of CXCL10; the association between the CXCL10 expression in tissues and overall survival of osteosarcoma patients need to be further studied. Genetic mutations of tumor suppressor such as TP53 and RB1 are highly associated with osteosarcoma development (, ). In addition, the mis-regulated expression of oncogene MYC is often found in osteosarcoma patients (). A study indicated that the overexpression of c-myc promoted osteosarcoma cells invasion through MEK-ERK pathway (). Consistent with above result, the expression of MYC was higher in metastatic osteosarcoma tissues compared with the primary tissues, and the low expression of MYC predicted a better overall survival in our study. In a word, the MYC expression may act as a biomarker in osteosarcoma metastasis. Chemokine receptor 4 (CXCR4), a seven-transmembrane G protein, has been implicated to mediate the metastasis of several tumors and has become a potential target for tumor therapy (). Several studies highlighted that the overexpression of CXCR4 potentiated osteosarcoma growth and lung metastasis (–). The CXCR4 antibody () and antagonist AMD3100 () were reported to suppress osteosarcoma cell invasion and lung metastasis, and it revealed that CXCR4 may act as a therapeutic agent to inhibit osteosarcoma progression.
Conclusion
In conclusion, by using integrated bioinformatics analysis such as RRA method, we identified the significant robust DEGs and gene modules in osteosarcoma. The enrichment analyses of DEGs showed that they were closely associated with osteosarcoma development and progression. We not only identified immune cell infiltration in osteosarcoma tissues, but we also screened 10 hub genes. After the above discussion, we found that genes EGR1, CXCL10, MYC, and CXCR4 may be considered as novel biomarkers of osteosarcoma, and more studies need to be done to illuminate their contribution in the diagnosis and prognosis of osteosarcoma.
Statements
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
TY and JN conceived the study and thoroughly revised the manuscript. JN, WG, WW, ZZ, TR, and YH acquired and analyzed the data. JN, HZ, YY, and XL wrote the manuscript. All authors approved the final version of the manuscript.
Funding
This work was supported by a grant from the National Natural Science Foundation of China (#81572635 and #81972513).
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.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2020.01628/full#supplementary-material
Footnotes
References
1.
XieBBLiYRZhaoRJXuYZWuYHWangJet alIdentification of key genes and miRNAs in osteosarcoma patients with chemoresistance by bioinformatics analysis.Biomed Res Int. (2018) 10:4761064. 10.1155/2018/4761064
2.
LuetkeAMeyersPALewisIJuergensH.Osteosarcoma. Treatment–Where do we stand? A state of the art review.Cancer Treat Rev. (2014) 40:523–32. 10.1016/j.ctrv.2013.11.006
3.
BielackSSKempf-BielackBDellingGExnerGUFlegeSHelmkeKet alPrognostic factors in high-grade osteosarcoma of the extremities or trunk: an analysis of 1,702 patients treated on neoadjuvant cooperative osteosarcoma study group protocols.J Clin Oncol. (2002) 20:776–90. 10.1200/jco.20.3.776
4.
LindseyBAMarkelJEKleinermanES.Osteosarcoma overview.Rheumatol Therapy. (2017) 4:25–43. 10.1007/s40744-016-0050-2
5.
DamronTAWardWGStewartA.Osteosarcoma, chondrosarcoma, and Ewing’s sarcoma.Clin Orthopaed Relat Res. (2007) 459:40–7. 10.1097/BLO.0b013e318059b8c9
6.
MillerBJCramPLynchCFBuckwalterJA.Risk factors for metastatic disease at presentation with osteosarcoma: an analysis of the SEER database.J Bone Joint Surgery Am. (2013) 95:e89–89. 10.2106/jbjs.l.01189
7.
HeXJGaoZZXuHWZhangZWFuP.A meta-analysis of randomized control trials of surgical methods with osteosarcoma outcomes.J Orthopaed Surgery Res. (2017) 12:6. 10.1186/s13018-016-0500-0
8.
MeyersPASchwartzCLKrailoMDHealeyJHBernsteinMLBetcherDet alOsteosarcoma: the addition of muramyl tripeptide to chemotherapy improves overall survival–a report from the children’s oncology group.J Clin Oncol. (2008) 26:633–8. 10.1200/jco.2008.14.0095
9.
ZhangKFGaoJWNiY.Screening of candidate key genes associated with human osteosarcoma using bioinformatics analysis.Oncol Lett. (2017) 14:2887–93. 10.3892/ol.2017.6519
10.
DiaoCYXiYXiaoT.Identification and analysis of key genes in osteosarcoma using bioinformatics.Oncol Lett. (2018) 15:2789–94. 10.3892/ol.2017.7649
11.
YangXZhuSMLiLZhangLXianSWangYQet alIdentification of differentially expressed genes and signaling pathways in ovarian cancer by integrated bioinformatics analysis.Oncotar Therapy. (2018) 11:1457–74. 10.2147/ott.s152238
12.
NiMWLiuXKWuJRZhangDTianJHWangTet alIdentification of candidate biomarkers correlated with the pathogenesis and prognosis of non-small cell lung cancer via integrated bioinformatics analysis.Front Genet. (2018) 9:14. 10.3389/fgene.2018.00469
13.
HouQBingZTHuCLiMYYangKHMoZet alRankProd combined with genetic algorithm optimized artificial neural network establishes a diagnostic and prognostic prediction model that revealed C1QTNF3 as a biomarker for prostate cancer.Ebiomedicine. (2018) 32:234–44. 10.1016/j.ebiom.2018.05.010
14.
YangYLuQShaoXJMoBHNieXQLiuWet alDevelopment of a three-gene prognostic signature for hepatitis B virus associated hepatocellular carcinoma based on integrated transcriptomic analysis.J Cancer. (2018) 9:1989–2002. 10.7150/jca.23762
15.
RitchieMEPhipsonBWuDHuYFLawCWShiWet allimma powers differential expression analyses for RNA-sequencing and microarray studies.Nucleic Acids Res. (2015) 43:13. 10.1093/nar/gkv007
16.
KoldeRLaurSAdlerPViloJ.Robust rank aggregation for gene list integration and meta-analysis.Bioinformatics. (2012) 28:573–80. 10.1093/bioinformatics/btr709
17.
YuGCWangLGHanYYHeQY.clusterProfiler: an R package for comparing biological themes among gene clusters.Omics J Integrat Biol. (2012) 16:284–7. 10.1089/omi.2011.0118
18.
ChenBKhodadoustMSLiuCLNewmanAMAlizadehAA.Profiling tumor infiltrating immune cells with CIBERSORT.Methods Mol Biol. (2018) 1711:243–59. 10.1007/978-1-4939-7493-1_12
19.
ChinCHChenSHWuHHHoCWKoMTLinCY.cytoHubba: identifying hub objects and sub-networks from complex interactome.BMC Systems Biol. (2014) 8:7. 10.1186/1752-0509-8-s4-s11
20.
WangJ-SDuanM-YZhongY-SLiX-DDuS-XXiePet alInvestigating age-induced differentially expressed genes and potential molecular mechanisms in osteosarcoma based on integrated bioinformatics analysis.Mol Med Rep. (2019) 19:2729–39. 10.3892/mmr.2019.9912
21.
DaiPHeYLuoGDengJJiangNFangTet alScreening candidate microRNA-mRNA network for predicting the response to chemoresistance in osteosarcoma by bioinformatics analysis.J Cell Biochem. (2019) 120:16798–810. 10.1002/jcb.28938
22.
MartinTA.The role of tight junctions in cancer metastasis.Semin Cell Dev Biol. (2014) 36:224–31. 10.1016/j.semcdb.2014.09.008
23.
ThieryJP.Epithelial-mesenchymal transitions in tumour progression.Nat Rev Cancer. (2002) 2:442–54. 10.1038/nrc822
24.
HanahanDWeinbergRA.Hallmarks of cancer: the next generation.Cell. (2011) 144:646–74. 10.1016/j.cell.2011.02.013
25.
YilmazMChristoforiG.EMT, the cytoskeleton, and cancer cell invasion.Cancer Metastas Rev. (2009) 28:15–33. 10.1007/s10555-008-9169-0
26.
YuYElbleRC.Homeostatic signaling by cell-cell junctions and its dysregulation during cancer progression.J Clin Med. (2016) 5:26. 10.3390/jcm5020026
27.
HamidiHIvaskaJ.Every step of the way: integrins in cancer progression and metastasis.Nat Rev Cancer. (2018) 18:532–47. 10.1038/s41568-018-0038-z
28.
YangRSWuCTLinKHHongRLLiuTKLinKS.Relation between histological intensity of transforming growth factor-beta isoforms in human osteosarcoma and the rate of lung metastasis.Tohoku J Exp Med. (1998) 184:133–42. 10.1620/tjem.184.133
29.
FranchiAArganiniLBaroniGCalzolariACapannaRCampanacciDet alExpression of transforming growth factor beta isoforms in osteosarcoma variants: association of TGF beta 1 with high-grade osteosarcomas.J Pathol. (1998) 185:284–9. 10.1002/(sici)1096-9896(199807)185:33.0.co;2-z
30.
GuoYZiXKoontzZKimAXieJGorlickRet alBlocking Wnt/LRP5 signaling by a soluble receptor modulates the epithelial to mesenchymal transition and suppresses met and metalloproteinases in osteosarcoma Saos-2 cells.J Orthopaed Res. (2007) 25:964–71. 10.1002/jor.20356
31.
YuanXWZhuXFLiangSGFanQChenZXHuangXFet alInterferon alpha enhances etoposide-induced apoptosis in human osteosarcoma U20S cells by a p53-dependent pathway.Life Sci. (2008) 82:393–401. 10.1016/j.lfs.2007.11.025
32.
YuanXWZhuXFHuangXFShengPYHeASYangZBet alInterferon-alpha enhances sensitivity of human osteosarcoma U2OS cells to doxorubicin by p53-dependent apoptosis.Acta Pharmacol Sinica. (2007) 28:1835–41. 10.1111/j.1745-7254.2007.00662.x
33.
PerissinottoECavalloniGLeoneFFonsatoVMitolaSGrignaniGet alInvolvement of chemokine receptor 4/stromal cell-derived factor 1 system during osteosarcoma tumor progression.Clin Cancer Res. (2005) 11:490–7.
34.
PradelliEKarimdjee-SoilihiBMichielsJFRicciJEMilletMAVandenbosFet alAntagonism of chemokine receptor CXCR3 inhibits osteosarcoma metastasis to lungs.Int J Cancer. (2009) 125:2586–94. 10.1002/ijc.24665
35.
Goguet-SurmenianERichard-FiardoPGuillemotEBenchetritMGomez-BrouchetABuzzoPet alCXCR7-mediated progression of osteosarcoma in the lungs.Br J Cancer. (2013) 109:1579–85. 10.1038/bjc.2013.482
36.
BurnsJMSummersBCWangYMelikianABerahovichRMiaoZHet alA novel chemokine receptor for SDF-1 and I-TAC involved in cell survival, cell adhesion, and tumor development.J Exp Med. (2006) 203:2201–13. 10.1084/jem.20052144
37.
NaitoYYoshiokaYYamamotoYOchiyaT.How cancer cells dictate their microenvironment: present roles of extracellular vesicles.Cell Mol Life Sci. (2017) 74:697–713. 10.1007/s00018-016-2346-3
38.
WhitesideTL.Exosome and mesenchymal stem cell cross-talk in the tumor microenvironment.Semin Immunol. (2018) 35:69–79. 10.1016/j.smim.2017.12.003
39.
LuoZWangQLauWBLauBXuLZhaoLet alTumor microenvironment: the culprit for ovarian cancer metastasis?Cancer Lett. (2016) 377:174–82. 10.1016/j.canlet.2016.04.038
40.
ArasSZaidiMR.TAMeless traitors: macrophages in cancer progression and metastasis.Br J Cancer. (2017) 117:1583–91. 10.1038/bjc.2017.356
41.
SegalinyAIMohamadiADizierBLokajczykABrionRLanelRet alInterleukin-34 promotes tumor progression and metastatic process in osteosarcoma through induction of angiogenesis and macrophage recruitment.Int J Cancer. (2015) 137:73–85. 10.1002/ijc.29376
42.
ChanHJLiHLiuZYuanY-CMortimerJChenS.SERPINA1 is a direct estrogen receptor target gene and a predictor of survival in breast cancer patients.Oncotarget. (2015) 6:25815–27. 10.18632/oncotarget.4441
43.
KwonCHParkHJChoiJHLeeJRKimHKJoHJet alSnail and serpinA1 promote tumor progression and predict prognosis in colorectal cancer.Oncotarget. (2015) 6:20312–26. 10.18632/oncotarget.3964
44.
ZhangHYeYLLiMXYeSBHuangWRCaiTTet alCXCL2/MIF-CXCR2 signaling promotes the recruitment of myeloid-derived suppressor cells and is correlated with prognosis in bladder cancer.Oncogene. (2017) 36:2095–104. 10.1038/onc.2016.367
45.
WangDGuanM-PZhengZ-JLiW-QLyvF-PPangR-Yet alTranscription factor Egr1 is involved in high glucose-induced proliferation and fibrosis in rat glomerular mesangial cells.Cell Physiol Biochem. (2015) 36:2093–107. 10.1159/000430177
46.
TarcicGAvrahamRPinesGAmitIShayTLuYet alEGR1 and the ERK-ERF axis drive mammary cell migration in response to EGF.Faseb J. (2012) 26:1582–92. 10.1096/fj.11-194654
47.
LiuH-TLiuSLiuLMaR-RGaoP.EGR1-Mediated transcription of lncRNA-HNF1A-AS1 promotes cell-cycle progression in gastric cancer.Cancer Res. (2018) 78:5877–90. 10.1158/0008-5472.can-18-1011
48.
LiLAmeriAHWangSJanssonKHCaseyOMYangQet alEGR1 regulates angiogenic and osteoclastogenic factors in prostate cancer and promotes metastasis.Oncogene. (2019) 38:6241–55. 10.1038/s41388-019-0873-8
49.
SinghUPSinghSIqbalNWeaverCTMcGheeJRLillardJW.IFN-gamma-inducible chemokines enhance adaptive immunity and colitis.J Interfer Cytokine Res. (2003) 23:591–600. 10.1089/107999003322485099
50.
NevilleLFMathiakGBagasraO.The immunobiology of interferon-gamma inducible protein 10kD (IP-10): a novel, pleiotropic member of the C-X-C chemokine superfamily.Cytokine Growth Fact Rev. (1997) 8:207–19. 10.1016/s1359-6101(97)00015-4
51.
ToiyamaYFujikawaHKawamuraMMatsushitaKSaigusaSTanakaKet alEvaluation of CXCL10 as a novel serum marker for predicting liver metastasis and prognosis in colorectal cancer.Int J Oncol. (2012) 40:560–6. 10.3892/ijo.2011.1247
52.
FloresRJKellyAJLiYNakkaMBarkauskasDAKrailoMet alA Novel prognostic model for osteosarcoma using circulating CXCL10 and FLT3LG.Cancer. (2017) 123:144–54. 10.1002/cncr.30272
53.
Kempf-BielackBBielackSSJurgensHBranscheidDBerdelWEExnerGUet alOsteosarcoma relapse after combined modality therapy: an analysis of unselected patients in the cooperative osteosarcoma study group (COSS).J Clin Oncol. (2005) 23:559–68. 10.1200/jco.2005.04.063
54.
KwiatkowskiNZhangTRahlPBAbrahamBJReddyJFicarroSBet alTargeting transcription regulation in cancer with a covalent CDK7 inhibitor.Nature. (2014) 511:616–20. 10.1038/nature13393
55.
GamberiGBenassiMSBohlingTRagazziniPMolendiniLSollazzoMRet alC-myc and c-fos in human osteosarcoma: Prognostic value of mRNA and protein expression.Oncology. (1998) 55:556–63. 10.1159/000011912
56.
HanGWangYBiW.c-Myc overexpression promotes osteosarcoma cell invasion via activation of MEK-ERK pathway.Oncol Res. (2012) 20:149–56. 10.3727/096504012x13522227232237
57.
WongDKorzW.Translating an antagonist of chemokine receptor CXCR4: from bench to bedside.Clin Cancer Res. (2008) 14:7975–80. 10.1158/1078-0432.ccr-07-4846
58.
ZhuYTangLZhaoSSunBChengLTangYet alCXCR4-mediated osteosarcoma growth and pulmonary metastasis is suppressed by MicroRNA-613.Cancer Sci. (2018) 109:2412–22. 10.1111/cas.13653
59.
RenZLiangSYangJHanXShanLWangBet alCoexpression of CXCR4 and MMP9 predicts lung metastasis and poor prognosis in resected osteosarcoma.Tumor Biol. (2016) 37:5089–96. 10.1007/s13277-015-4352-8
60.
XiYQiZMaJChenY.PTEN loss activates a functional AKT/CXCR4 signaling axis to potentiate tumor growth and lung metastasis in human osteosarcoma cells.Clin Exp Metastas. (2020) 37):173–85. 10.1007/s10585-019-09998-7
61.
BrenneckePArltMJECampanileCHusmannKGvozdenovicAApuzzoTet alCXCR4 antibody treatment suppresses metastatic spread to the lung of intratibial human osteosarcoma xenografts in mice.Clin Exp Metastas. (2014) 31:339–49. 10.1007/s10585-013-9632-3
62.
FontanellaRPelagalliANardelliAD’AlterioCIeranoCCerchiaLet alA novel antagonist of CXCR4 prevents bone marrow-derived mesenchymal stem cell-mediated osteosarcoma and hepatocellular carcinoma cell migration and invasion.Cancer Lett. (2016) 370:100–7. 10.1016/j.canlet.2015.10.018
Summary
Keywords
bioinformatics analysis, osteosarcoma, robust rank aggregation, immune cell infiltration, survival analysis
Citation
Niu J, Yan T, Guo W, Wang W, Zhao Z, Ren T, Huang Y, Zhang H, Yu Y and Liang X (2020) Identification of Potential Therapeutic Targets and Immune Cell Infiltration Characteristics in Osteosarcoma Using Bioinformatics Strategy. Front. Oncol. 10:1628. doi: 10.3389/fonc.2020.01628
Received
14 June 2020
Accepted
27 July 2020
Published
21 August 2020
Volume
10 - 2020
Edited by
Zhenyu Jia, University of California, Riverside, United States
Reviewed by
Md. Rakibul Islam, Daffodil International University, Bangladesh; Bikash Kumar Paul, Mawlana Bhashani Science and Technology University, Bangladesh
Updates
Copyright
© 2020 Niu, Yan, Guo, Wang, Zhao, Ren, Huang, Zhang, Yu and Liang.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Taiqiang Yan, yantqzh@163.com
This article was submitted to Molecular and Cellular Oncology, a section of the journal Frontiers in Oncology
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.