3.1 Four differential expressed CRSs between GC and adjacent tissue samples
In order to identify whether cuproptosis has a high correlation with cancers progression, we obtained differential expression information of 10 CRGs including CDKN2A, DLAT, DLD, FDX1, GLS, LIAS, LIPT1, MTF1, PDHA1, and PDHB from GSCA database for differential expression analyses. The result shows that CDKN2A and LIAS have a significant difference in at least seven cancers (Fig. 1A). After, the expression level of 10 CRGs between the normal and tumor samples was acquired from the TCGA database. The box plots demonstrated that the expression level of CDKN2A and GLS were significantly up-regulated and LIAS and PDHB were significantly down-regulated in GC samples (Fig. 1B). Of note, expression level of DLD significantly changed in Bladder urothelial carcinoma, Breast invasive carcinoma, and Thyroid carcinoma (Fig. 1C-E). Besides, the box plots visually reveal the normalized expression level of each CRG between the normal stomach and GC samples that with pairing correlation (Fig. 1F). The heatmap shows that there is no significant correlation among CRGs (Fig. 1G). Moreover, tumor-specific CNVs and SNVs are helpful for exploring the molecular mechanisms of GC progression. We analyzed the mutation frequencies of CRGs of GC samples, and the result shows that the LIPT1 and CDKN2A have the highest mutation frequency of 6% and 4% respectively (Fig. 1H). After analysis of CNV was performed and we found that the CDKN2A, LIAS, and PHDB frequently occur the copy number deletion (DEL), and DLD frequently occurs the cope number amplification (AMP) (Fig. 1I).
3.2 Survival analyses between high- and low-expression groups
A series of survival analyses between high- and low-expression groups were performed to identify the correlation between the expression level of CRGs and GC prognosis. The OS curves show that the expression level of MTF1, DLAT, and LIAS are positively correlated, while the expression level of GLS is negatively correlated with the prognosis of GC (Additional File 2). After, the disease-special survival analysis results demonstrate that the expression level of MTF1, DLAT, PDHA1, DLD, and LIAS are positively correlated, while the expression level of CDKN2A is negatively correlated with the prognosis of GC (Additional File 3). Moreover, disease-free survival analysis results illustrate that the expression level of three CRGs including MTF1, PDHA1, and DLD are positively correlated with the prognosis of GC (Additional File 4). Finally, progressive-free survival curves show that the expression level of five CRGs including MTF1, PDHA1, DLD, LIAS, and DLAT is positively correlated with the prognosis of GC (Additional File 5). As shown in Table 1, there were 4 differentials expressed CRGs had a significant survival relationship.
Table 1
The correlation between survival and CRGs in GC patients.
gene
|
OS
|
DSS
|
DFS
|
PFS
|
*GLS
|
TRUE
|
FALSE
|
FALSE
|
FALSE
|
*LIAS
|
TRUE
|
TRUE
|
FALSE
|
TRUE
|
*CDKN2A
|
FALSE
|
TRUE
|
FALSE
|
FALSE
|
*DLD
|
FALSE
|
TRUE
|
TRUE
|
TRUE
|
*PDHB
|
FALSE
|
FALSE
|
FALSE
|
FALSE
|
DLAT
|
TRUE
|
TRUE
|
FALSE
|
TRUE
|
MTF1
|
TRUE
|
TRUE
|
TRUE
|
TRUE
|
FDX1
|
FALSE
|
FALSE
|
FALSE
|
FALSE
|
LIPT1
|
FALSE
|
FALSE
|
FALSE
|
FALSE
|
PDHA1
|
FALSE
|
TRUE
|
TRUE
|
TRUE
|
* represents differentially expressed genes in GC; TRUE represents survival significance. OS: overall survival; DSS: disease-special survival; DFS: disease-free survival; PFS: progressive-free survival.
3.3 Three biomarkers of clinical prognosis prediction were obtained
We perform single-gene ROC analysis to evaluate the value of CRGs for GC prognosis, and the ROC curves demonstrate that PDHB, CDKN2A, LIAS, and GLS have great prognosis prediction abilities (Fig. 2A-D). Then, considering the differential expression between normal and GC samples, and the survival correlation with GC, three CRGs including CDKN2A, LIAS, and GLS were chosen as biomarkers of GC prognosis. After, the expression level of biomarkers between high- and low-risk groups in different clinical subgroups including age, gender, vital, state, grade, T stage, N stage, and M stage was compared. We found that the expression level of CDKN2A between G2 and G3 is significantly different. The expression level of GLS in T3 and T4 are significantly higher than in T1, while the expression level of LIAS in T2 and T3 are significantly lower than in T1 (Fig. 2E-G). The total information is provided in Additional File 1. Finally, clinical features were combined with the expression of biomarkers to generate a nomogram of GC prognosis prediction (Fig. 2H). The calibration curve demonstrates that the nomogram has a great prediction accuracy of 1-year survival (Fig. 2I).
3.4 Three methylation sites of CRGs significantly related to survival
Gene methylation significantly affects gene function. In this study, the methylation site of biomarkers was detected through the UCSC Xena website, and visualized in a heatmap (Fig. 3A). The scatter plots show that cg07562918 and cg13601799 were significantly related with CDKN2A, cg07253264 and cg09390371 were significantly associated with LIAS, and cg03962451 and cg19300307 were significantly correlated with GLS (Fig. 3B-G). Total information of the correlation of methylation sites were provided in Additional File 1. Finally, K-M survival analysis between high- and low-expressed groups was utilized, and we found that cg13601799 and cg07562918 of CDKN2A ,and cg07253264 of LIAS have a strongly significant correlation with GC survival (Fig. 3H-J).
3.5 Functional enrichment of CDKN2A, GLS, and LIAS
In order to identify the biological functions of biomarkers, differential expression analysis was implemented to screen DEGs between high- and low-expressed groups based on the expression level of biomarkers, and then DEGs were used for enrichment analysis. Enrichment results of biomarker CDKN2A show that there are 82 functions (such as digestion, antimicrobial humoral response, and thyroid hormone generation) in GO-BP, 35 functions (such as apical plasma membrane, apical part of cell, and chylomicron) in GO-CC, 33 functions (such as endopeptidase activity, cholesterol transfer activity, and sterol transfer activity) in GO-MF (Fig. 4A), and 9 KEGG pathways such as Fat digestion and absorption, Cholesterol metabolism, and Proximal tubule bicarbonate reclamation (Fig. 4B). Enrichment results of biomarker GLS demonstrate that there are 90 functions (such as digestion, digestive system process, and maintenance of gastrointestinal epithelium) in GO-BP, 12 functions (such as digestion, digestive system process, and maintenance of gastrointestinal epithelium) in GO-CC, 42 functions (such as aspartic-type endopeptidase activity, aspartic-type peptidase activity, and carboxylic acid transmembrane transporter activity) in GO-MF (Fig. 4C), and 4 KEGG pathways include Fat digestion and absorption and Protein digestion, absorption, Pancreatic secretion, and Mineral absorption (Fig. 4D). Enrichment results of biomarker LIAS illustrate that there are 171 functions (such as uterus development, skin development, and negative regulation of cardiac muscle tissue development) in GO-BP, 51 functions (such as serine-type endopeptidase inhibitor activity, co-receptor binding, and endopeptidase inhibitor activity) in GO-MF (Fig. 4E), and 3 KEGG pathways include Cell adhesion molecules, Wnt signaling pathway, and Fat digestion and absorption (Fig. 4F).
3.6 Expression of biomarkers is correlated with immune cells
Immune cells are an important part of the tumor microenvironment (TME) and profoundly influence the progression of GC. We calculated the percentage of the immune cells in GC samples, and the results were visualized as a stacked bar chart (Fig. 5A). The bubble charts demonstrate that biomarker GLS has a negative correlation with T cells regulatory (Tregs), NK cells resting, and Neutrophils, biomarker LIAS has a positive correlation with T cells CD4 memory activated, and has a negative correlation with Macrophages M0 and Macrophages M1, and biomarker CDKN2A has a negative correlation with B cells naive and Plasma cells (Fig. 5B). The results of differential analysis demonstrated that the percentage of B cells native in high-expressed group of CDKN2A was significant different with low-expression group (Fig. 5C), And the proportion of B cells naive, T cells regulatory (Tregs), Mast cells resting, Mast cells activated, and Dendritic cells resting in high-expressed group of GLS were significant different with the proportion in low-expressed group(Fig. 5D), and the percentage of B cells naive, T cells regulatory (Tregs), Mast cells resting, Mast cells activated, and Dendritic cells resting in high-expressed group of LIAs is significant different with low-expressed group (Fig. 5E). The correlations between expression level of biomarkers and immune cells were visualized in lollipop charts (Fig. 5F-H). Besides, we detected the correlation between biomarkers and Immune checkpoints, and the result shows that CDKN2A has a negative correlation with IDO1, TIGIT, CD274, PDCD1LG2, ICOS, and HAVCR2, LIAS is negatively related with CD27, and GLS is positively related with ICOS (Fig. 5I).
3.7 MiRNA-mRNA network of biomarkers
Mi-RNAs can affect the function of biomarkers by degrading mRNA. The volcano plot indicates that there are 724 DE-mRNAs between normal and GC samples (Fig. 6A). Then, a miRNA-mRNA network of biomarkers including 22 miRNAs was generated on the basis of correlation of miRNA-mRNA which identified via miRWalk database (Fig. 6B). The heatmap of correlation biomarkers and miRNAs indicates that CDKN2A significantly correlate with hsa-miR-330-5P, hsa-miR-181d-5P, hsa-miR-181b-5P, and hsa-miR-1301-3P, and GLS significantly correlates with hsa-miR-34a-5P, hsa-miR-181d-5P, and hsa-miR-1301-3P, and LIAS significantly correlates with 6 miRNAs including hsa-miR-942-5P, hsa-miR-432-5P, hsa-miR-181d-5P, hsa-miR-181b-5P, hsa-miR-125b-5P, and hsa-miR-125a-5P (Fig. 6C). Besides, there are 8 miRNA are significantly correlated with OS of patients, which includes hsa-miR-185-5P, hsa-miR-320b, hsa-miR-181b-5P, hsa-miR-34b-5P, hsa-miR-370-3P, hsa-miR-490-3P, hsa-miR-942-5P, and hsa-miR-432-5P (Fig. 6D-K).
3.8 TFs-mRNA network of biomarkers
TFs can affect the function of biomarkers by affecting the transcription process of mRNA. Differentially expressed gene analysis was implemented to obtain 3489 DEGs between normal and GC samples, and then those DEGs were intersected with TFs identified through the NetworkAnalyst database (Fig. 7A, B). There are 5 differentially expressed TFs in the intersection and they were used for constructing a TFs-mRNA network (Fig. 7C). The heatmap demonstrates that CDKN2A significantly correlates with RCOR2, GLS significantly correlates with SOX5, RCOR2, EZH2, and ELF3, and LIAS is significantly correlated with SOX5, KLF9, EZH2, and ELF3 (Fig. 7D). The survival analysis of 5 TFs indicated that 3 TFs including EZH2, SOX5, and KLF9 have a significant correlation with OS of GC patients (Fig. 7E-G).
3.9 Co-expression Network of biomarkers
We constructed a Co-expression Network to indicate reciprocity between DEGs and biomarkers (Fig. 8A), and there are 82 genes co-express with GLS, 21 genes co-express with LIAS, and 19 genes co-express with CDKN2A in the network. Then, we performed enrichment analysis on the basis of the network, and we found that there are 20 functions of biomarkers and co-expressed genes, such as cardiac septum development, collagen fibril organization, and sensory organ morphogenesis (Fig. 8B).