- Identification and PCA analysis of the DEGs.
After standardization of the microarray results, DEGs (5,744 in GSE68591, 1,859 in GSE45544 and 2,066 in GSE68776) were identified. The overlap among the 3 datasets contained 142 genes as shown in Fig. 1A. While the PCA dimensionality reduction showed a good distinction between the expression of genes among the tumor sarcoma tissue, normal tissue, and skeletal muscle tissue (Fig. 1B).
- GO functional analysis and pathway enrichment of DEGs
The biological classification, functional and pathway enrichment analyses of DEGs are shown in Fig. 2 (details are shown in Table 1). In biological processes (BP) ontology, DEGs were significantly enriched in positive regulation of cell motility, leukocyte migration, and neutrophil mediated immunity (Fig. 2A). In molecular function (MF) ontology, DEGs were mainly enriched in growth factor binding, extracellular matrix structural constituent, and transmembrane receptor protein kinase activity (Fig. 2B). In cell component (CC) ontology, DEGs were mainly enriched in the extracellular matrix, secretory granule membrane, and specific granule (Fig. 2C). Besides that, KEGG pathway analysis revealed that the terms were related to the Rap1 signaling pathway, Ras signaling pathway, and Glycine, serine, and threonine metabolism (Fig. 2D).
- Construction of the PPI network to select hub genes
The PPI network of the DEGs was constructed by Cytoscape as shown in Fig. 1C, while the most significant module with 10 nodes and 17 edges as shown in Fig. 1D. These 10 genes were identified as hub genes and the detailed information of hub genes is provided in Table 2. The hub genes may be an important regulatory network in the development of ES.
- Regulatory network of the hub genes
The DEG-miRNA (Fig. 3A), TF-DEG (Fig. 3B), hub gene-miRNA (Fig. 3C), and TF-hub gene (Fig. 4D) interaction regulatory networks were constructed and visualized by Cytoscape. The top 5 miRNAs that significantly regulated the DEGs were hsa-mir-124-3p, hsa-mir-335-5p, hsa-mir-16-5p, hsa-mir-26b-5p, and hsa-mir-1-3p While the top 5 TFs significantly interacted with the DEGs were TFDP1, ZNF580, ZBTB11, PHF8, and KDM5B. Furthermore, the top 5 miRNAs that significantly regulated the hub genes were hsa-mir-218-5p, hsa-mir-122-5p, hsa-mir-335-5p, hsa-mir-128-3p, and hsa-mir-124-3p. While the top 10 TFs that significantly interacted with the hub genes were ZNF580, MXD3, SAP30, KDM5B, and PHF8. The topological properties of the miRNAs and TFs are displayed in Table 3 and Table 4. These DEG-miRNA, TF-DEG, hub gene-miRNA and TF-hub gene networks provide evidence for further experimental verification.
- Immune infiltration analysis of the hub genes
The expression levels of CTSB, CTSD, HMOX1, BCL2L1, CD68, TEK, IL1R1, and VEGFC were negatively associated with tumor purity while FGF7 and PTGS2 were positively associated with tumor purity (Fig. 4). This result indicated that most hub genes are probably expressed in the microenvironment, while FGF7 and PTGS2 are probably expressed in the tumor cells.
- DNA methylation analysis and Isomer expression analysis of the hub genes
Furthermore, the relationship between DNA methylation level in the promoter region of the hub genes and the occurrence of SARC was explored. It was found that only the DNA methylation levels of CTSB in SARC were significantly lower than those in corresponding normal tissues, while there was no significant difference in other genes (Fig. 5). Meanwhile, the expression levels of isomers of the hub genes were also analyzed. For BCL2L1, BCL2L1-004 are highly expressed in SARC. For CD68, CD68-001 and CD68-002 are highly expressed in SARC. For CTSB, CTSB-006 and CTSB-001 are highly expressed. For CTSD, CTSD-001 and CTSD-006 are highly expressed. For FGF7, FGF7-001 are highly expressed. For HMOX1, HMOX1-001 are highly expressed. For PTGS2, PTGS2-001 are highly expressed. For TEK, TEK-002 is highly expressed. For IL1R1, IL1R1-005 and IL1R1-001 are highly expressed. For VEGFC, VEGFC-001 is highly expressed (Fig. 6).
- Survival analysis of the hub genes
The overall survival plots of hub genes could be obtained as shown in Fig. 7. Unfortunately, most of them, namely CTSB, CTSD, HMOX1, BCL2L1, CD68, TEK, IL1R1, PTGS2, and VEGFC, had no statistically significant difference in the survival curve. While only FGF7 had a p-value < 0.05. High expression of FGF7 had a better overall survival (Fig. 7).
Together with all the results above, FGF7 is one of the hub genes which could probably be expressed in the tumor cells and its expression had a relationship with the prognosis of ES patients. Therefore, FGF7 was considered as a candidate gene for further analyses.
- Prediction of functions of the candidate genes
A functional network of FGF7 and interactive association genes was constructed by GeneMANIA. The function of this gene set was mainly enriched to fibroblast growth factor receptor signaling pathway, mesenchymal cell proliferation, and renal system development (Fig. 8A). Besides that, the expression profiles of the FGF7 in human tissue are displayed. As shown, FGF71 mRNA in thyroid, lung, and colon normal tissues displayed higher levels than in the matched cancer tissues (Fig. 8B). Subsequently, the gene expression correlation heatmaps indicated that ATP8B4, IL33, A2M, NPR1, and PPAP2A were the most positively correlated significant genes of FGF7 (Fig. 8C), and WNT7B, PRR7, NOTUM, ETV4, and RPSAP52 were the most negatively correlated significant genes of FGF7 (Fig. 8D). All these findings predict the function of FGF7 and can be helpful for further research.