Identification of DEGs
Firstly, all the genes from GSE15641, GSE72304, GSE71963, GSE53757, and GSE36895 were shown in the volcano plots (Figure 1A). In this study, 173 KIRC tissues and 138 normal kidney tissues from five GEO data sets were involved. A total of 1760 upregulated and 421 downregulated genes were filtered from GSE15641; 962 upregulated and 396 downregulated genes from GSE72304; 1580 upregulated and 1484 downregulated genes from GSE71963; 1593 upregulated and 1985 downregulated genes from GSE53757; 910 upregulated and 1082 downregulated genes from GSE36895 based on the threshold of P<0.05 and |logFC| > 1, respectively. Subsequently, 99 common upregulated genes (Figure 1B) and 42 common downregulated genes (Figure 1C) were selected among the DEGs with Venn software online. The top 25 upregulated and downregulated genes were shown in Heat map with R package “ComplexHeatmap” (Figure 2).
GO and KEGG pathway analyses
GO analysis of the 141 DEGs (99 upregulated genes and 42 downregulated genes) was performed with the DAVID online analysis tool. In biological process (Figure 3A), the DEGs were significantly enriched in oxidation-reduction process (GO:0055114), response to hypoxia (GO:0001666), gluconeogenesis (GO:0006094), cellular response to hypoxia (GO:0071456), response to drug (GO:0042493), and lung development (GO:0030324). In cellular component terms (Figure 3B), the DEGs were mainly involved in the extracellular exosome (GO:0070062), extracellular space (GO:0005615), extracellular region (GO:0005576), protein complex (GO:0043234), and extracellular matrix (GO:0031012). For molecular functions (Figure 3C), the DEGs were mainly enriched in iron ion binding (GO:0005506), transporter activity (GO:0005215), fatty acid binding (GO:0005504), L-ascorbic acid binding (GO:00314180), and receptor binding (GO:0005102). The KEGG pathway analysis results were shown in Figure 3D, and it demonstrated that these targets were particularly enriched in PPAR signaling pathway, carbon signaling pathway, HIF-1 signaling pathway, valine, leucine ang isoleucine degradation, and gluconeogenesis.
PPI network and hub gene selection
The PPI network was constructed for further investigation of the interaction among the DEGs used Cytoscape software, based on the STRING database. The PPI network consisted of 104 nodes and 220 edges. 5 central nodes modules of the PPI networks were selected with Cytotype MCODE (Figure 4). The identified hub genes were AOX1, ALDH6A1, ABAT, HADH, and PCCA.
Hub gene validation
To further evaluate the role of the five candidate hub genes in KIRC prediction, TCGA_KIRC database was selected for validation. The results showed that the expression of these five hub genes were significantly downregulated in KIRC, which consistent with the mRNA level from GEO dataset (Figure 5A). To further clarify the clinical significance of AOX1, ALDH6A1, ABAT, HADH, and PCCA expression in KIRC, we analyzed their expression in different stages of KIRC. The low expression of hub genes was closely related with advanced tumor histology stage (Figure 5B), pathological stage (Figure 5C), distant metastasis (Figure 6A), and lymph node metastasis (Figure 6B), tumor status (Figure 6C) in KIRC. These results demonstrated that the low expression of AOX1, ALDH6A1, ABAT, HADH, and PCCA might be associated with KIRC progression and have prognostic significance for KIRC patients.
Prognostic significance of hub genes in patients with KIRC
We investigated the R package “survival” for the prognostic significance (overall survival and disease specific survival) of AOX1, ALDH6A1, ABAT, HADH, and PCCA expression in KIRC. These findings indicated that low expression of AOX1, ALDH6A1, ABAT, HADH, and PCCA predicted poor prognostic in in patients with KIRC (Figures 7A and 7B), which suggested that AOX1, ALDH6A1, ABAT, HADH, and PCCA were significant factors for predicting worse prognostic in KIRC.
ROC curve analysis of hub genes
ROC curve analysis was performed on five hub genes with the R package of “pROC”. AUC > 0.90 was selected as the cutoff value, and it was found that three of the five hub genes with AUC >90% included ALDH6A1 (AUC=0.939), HADH (AUC=0.979), and PCCA (AUC=0.931), respectively. The expression of ALDH6A1, HADH and PCCA have high accuracy in distinguishing normal kidney tissue from KIRC tissue, which has important significance for the accurate diagnosis of KIRC (Figure 8).
To further identify the possible mechanism of the hub genes in KIRC, GSEA was conducted to obtain the biological pathway. The TCGA_KIRC data set was divided into high and low‐expression groups based on the median value of AOX1, ALDH6A1, ABAT, HADH, and PCCA, respectively. As shown in Figure 9, the result revealed that the low expression of AOX1 was enriched in the “REACTOME_FCERI_MEDIATED_NF-KB activation”, ALDH6A1 was enriched in the “KEGG_CYTOKINE_CYTOKINE_RECEPTOR_INTERACTION”, ABAT was enriched in the “KEGG_FOCAL_ADHESION”, HADH was enriched in the “REACTOME_SIGNALING_B_CELL_RECEPTOR”, and PCCA was enriched in the “KEGG_JAK-STAT_SIGNALING_PATHWAY”, respectively.
The hub genes expression in different cancer types
We applied the Oncomine database to analyze the expression of AOX1, ALDH6A1, ABAT, HADH, and PCCA. The results showed that these genes expression was significantly different in various cancer types (Figure 10A). And AOX1, ALDH6A1, ABAT, HADH, and PCCA expression were significantly low in KIRC, compared with the corresponding normal tissues (Figure 10B).