Skip to main content

Constructing a novel mitochondrial-related gene signature for evaluating the tumor immune microenvironment and predicting survival in stomach adenocarcinoma

Abstract

Background

The incidence and mortality of gastric cancer ranks fifth and fourth worldwide among all malignancies, respectively. Accumulating evidences have revealed the close relationship between mitochondrial dysfunction and the initiation and progression of stomach cancer. However, rare prognostic models for mitochondrial-related gene risk have been built up in stomach cancer.

Methods

In current study, the expression and prognostic value of mitochondrial-related genes in stomach adenocarcinoma (STAD) patients were systematically analyzed to establish a mitochondrial-related risk model based on available TCGA and GEO databases. The tumor microenvironment (TME), immune cell infiltration, tumor mutation burden, and drug sensitivity of gastric adenocarcinoma patients were also investigated using R language, GraphPad Prism 8 and online databases.

Results

We established a mitochondrial-related risk prognostic model including NOX4, ALDH3A2, FKBP10 and MAOA and validated its predictive power. This risk model indicated that the immune cell infiltration in high-risk group was significantly different from that in the low-risk group. Besides, the risk score was closely related to TME signature genes and immune checkpoint molecules, suggesting that the immunosuppressive tumor microenvironment might lead to poor prognosis in high-risk groups. Moreover, TIDE analysis demonstrated that combined analysis of risk score and immune score, or stromal score, or microsatellite status could more effectively predict the benefit of immunotherapy in STAD patients with different stratifications. Finally, rapamycin, PD-0325901 and dasatinib were found to be more effective for patients in the high-risk group, whereas AZD7762, CEP-701 and methotrexate were predicted to be more effective for patients in the low-risk group.

Conclusions

Our results suggest that the mitochondrial-related risk model could be a reliable prognostic biomarker for personalized treatment of STAD patients.

Introduction

The incidence (accounting for 5.6% of all cancer cases) and mortality (accounting for 7.7% of all cancer deaths) of gastric cancer (GC) ranks fifth and fourth worldwide among all malignancies, which critically threatens human health. It was estimated that over one million new cases of GC were reported worldwide in 2020 [1]. Gastric adenocarcinomas derived from gastric glandular epithelial cells accounts for more than 90% of all GCs [2]. GC is a multifactorial disease, which is contributed by both environmental and genetic factors [3], such as smoking, family history, Epstein–Barr virus (EBV) infection, alcohol consumption, and diet [4]. Most GCs are diagnosed at late stage of progression due to limited premalignant indications and symptoms [2]. Nowadays, the application of the endoscopic examination largely improved the survival rate of GC patients, and a 30% reduction in GC mortality using endoscopic screening [3]. However, stomach cancer is still one of the most lethal malignant tumors, with a 5 year survival rate of around 20% [5]. Therefore, it is necessary to identify more reliable biomarkers for predicting the prognosis and exploring more potential therapeutic targets in GC.

Accumulating evidences indicated that mitochondria plays essential roles in regulation of cell growth, cell death, and cell metabolism during the whole process of tumor progression [6]. Mitochondria are involved in bioenergetics metabolism, such as ATP production, reactive oxygen species (ROS) production, apoptosis, and calcium homeostasis [7]. Moreover, mitochondrial dysfunction may contribute to the chemoresistance [8]. Therefore, mitochondrial-targeting therapies may be applied for the treatment of GC, including ROS production and elimination, mitochondrial fission and fusion, ATP production, and apoptosis [6, 9, 10]. For instance, nanohybrid-induced oxidative stress triggered mitochondria-mediated autophagy, which inhibited cell growth in cancer cell [11].

Considering that mitochondrial dysfunction was a risk factor for the tumorigenesis of GC, identifying effective mitochondrial-related biomarkers for the prognosis of GC patients should be an encouraging direction of research. Several studies have constructed GC prognosis-related models to predict patient survival [12,13,14]. However, rare studies have been applied for the establishment of prognostic models for GC associated with mitochondria.

Tumor microenvironment (TME) was mainly composed of the stromal cells, immune cells and cytokines [15]. The components of TME affected the immune cell evasion or inhibition, and drug resistance in malignancies. For example, immune progenitors in the complex microenvironment of the TME were more likely to differentiate into M2 macrophages and Treg cells, but not to play their tumor-inhibiting functions as fully mature immune cells [16]. In addition, the response to immune checkpoint blockade (ICB) was closely related to the constitution of the TME. ICB revived an effective anti-tumor immune response [17]. It was reported that PD-1/PD-L1 inhibitors immunotherapy has an impact on the therapy of patients with advanced gastrointestinal malignancy. Moreover, a study reported that tumors with a higher tumor mutation burden (TMB), had a better immunotherapy response, especially with PD-1/PD-L1 blockade [18]. Thus, to figure out the correlation between the risk score and TME, we explored the TME signatures, the immune cell infiltration in TME, the expression level of immune checkpoints, and the response to immunotherapy.

In summary, a novel mitochondrial-related risk model was constructed using NOX4, FKBP10, ALDH3A2, and MAOA gene set, which could effectively predict the prognosis and immunotherapy responsiveness for patients with STAD. In addition, we estimated the drug sensitivity of STAD patients to 138 drugs, including chemotherapy drugs, immunotherapy drugs, and targeted drugs, et al., and found that patients in high-risk group was more sensitive to rapamycin, PD-0325901 and dasatinib, whereas patients in low-risk group was more sensitive to AZD7762, CEP-701 and methotrexate. Taken together, our mitochondrial-related risk model could be a reliable prognostic biomarker for personalized treatment of STAD patients.

Methods

Data collection

RNA-seq data and microsatellite status information for 407 STAD samples were downloaded from the TCGA database (https://www.cancer.gov/tcga). The clinical information was extracted from the UCSC Xena (http://xena.ucsc.edu) [19]. 61 samples were excluded due to incomplete clinical information or survival less than 30 days. In total, 346 samples, comprising 317 tumor samples and 29 healthy samples, were analyzed in the present study. Two validation cohorts, GSE66229 and GSE15459, were applied in the present study (https://www.ncbi.nlm.nih.gov/geo/). In GSE15459, 10 samples were excluded due to survival less than 30 days. Altogether, 300 and 182 samples were analyzed in GSE66229 and GSE15459, respectively. The list of mitochondrial-related genes was collected from MitoCarta 3.0 database (https://www.broadinstitute.org/mitocarta/mitocarta30-inventory-mammalian-mitochondrial-proteins-and-pathways) [20] and the Gene set enrichment analyses (GSEA, http://www.gsea-msigdb.org/gsea/index.jsp) [21, 22] (Additional file 10: Table S1).

Identification of differentially expressed genes (DEGs)

The “limma” package of R (version 3.5.1) was applied to produce DEGs between normal and tumor samples, or between high-risk and low-risk groups from the training set. |Log (2) fold change|> 2 and adjusted P < 0.01 were the criteria for defining DEGs. “GdcVolcanoPlot” packages in R were employed to generate volcano map to visualize the DEGs, and a Venn plot was exploited to display the common DEGs in both DEGs groups and mitochondrial-related genes.

Construction and validation of prognostic mitochondrial-related risk score signature

The mitochondrial-related genes were screened by univariate Cox regression, Lasso regression analysis and multivariable Cox regression analysis to construct a novel prognostic gene signature. Each sample’s risk score was calculated using the following formula:

$${\text{Risk score }} = \, \Sigma {\text{expgenei}} * {\beta i}$$

where expgene, i, and βi represent the expression level of gene, the number of signature genes, and the coefficient index, respectively. In all participated cohorts, the samples were divided into low-risk and high-risk groups based on the risk score (median cut-off value). To analyze the survival conditions for the prognosis signature, the optimized cutoff and the Kaplan–Meier (K–M) survival curve were conducted by R package “survival” and “survminer”. The predictive performance was presented by ROC curve, risk plot and concordance index (C-index). Detailed information for prognostic genes was obtained from The Human Protein Atlas (HPA, https://www.proteinatlas.org/) and National Center for Biotechnology Information (NCBI, https://www.ncbi.nlm.nih.gov/).

Construction and valuation of nomogram

Risk score and clinical factors including age, gender, T stage, N stage, M stage, tumor stage, family history, H pylori infection, grade, reflux history, and disease types were analyzed using univariate Cox regression analysis to screen the factors significantly related to survival (P < 0.1). Then, multivariate COX regression analysis was applied to identify the candidate predictors significantly related to survival (P < 0.05). Based on this, nomograms were constructed using these predictors, and scores in nomogram model were assigned for these variables. By adding the scores of the predictors enrolled in nomogram model, the total score of each patient was obtained. Finally, the patient's survival outcome in 1, 3 and 5 years can be calculated using the total score and the probability of survival outcome. ROC curve, calibration curves and decision curve analysis (DCA) were applied to estimate the discrimination and accuracy of the nomogram model.

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analyses

In the present study, R “clusterProfiler”, “org.Hs.eg.db”, “enrichplot” and “ggplot2” package (R version: 3.5.1) were employed to analyze the function of mitochondrial-related DEGs, or the DEGs between high-risk and low-risk groups. Furthermore, adjusted P < 0.05 was used to filter the functional candidates.

Gene set enrichment analyses (GSEA)

Curated sets v7.4 collections were obtained from the Molecular Signatures Database as the target sets with which GSEA was performed by using GSEA 4.2.1 software. The total transcriptome of tumor samples was used for the GSEA, and only gene sets with P < 0.001 and FDR, q < 0.001 were regarded to be statistically significant.

Tumor microenvironment

Stromal scores were calculated using the ESTIMATE algorithm by R (version 3.5.1) package “estimate”. The list of TME-related biomarkers was extracted from the Gene set enrichment analyses (GSEA, http://www.gsea-msigdb.org/gsea/index.jsp) [21, 22] (Additional file 10: Table S2–S5).

Calculation of relative abundance of 22 immune cell subtypes

The abundance of 22 tumor-infiltrating immune cells (TIICs) in STAD samples was calculated using the CIBERSORT algorithm by R package. CIBERSORT is a deconvolution algorithm that can infer 22 kinds of TIICs and harnesses the ability to predict the relative abundance of each immune cell population by calculating the expression of specific marker [23]. The relative abundance of the TIICs between high-risk and low-risk groups was compared using the Wilcox text. The list of immune check points was referenced from a published study [24]. Immune score and tumor purity were also calculated using the ESTIMATE algorithm by R (version 3.5.1) package “estimate”. The list of immune cell signatures was downloaded from TISIDB (http://cis.hku.hk/TISIDB/download.php) [25].

Prediction of therapeutic sensitivity in patients with different risk scores

The capability of risk score in predicting the response to immunotherapy or 138 drugs for chemotherapies/targeted therapies was explored in the present study. The 50% inhibiting concentration (IC50) values of the 138 drugs were calculated using the “pRRophetic” package of R (version 3.5.1) and the value was normally transformed. The detailed information of 138 drugs was acquired from Genomics of Drug Sensitivity in Cancer (GDSC, https://www.cancerrxgene.org/). The potential response to immunotherapy was inferred by the tumor immune dysfunction and exclusion (TIDE, http://tide.dfci.harvard.edu) score.

Mutation analysis

The somatic mutation data were downloaded from cBioportal database (https://www.cbioportal.org/) [26, 27]. The R (version 3.5.1) package “maftools” was then used to draw a waterfall plot to illustrate the mutation landscape in STAD patients with the high- and low-risk group and calculate the TMB score for each sample.

Cell culture and patient sample collection

The normal gastric epithelial cell line GES-1 and human gastric cancer cell lines SGC-7901 and HGC-27 were purchased from the American Type Culture Collection (ATCC, Manassas, VA, USA). Cells were cultured in DMEM medium (Gibco, Thermo Fisher Scientific, Inc., Waltham, MA, USA) supplemented with 10% fetal bovine serum (FBS, Hyclone). Cells were routinely cultured in a humidified atmosphere containg 5% CO2 at 37 ℃. A total of 10 fresh tumor and paired adjacent normal tissues from patients with STAD were collected in the First Affiliated Hospital of Shanxi Medical University (Taiyuan, Shanxi, China). All patients provided written informed consent, and this study was approved by the ethics committee of Shanxi Medical University.

RNA extraction and qRT-PCR assays

RNAs in tissues and cell lines of STAD were extracted with a RNAiso Plus (Takara,Tokyo, Japan) and were reversely transcribed into cDNA usin PrimeScriptTM RT Master Mix (Takara). Quantitative real-time PCR (qRT-PCR) was performed by TB Green ®Premix Ex TaqTM II (Takara). β-actin was used as reference genes. The primers were listed in Additional file 10: Table S6.

Statistical analysis

R (version 3.5.1) and the GraphPad Prism 8 software were applied for statistical analysis. A Student’s t-test was used to analyze the expression and the distribution of risk score, stromal score, immune score, tumor purity and TMB in different groups. Chi-square test is applied to evaluate the difference in immunotherapy response, status of top 5 mutant genes and clinical factors in different groups. The correlation was evaluated using the Spearman method. C-index was used to estimate the predictive power of age and risk score to OS. P < 0.05 was defined as statistically significant.

Results

Identification of DEGs related to mitochondrion and functional enrichment analysis in STAD

The general workflow of our current study was presented in Fig. 1. As shown in Additional file 10: Table S7, 2381 DEGs, including 2145 protein-coding genes, were screened and visualized via volcano maps between normal and tumor groups (Fig. 2A, Additional file 10: Table S8). Next, combined analysis for selected mitochondrial-related genes from the GSEA and 2145 DEGs from our study were performed to filter out 183 candidate mitochondrial-related DEGs in STAD (Fig. 2B, Additional file 10: Table S9).

Fig. 1
figure 1

Workflow diagram. The flowchart graph of this study

Fig. 2
figure 2

Identification of DEGs related to mitochondrion and construction of prognostic risk model using TCGA-STAD cohort. A Volcano plot of 2381 DEGs in STAD tumor and normal groups. B Venn diagram showed that the overlap of 2381 DEGs and 2030 mitochondrial genes led to 183 hub genes being identified. C Univariate Cox regression analysis revealed 19 genes were associated with prognosis of patients with STAD. D LASSO regression of the 19 OS-related genes. Cross-validation in the LASSO regression model to select the tuning parameter. The abscissa shows the log (λ) value, and the ordinate shows partial likelihood deviance. The red dots in the figure show partial likelihood deviations ± standard error for diverse tuning parameters. E Multivariable Cox regression analysis revealed 4 genes were associated with prognosis of patients with STAD. F Gene expressions of the 4 prognosis-related genes in TCGA-STAD. P values were showed as: ***P < 0.001

GO enrichment analysis were then carried out to uncover important roles of mitochondrial-related DEGs in STAD. These DEGs potentially participated in small molecule catabolic process, regulation of mitochondrial organization, et al. Regarding the cellular component, they were mainly related to mitochondrial matrix, mitochondrial outer membrane, et al. In terms of molecular function, these DEGs were involved in tubulin binding, and ubiquitin-like protein ligase binding, et al. (Additional file 1: Fig. S1A, Additional file 10: Table S10). Moreover, KEGG pathway analysis was also applied to demonstrate important pathways being involved by these DEGs, such as lipid and atherosclerosis, Hepatitis B infection, Diabetic cardiomyopathy, and apoptosis, et al. (Additional file 1: Fig. S1B, Additional file 10: Table S11).

Construction and validation of a mitochondrial-related risk signature

Based on above 183 mitochondrial-related DEGs, 19 genes were further selected as potential risk factors for the prognosis of patients with STAD through univariate Cox regression analysis (P < 0.05, Fig. 2C). The gene number was further narrowed down to 9 according to LASSO regression analysis and to 4 by multivariable Cox regression analysis (Fig. 2D, E). Finally, 4 mitochondrial-related DEGs, including NOX4, FKBP10, ALDH3A2 and MAOA, were utilized to establish a prognostic model for patients with STAD (Table 1).

Table 1 The information of 4 prognosis-related genes

As shown in Fig. 2F, higher expressions of NOX4 and FKBP10, and lower expressions of ALDH3A2 and MAOA were observed in tumor samples compared with the normal tissues, respectively. Moreover, the immunohistochemistry results from HPA database showed that FKBP10 was upregulated in gastric cancer tissues, while ALDH3A2 and MAOA were downregulated in gastric cancer tissues, when compared with corresponding non-cancerous tissues (Additional file 1: Fig. S1C). Further K-M analysis demonstrated that patients with higher expression of NOX4 (P = 0.030), FKBP10 (P = 0.040), and MAOA (P = 0.018) had a shorter OS than those with lower expression, respectively (Additional file 1: Fig. S1D). However, the patients with higher level of ALDH3A2 had a better OS than those with lower expression, even though it is a bit beyond statistically significant difference (P = 0.052, Additional file 1: Fig. S1D).

Then, the risk score for each patient with STAD in both training and validation cohorts was computed based on the following formula:

$${\text{Risk score}}\, = \,0.{157}*{\text{MAOA}} - 0.{198}*{\text{ALDH3A2}}\, + \,0.{133}*{\text{FKBP1}}0\, + \,0.{146}*{\text{NOX4}}{.}$$

Patients were divided into high-risk and low-risk subgroups based on the median risk score. K–M curves showed that patients in high-risk group had worse OS (P = 0.0009, Fig. 3A). To assess the accuracy of prognostic risk models in predicting 1-, 3-, and 5-year OS, ROC curves were plotted with AUC values of 0.635, 0.640, and 0.793, respectively (Fig. 3B). The relationship between the risk score and the survival time, survival status, and risk ranking, and a heatmap of the expressions of the 4 genes were shown in Fig. 3C. Taken together, these results demonstrated the robustness of our risk model in predicting the prognosis of patients with STAD.

Fig. 3
figure 3

Assessing the performance of the prognostic risk model in the training and validation cohort. A, D Kaplan–Meier curves of the OS of patients in the high- and low-risk groups in the TCGA-STAD training cohort (A), and GSE66229 cohort (D). B, E ROC curves for predicting 1-, 3-, and 5-year OS in the TCGA-STAD training cohort (B), and GSE66229 cohort (E). C, F Distribution of risk score, survival status (red dots indicate dead, blue dots indicate alive) and the gene expression of 4 model genes in the TCGA-STAD training cohort (C), and GSE66229 cohort (F)

The robustness of the prognostic risk model was further validated in GSE66229 and GSE15459 datasets. In line with that of the training cohort (TCGA-STAD), patients in high-risk group also had worse prognosis in the validation cohorts (Fig. 3D, Additional file 2: Fig. S2A). In GSE66229 dataset, the AUC values of the ROC curve for 1-year, 3-year and 5-year survival were 0.620, 0.625, and 0.601, respectively (Fig. 3E). Corresponding AUC values of 0.620, 0.647, and 0.657 were observed in GSE15459 dataset (Additional file 2: Fig. S2B). The higher the risk score, the worse the survival (Fig. 3F, Additional file 2: Fig. S2C). The heatmaps of the expressions of the 4 genes were shown in Fig. 3F, Additional file 2: Fig. S2C. Consistent with the TCGA-STAD training cohort, the expressions of the NOX4 and FKBP10 were significantly up-regulated, while the expressions of the MAOA and ALDH3A2 were significantly down-regulated in STAD in GSE66229 validation cohort (Additional file 2: Fig. S2D). Next, we systematically analyzed the relationship between the risk score and clinical characteristics in STAD. The risk scores were remarkably higher in patients with H pylori infection, and cystic, mucinous and serous neoplasms (Additional file 3: Fig. S3). Nevertheless, no differences were observed in the mean of risk score among the groups of age, gender, T stage, N stage, M stage, tumor stage, family history of GC, grade, and reflux history (Additional file 3: Fig. S3). The clinical characteristics of the low-risk and high-risk subgroups were then compared, and the difference of Gender (P = 0.006), N stage (P = 0.043), H pylori infection (P = 0.035), disease type (P = 0.046) and survival status (P = 0.008) among the two risk subgroups reached statistical significance (Table 2).

Table 2 Clinical characteristics between low- and high-risk groups

Construction of nomogram

The nomogram integrated the risk score and all important clinical features, which can be used to quantitatively predict the prognosis of patients and provide a reference for clinical decision making. In our study, risk score (P = 0.0005) and age (P = 0.020) were finally identified as prognostic indicators by using univariate and multivariate Cox regression analysis to construct nomogram (Table 3). As a result, a predictive nomogram integrating risk score (a score of 100) and age (a score of 67.5) for prognosis was constructed (Additional file 4: Fig. S4A). ROC curves showed that the AUC values of the nomogram were 0.651, 0.664, and 0.749 for 1-, 3-, and 5-years OS, respectively (Additional file 4: Fig. S4B). The calibration curve showed that the actual survival probabilities at 1-, 3- and 5-year were almost in accordance with the survival probabilities predicted by the nomogram model (Additional file 4: Fig. S4C). The decision curves showed that the nomogram model was better than other factors in predicting the prognosis in STAD (Additional file 4: Fig. S4D).

Table 3 Univariate and multivariate Cox regression analysis of various prognostic parameters in STAD patients

Functional enrichment analysis of the DEGs in high-risk and low-risk groups

We further conducted functional enrichment analyses of the 298 DEGs in high-risk and low-risk groups (Additional file 10: Table S12). GO enrichment analysis indicated that the differential genes annotated to biological processes were involved in extracellular matrix (ECM) organization and extracellular structure organization. Differential genes annotated to cellular component categories were mainly enriched in collagen-containing ECM and collagen trimer. Differential genes annotated to molecular function categories were mainly enriched in ECM structural constituent and collagen binding (Fig. 4A, Additional file 10: Table S13). The top 10 pathways obtained by KEGG analysis were: protein digestion and absorption, proteoglycans in cancer, focal adhesion, human papillomavirus infection, PI3K-Akt signaling pathway, ECM-receptor interaction, cell adhesion molecules, axon guidance, cAMP signaling pathway, and vascular smooth muscle contraction (Fig. 4B, Additional file 10: Table S14). GSEA results showed that risk score was significantly associated with ECM glycoproteins, core matrisome, ECM organization in high-risk group (Fig. 4C). The detailed GSEA results for the high-risk and low-risk groups were presented in Additional file 5: Fig. S5.

Fig. 4
figure 4

Enrichment analysis in high-risk group and the low-risk group. A Circle map. Bands with different colors in the right half circle symbolized 6 significant GO pathways, including biological process (BP), cellular component (CC), and molecular function (MF). The 6 pathways were enriched by genes listed in the left half circle. B Circle map. Bands with different colors in the right half circle symbolized top 10 significant KEGG pathways. The 10 top pathways were enriched by genes listed in the left half circle. C GSEA recognized different gene sets in the high-risk groups

Mitochondrial-related risk score was associated with TME signatures in STAD

Given the TME-associated signal pathway was enriched through functional enrichment analyses, we explored the relationship between the risk score and the TME signatures. As shown in Fig. 5A, risk score was highly positively correlated with stromal score in STAD, and stromal score was higher in high-risk group compared to that in the low-risk group. We further investigated the relationship between risk score and TME components. Our results indicated that risk score had a significant and positive correlation with the expressions of the majority of carcinoma associated fibroblast (CAF) signatures (Fig. 5B), as well as ECM-collagen and matrisome signatures (Fig. 5C, D). Taken together, these results suggested a close relationship between mitochondrial-related risk score and TME signatures in STAD.

Fig. 5
figure 5

Risk score was associated with TME signatures in STAD. A Association between stromal score and risk score and its distribution in the low- and high-risk groups. B Correlation analysis for risk score and the expressions of carcinoma associated fibroblast (CAF) up and down signatures. C Correlation analysis for risk score and the expressions of ECM and collagen signatures. D Correlation analysis for risk score and the expressions of matrisome signatures. P values were showed as: ***P < 0.001

Mitochondrial-related risk score was associated with immune signatures and immunotherapy responses in STAD

The tumor immune microenvironment was closely related with the therapeutic effects and prognosis of patients with malignant tumor. It is reasonable to check the relationship between the risk score and the immune cell infiltration in STAD. The contents of naive B cells, regulatory T cell (Tregs), M0 macrophage, and M2 macrophage were remarkably higher in the high-risk group. In contrast, CD8+ T cells and resting CD4+ T cells were higher in low-risk group (Fig. 6A, Additional file 6: Fig. S6A, B). Consistent with the above results, risk score was positively correlated with the expressions of the majority signatures of M2 macrophage, while negatively correlated with the expressions of the majority signatures of activated CD8+ T cell (Fig. 6B). The correlations between the risk score and the signatures of other immune cells were presented in Additional file 8: Fig. S8. Consistent with the previous studies [28], the ESTIMATE results showed that the patients in high-risk had higher immune score, and significantly lower tumor purity (P = 0.370), than those in the low-risk groups (Fig. 6C, D). Due to the positive correlation between the risk score, and matrisome and CAF signatures, as well as the negative correlation between the risk score and activated CD8+ T cells signatures, we speculated activated CD8+ T cells signatures was negatively correlated with the matrisome and CAF signatures. Interestingly, our results showed that the expression of activated CD8+ T cell signatures were negatively correlated with both matrisome and CAF signatures (Additional file 6: Fig. S6C). Taken together, these results strongly suggested the tumor immunosuppressive microenvironment might contributed to the worse prognosis of the patients with STAD in high-risk group, which needs to be validated in further study.

Fig. 6
figure 6

The different immune profiles between the low- and high-risk groups in the TCGA-STAD dataset. Two risk groups were divided based on the median risk score. A CIBERSORT analysis. B Correlation between risk score and the expressions of activated CD8+ T cell and M2 macrophages signatures. C, D ESTIMATE algorithm. E Expression variation of immune checkpoint. P values were showed as: ns not significant; *P < 0.05; **P < 0.01; ***P < 0.001

Nowadays, immune checkpoint inhibitors were studied and well applied in cancer immunotherapy. In the present study. Our results indicated that 43 immune checkpoints were considerably modulated in high-risk group (Fig. 6E). In addition, risk score was significantly positively correlated with the expression level of 7 immune checkpoints, including CD200, NRP1, TNFSF4, B7-H3, TNFSF18, LAIR1 and OX40 (r > 0.2, Additional file 6: Fig. S6D). Currently, the inhibitors for PD-1 and CTLA-4 are research hotspots in the treatment of advanced STAD. As shown in Fig. 6C, the expressions of PD-1, PD-L1 and CTLA-4 were significantly down-regulated in the high-risk group. Consistently, risk score was significantly negatively correlated with the expressions of PD1, PD-L1 and CTLA-4, respectively (Additional file 6: Fig. S6D).

Given the above results, we further used the TIDE algorithm to evaluate the ability of risk score in predicting the responses to immunotherapy in STAD. Our results showed that risk score had a significantly positive correlation with TIDE score (Additional file 7: Fig. S7), indicating that patients in low-risk group received better response to immunotherapy. The immunotherapy response rate in high-risk group (32.28%) was significantly lower than that in low-risk group (66.67%) (Fig. 7A).

Fig. 7
figure 7

Risk score is a potential biomarker to predict benefits from immune therapies in STAD. A TIDE predicted the proportion of patients with response to immunotherapy in low-risk and high-risk groups. B The proportion of patients with response to immunotherapy in low-risk and high-risk groups in the PRJEB25780 immunotherapy cohort (45 patients with advanced gastric cancer who had received PD-L1 inhibitor treatment). C TIDE predicted the proportion of patients with response to immunotherapy in low-immune score and high-immune score groups. D TIDE predicted the proportion of patients of four groups based on the risk score and immune score with response to immunotherapy. E TIDE predicted the proportion of patients with response to immunotherapy in low-stromal score and high-stromal score groups. F TIDE predicted the proportion of patients of four groups based on the risk score and immune score with response to immunotherapy. G TIDE predicted the proportion of patients with response to immunotherapy in MSS, MSI-L and MSI-H groups. H TIDE predicted the proportion of patients of six groups based on the risk score and microsatellite status with response to immunotherapy. MSS, Microsatellite stability; MSI-L, Microsatellite Instability-Low; MSI-H, Microsatellite Instability-High. P values were showed as: ns not significant; ***P < 0.001

Next, PRJEB25780 cohort (PD-L1 inhibitor treatment for 45 patients with advanced gastric cancer [29]) was applied to validate whether the mitochondrial-related risk signature could accurately predict the responses to immunotherapy for patients with STAD. Consistent with the prediction results by TIDE, the immunotherapy respond rate in high-risk group (13.04%) was significantly lower than that in low-risk group (40.91%) in PRJEB25780 validation cohort (Fig. 7B).

As shown in Fig. 7C, the immunotherapy response rate in the low-immune subgroup (34.59%) was remarkably lower than that in the high-immune subgroup (64.56%). Interestingly, the immunotherapy response rate in low-risk group (45.56%) was remarkably higher than that in high-risk group (20.29%) in the subgroup with low-immune score. Moreover, the immunotherapy response rate in low-risk group (94.20%) was significantly higher than that in high-risk group (41.57%) in the subgroup with high-immune score, strongly suggesting that combined risk score and immune score was a robust indicator to predict the responses to immunotherapy in STAD (Fig. 7D). In addition, the immunotherapy response rate in the high stromal group (37.97%) was significantly lower than that in the low stromal group (61.01%) (Fig. 7E). In the low stromal subgroup, the immunotherapy response rates were 65.35% and 53.45% in low risk and high risk subgroup, respectively, which were similar to the low-stromal group (61.01%), indicating combined risk score and stromal score was not better than stromal score alone in predicting the response to immunotherapy in STAD patients with low-stromal score. However, the immunotherapy response rate in the low-risk + high-stromal group (68.97%) was significantly higher than that in high-stromal group (37.97%), whereas the immunotherapy response rate in the high-risk + high-stromal group (20.00%) was significantly lower than that in the high-stromal group (37.97%), strongly suggesting that combined risk score and stromal score can more accurately predict response to immunotherapy in STAD patients with high-stromal score (Fig. 7F).

The phenotype for microsatellite instability–high (MSI-H) is a distinct tumor subclass that is highly susceptible to immunotherapy. Consistent with the previous studies [30], the immunotherapy response rate MSI-H subgroup (72.73%) was remarkably higher than that in the MSS subgroup (46.30%) and MSI-L subgroup (36.96%) (all P = 0.0005, Fig. 7G). Furthermore, our results also showed that the immunotherapy response rate in low-risk group (65.35%) was remarkably higher than that in high-risk group (29.57%) in the subgroup with MSS. The immunotherapy response rate in low-risk group (57.89%) was also significantly higher than that in high-risk group (22.22%) in the subgroup with MSI-L. The above results strongly suggested that combination of risk-score and MSS/MSI-L can be used as a robust indicator to predict the response to immunotherapy in STAD (Fig. 7H).

Mutation status of STAD patients in high-risk and low-risk groups

Progressive accumulation of mutations throughout life can lead to cancer. Genome sequencing has revolutionized our understanding of somatic mutation in cancer, providing a detailed view of the mutational processes and genes that drive cancer [31]. Therefore, we mapped the mutation landscape of STAD in both the high risk and low risk groups, and analyzed the relationship between risk score and mutation profile. As shown in Fig. 8A, the top 20 high-frequency mutated genes in high-risk and low-risk group were presented. TP53, TTN, MUC16, LRP1B, SYNE1, CSMD3, FAT4, OBSCN, ARID1A, FLG, CSMD1, DNAH5, SPTA1, PCLO and RYR2 were the common high-frequency mutation genes in both groups. In top 5 mutant genes, the mutation rates of TTN and MUC16 were significantly decreased in high-risk group (Fig. 8B).

Fig. 8
figure 8

Mutation status in the high- and the low-risk groups in STAD. A The top 20 genes according to mutation frequency in low and high-risk groups, respectively. B Mutation rate of the top five mutant genes in high-risk and low-risk groups. C Relationship between the risk score and TMB. D Correlation between risk score and TMB score in STAD. E Kaplan–Meier curves of the OS of patients in the high- and low-TMB groups in the TCGA-STAD training cohort. P values were showed as: ns not significant; *P < 0.05; ***P < 0.001

Besides, accumulating evidences supported that TMB functioned as a potentially predictive biomarker for multiple applications, including the biomarker for response to immunotherapy in malignancies [32,33,34,35,36]. Our data showed that TMB in high-risk group was significantly lower than that in low-risk group, and the risk score had a significant negative correlation with TMB in STAD (Fig. 8C, D). Interestingly, the higher TMB tended to have a better OS compared with the lower TMB, but without a statistically significant difference (P = 0.153, Fig. 8E). Moreover, low-risk group tended to have a better OS compared with high-risk group in high TMB subgroup. The similar results were obtained in the low TMB subgroup, but without a statistically significant difference (P = 0.099, Fig. 8E). Taken together, these results strongly suggested that combination of risk score and TMB might be a valuable biomarker for predicting the prognosis for STAD patients (Fig. 8E).

Risk score predicts therapeutic benefits in STAD

To find the potency of risk score as an index for predicting the response to drugs (including chemotherapy, targeted therapy, and immunotherapy) in STAD, we inferred the IC50 value of the 138 drugs in TCGA-STAD patients. We found that patients in high-risk group might be more sensitive to rapamycin, PD-0325901, dasatinib, et al., whereas patients in low-risk group might be more sensitive to AZD7762, CEP-701, methotrexate, et al., which could provide a reliable reference for clinical treatment (Fig. 9). Detailed information on the top 10 sensitive drugs for the high-risk and low-risk subgroup was illustrated in Tables 4 and 5.

Fig. 9
figure 9

Risk score predicts drug therapeutic benefits in STAD. Proportion of normalized IC50 value of the 60 drugs between the low-risk and high-risk groups. P values were showed as: ns not significant; *P < 0.05; **P < 0.01; ***P < 0.001

Table 4 Detailed information of the top 10 sensitivity drugs in high-risk groups
Table 5 Detailed information of the top 10 sensitivity drugs in low-risk groups

Experimental verification of hub gene expression in STAD

To verify the expressions of hub genes in STAD samples, qRT-PCR was conducted on 10 pairs of STAD tumor and paired adjacent normal tissues. As shown in Fig. 10A, the expressions of NOX4 and FKBP10 were significantly higher in STAD tumor tissues than those in the paired adjacent normal tissues, respectively (all P < 0.05). The expressions of ALDH3A2 and MAOA were significantly lower in STAD tumor tissues than those in the paired adjacent normal tissues, respectively (all P < 0.05, Fig. 10A). We also verified the expressions of the four hub genes in a human normal gastric epithelial cell line GES-1 and human gastric cancer cell lines through qRT-PCR. Our results showed that the expressions of NOX4 and FKBP10 were significantly higher than in the normal gastric epithelial cell line GSE-1 (all P < 0.05, Fig. 10B), while the expressions of ALDH3A2 and MAOA were significantly lower than in the normal gastric epithelial cell line GSE-1 (all P < 0.05, Fig. 10B). These results supported our hypothesis and provide solid evidence for the rationality of choosing these four genes for prognostic model construction (see Fig. 11).

Fig. 10
figure 10

Experimental verification of 4 genes expression in STAD. A Expression of 4 genes in 10 paired STAD tissues and normal tissues was evaluated by qRT-PCR. B Expression of 4 genes in a human normal gastric epithelial cell line GSE-1 and human gastric cancer cell lines through qRT-PCR

Fig. 11
figure 11

Graph summarization. The work summary graph of this study

Discussion

GC remains a frequent cancer worldwide with high incidence and mortality globally [1]. Effective biomarkers are still missing even though 3 biomarkers (HER2, MSI-H and PD-L1) have been proven to predict the responses of targeted therapy in GC [37]. Therefore, identifying more effective biomarkers for targeted therapy and prognosis prediction is highly demand. Mitochondria were important pharmacological targets due to their critical roles in cell proliferation and death. The mitochondrial energy metabolisms are now known to be reprogrammed to meet the challenges of high energy demand, with better use of available fuels for malignant cell growth and migration [38]. Thus, mitochondria play a vital and multifunctional role in tumor occurrence and development, and targeting mitochondria provides therapeutic opportunities [39]. A growing body of research showed that mitochondrial-related genes can be used as biomarkers for the diagnosis and treatment of malignancies.

NOX4 has been identified as a biomarker and therapeutic target for a variety of human cancers. NOX4 was upregulated in pancreatic cancer and was involved in the development of pancreatic cancer by promoting cell proliferation, regulating cell metabolism, and mediating angiogenesis, suggesting NOX4 was a potential therapeutic target for pancreatic cancer [40]. Up-regulation of NOX4 predicted worse prognosis and accelerated tumor growth in colorectal carcinoma [41]. Moreover, it has been shown that NOX4 recruited M2-macrophages via ROS/PI3K signaling pathway-dependent cytokines production, thus contributing to the cell division in NSCLC [42]. Consistent with the previous studies, NOX4 in our model was significantly up-regulated in STAD, and high level of NOX4 was associated with worse prognosis in patients with STAD.

It has been reported that a cancer-specific molecular mechanism for NSCLC was related with FKBP10-dependent protein translation. The expression of FKBP10 was positive in cancer lesions [43]. Li et al. reported that FKBP10 silencing decreased the expression of integrin αV and integrin α6, and P-AKT, suggesting that FKBP10 might promote metastasis [44]. FKBP10 was up-regulated in GC tissues and might be a reliable therapeutic target in GC [45]. Consistent with the above studies, our results showed that the high expression of FKBP10 was related with high risk, and predicated poor prognosis of STAD patients. Consistent with the previous studies, FKBP10 in our model was significantly up-regulated in STAD, and high level of FKBP10 was associated with worse prognosis in patients with STAD.

It has been shown that ALDH3A2 was overexpressed in low-grade GC compared with high-grade GC, and patients with low expression of ALDH3A2 had worse OS than those with high ALDH3A2 expression. ALDH3A2 was reported as a reliable biomarker for the immunotherapy, as well as an independent predictor for the prognosis of GC [46]. In renal clear cell carcinoma, low level of ALDH3A2 was related with shorter survival [47]. Consistently, ALDH3A2 was significantly down-regulated in STAD. However, ALDH3A2 couldn’t effectively predict the prognosis for patients with STAD, and this inconsistent effects of ALDH3A2 on prognosis may be caused by the different inclusion and exclusion criteria across different study, which needs more comprehensive investigation.

MAOA exerted different biological effects in different tumors. MAOA was found to be involved in mitochondrial dysfunction, and promoted malignant growth and metastasis in gastric cancer [48]. MAOA promoted prostate cancer progression by increasing cell growth and cancer stem cells, which suggested that MAOA might be a potential therapeutic target for the treatment of prostate cancer [49]. In the present study, MAOA was significantly down-regulated in STAD. However, the high level of MAOA was associated with worse prognosis in patients with STAD, which need to be validated in more samples in further study.

Currently, many biomarkers were applied for prognostic prediction of GC, such as NOX4, FKBP10 and ALDH3A2, but most of them are studied for a single biomarker [44,45,46, 50]. Increasing evidences indicated that prognostic model constructed by multi-genes as a prognostic index was more comprehensive and effective than single gene in kinds of malignancies. For instance, Nie et al. constructed a GC prognosis model based on metabolic signature, which were mainly related to the dysregulation of the metabolic microenvironment [51]. Wu et al. constructed a immune-related prognostic model [52]. As the dysfunction and dysregulation of mitochondria have been associated with cancer, we constructed a STAD prognostic model based on mitochondrial-related genes, which could effectively predict the prognosis for patients with STAD.

The DEGs between the high-risk and low-risk groups were mainly enriched to the extracellular matrix (ECM) and focal adhesion pathway. ECM accumulation was a classical characteristic feature of tumors, and a higher ECM content predicted a poorer prognosis in a broad range of cancer types [53]. The TME is a composition of cancer cells, non-cancerous stromal cells, soluble growth factors, cytokines, proteases, and ECM, which provides essential signals for tumor survival, growth, and acquisition of invasiveness, while hindering antitumor immunity [24, 54, 55]. Fibroblasts constitute one of the most vital cells in the stroma and turn into cancer-associated fibroblasts (CAFs) in TME. CAFs not only play active roles in tumorigenesis and progression both by soluble factors and direct cell-to-cell contact, but also sculpt TME by suppressing anti-tumor immune responses or by recruiting immunosuppressive cells [54]. Matrisome defined as the compendium of genes encoding core ECM proteins, or structural component of the ECM [53]. Consistently, our results indicated that risk score was positively correlated with the CAF signature, ECM signature, and Matrisome signatures. Moreover, risk score had a positive correlation with stromal score and a negative correlation with tumor purity, which could be on behalf of higher infiltration degrees of stromal cells in the TME of the high-risk group.

The immune cells play essential roles in TME. The success of cancer immunotherapy relies on the comprehensive understanding of the tumor microenvironment and immune evasion mechanisms in which the tumor, stroma, and infiltrating immune cells coordinated in a complex network. The main benefit of immunotherapy is to generate memory CD8+ T cells for sustained protection against metastasis and preventing recurrence of the disease [56]. Active immune cells could enter into the tumor parenchyma and perform their anti-tumor function [55]. Therefore, the ultimate goal of immunotherapy is to convert an immunodorminant TME into an immunostimulatory TME, which allows the immune system to clear tumor lesions [56]. Treg cells are involved in tumor progression by inhibiting antitumor immunity. High Treg cell infiltration in the TME was involved in unfavorable prognosis in patients with various types of cancer [53]. M2-polarized macrophages, commonly deemed tumor-associated macrophages (TAMs), were contributors to many pro-tumorigenic outcomes in cancer [57]. Macrophage type 2 (M2) cells, and Tregs cells could make immunologic barriers against CD8+ T cell‐mediated antitumor immune responses [58]. Our results indicated that naive B cells, regulatory T cell (Tregs), M0 macrophage and M2 macrophage were significantly enriched in high-risk group, whereas CD8+ T cells and resting CD4+ T cells were remarkably enriched in low-risk group. Taken together, these results suggested that high-risk group was suffused with immunosuppressive cells such as Tregs, M2 macrophages, producing the immunosuppressive microenvironment to hamper CD8+ T cells-mediated eradication for tumor cells.

Monoclonal antibodies against immunological checkpoint molecules provided a vast breakthrough in cancer therapeutics. For instance, PD-1/PD-L1 and CTLA-4 inhibitors showed promising therapeutic outcomes [59]. In our study, high-risk group had a considerably lower rate of immunotherapy response than that in the low-risk group, which was consistent with the expression levels of PD-1/PD-L1 and CTLA-4 in the high-risk and low-risk groups. Interestingly, combination of risk score and immune score, or stromal score can more accurately predict the responsiveness to immunotherapy of patients with STAD. These findings further demonstrated the effectiveness of the risk score as a biomarker in predicating the response to immunotherapy. The TMB was associated with the formation of neoantigens which activated antitumor immunity, which was a reliable biomarker to predict the response to PD-L1 therapy [60]. In our study, the TMB in high-risk group was lower than that in low-risk group, which strongly suggested that the lower response rate to immunotherapy in high-risk group might be due to the lower TMB.

In addition, risk score might be helpful in screening the therapeutics drugs for patients with STAD. For instance, an independent study showed that high-risk STAD patients showed higher sensitivity to the chemotherapy agents, including rapamycin [61]. Another study found that methotrexate is suitable to inhibit the function of Early B-cell factors (EBFs) in gastric cancer [62]. In the present study, rapamycin, PD-0325901 and dasatinib were found to be more effective for patients in the high-risk group, whereas AZD7762, CEP-701 and methotrexate were predicted to be more effective for patients in the low-risk group. However, the toxicities of the screened drugs was uncertain. For instance, due to the unpredictable cardiac toxicity, the development of AZD7762 was not going forward in patients with advanced solid tumors [63]. 5-fluorouracil, doxorubicin, high-dose methotrexate (FAMTX) schedule was reported to be active in advanced gastric cancer, and the main toxicity was myelosuppression [64]. Oxidative stress is a component of many diseases, including cancer. Although numerous small molecule drugs evaluated as antioxidants have exhibited potential therapeutic ability in preclinical studies, results from clinical trial was disappointed. A greater understanding of the pharmacological mechanisms through which anti oxidative drug act might provide a rational usage would lead to greater therapeutic success in malignancies [65, 66].

Our research has some unique advantages. First, the combination of multigene has robust predictive capability for cancer prognosis than single genes. An integrated mitochondria-related gene prognostic risk model would play more vital roles in the diagnosis and prognosis of STAD patients. Second, the results of the study provide us with a more accessible method to determine whether patients belong to the high- or low-risk group, which is simple and feasible. In addition, we evaluated its predictive value, chemotherapy efficacy, immunotherapy efficacy and immune cell infiltration for patients with STAD, which could provide individualized management and treatment for STAD patients.

Limitations

Although our findings in this study have important clinical consequences, there are still some limitations. Firstly, this is a retrospective study, and an independent prospective cohort is needed to verify the risk model constructed in this study. Secondly, this study heavily relied on datasets and computational predictions while validation component is poor. More experimental studies are needed to validate in further studies. Finally, the carcinogenic effects of the prognostic genes in the model and the mechanisms of interaction between prognostic genes and mitochondrial dysfunction in STAD are mainly unknown and need to be further explored. Based on above information, our future direction will focus on three aspects: (1) Applying mouse model to verify our current hypothesis; (2) Collecting gastric cancer cases and clinical information to validate the risk score model, and compare it with the gold standard of clinical diagnosis, and ensure that the risk model constructed in the present study can be applicable for clinical practice; (3) To screen more novel genes associated with mitochondria in more datasets.

Conclusions

We established a STAD patient risk score model including NOX4, FKBP10, ALDH3A2, and MAOA. Functionally, the risk score was highly correlated to the TME and immune cell infiltration of STAD patients. Combined analysis for risk score and stromal score, or immune score, or MSS/MSI can predict the response to immunotherapy more accurately than single index in STAD. Regarding drug sensitivity, patients in high-risk group was more sensitive to rapamycin, PD-0325901 and dasatinib, whereas patients in low-risk group was more sensitive to AZD7762, CEP-701 and methotrexate. Taken together, our mitochondrial-related risk model could be a reliable prognostic biomarker for personalized treatment of STAD patients.

Availability of data and materials

All data are available in a public, open access repository. R and other custom scripts for analyzing data are available upon reasonable request.

Abbreviations

STAD:

Stomach adenocarcinoma

TME:

Tumor microenvironment

GC:

Gastric cancer

EBV:

Epstein–Barr virus

OS:

Overall survival

ROS:

Reactive oxygen species

ICB:

Immune checkpoint blockade

TMB:

Tumor mutation burden

GSEA:

Gene set enrichment analyses

DEGs:

Differentially expressed genes

K-M:

Kaplan–Meier

C-index:

Concordance index

HPA:

The Human Protein Atlas

NCBI:

National Center for Biotechnology Information

DCA:

Decision curve analysis

GO:

Gene Ontology

KEGG:

Kyoto Encyclopedia of Genes and Genomes

TIICs:

Tumor-infiltrating immune cells

IC50:

50% Inhibiting concentration

GDSC:

Genomics of Drug Sensitivity in Cancer

TIDE:

Tumor immune dysfunction and exclusion

ECM:

Extracellular matrix

CAF:

Carcinoma associated fibroblast

Tregs:

Regulatory T cell

MSS:

Microsatellite stability

MSI-L:

Microsatellite instability-low

MSI-H:

Microsatellite instability-high

HER2:

Human epidermal growth factor receptor 2

PD-L1:

Programmed death-ligand 1

TAMs:

Tumor-associated macrophages

M2:

Macrophage type 2

BP:

Biological process

CC:

Cellular component

MF:

Molecular function

References

  1. Sung H, Ferlay J, Siegel RL, Laversanne M, Soerjomataram I, Jemal A, Bray F. Global cancer statistics 2020 GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin. 2021;71:209–49.

    Article  PubMed  Google Scholar 

  2. Hoft SG, Noto CN, DiPaolo RJ. Two distinct etiologies of gastric cancer: infection and autoimmunity. Front Cell Dev Biol. 2021;9:752346.

    Article  PubMed  PubMed Central  Google Scholar 

  3. Machlowska J, Baj J, Sitarz M, Maciejewski R, Sitarz R. Gastric cancer: epidemiology, risk factors, classification, genomic characteristics and treatment strategies. Int J Mol Sci. 2020;21:4012.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  4. Poorolajal J, Moradi L, Mohammadi Y, Cheraghi Z, Gohari-Ensaf F. Risk factors for stomach cancer: a systematic review and meta-analysis. Epidemiol Health. 2020;42: e2020004.

    Article  PubMed  PubMed Central  Google Scholar 

  5. Ilic M, Ilic I. Epidemiology of stomach cancer. World J Gastroenterol. 2022;28:1187–203.

    Article  PubMed  PubMed Central  Google Scholar 

  6. Tanprasert P, Limpakan S, Chattipakorn SC, Chattipakorn N, Shinlapawittayatorn K. Targeting mitochondria as a therapeutic anti-gastric cancer approach. Apoptosis. 2022;27:163–83.

    Article  CAS  PubMed  Google Scholar 

  7. Lee H. Somatic alterations in mitochondrial DNA and mitochondrial dysfunction in gastric cancer progression. World J Gastroenterol. 2014;20:3950.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  8. Rodrigues T, Ferraz LS. Therapeutic potential of targeting mitochondrial dynamics in cancer. Biochem Pharmacol. 2020;182: 114282.

    Article  CAS  PubMed  Google Scholar 

  9. Ahmadian E, Babaei H, Mohajjel NA, Eftekhari A, Eghbal MA. Venlafaxine-induced cytotoxicity towards isolated rat hepatocytes involves oxidative stress and mitochondrial/lysosomal dysfunction. Adv Pharm Bull. 2016;6:521–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Chodari L, Dilsiz AM, Vahedi P, Alipour M, Vahed SZ, Khatibi S, Ahmadian E, Ardalan M, Eftekhari A. Targeting mitochondrial biogenesis with polyphenol compounds. Oxid Med Cell Longev. 2021;2021:4946711.

    Article  PubMed  PubMed Central  Google Scholar 

  11. Ma X, Gong N, Zhong L, Sun J, Liang XJ. Future of nanotherapeutics: targeting the cellular sub-organelles. Biomaterials. 2016;97:10–21.

    Article  CAS  PubMed  Google Scholar 

  12. Li Z, Liu ZM, Xu BH. A meta-analysis of the effect of microRNA-34a on the progression and prognosis of gastric cancer. Eur Rev Med Pharmacol Sci. 2018;22:8281–7.

    CAS  PubMed  Google Scholar 

  13. Chang W. Inflammation-related factors predicting prognosis of gastric cancer. World J Gastroenterol. 2014;20:4586.

    Article  PubMed  PubMed Central  Google Scholar 

  14. Shao W, Yang Z, Fu Y, Zheng L, Liu F, Chai L, Jia J. The pyroptosis-related signature predicts prognosis and indicates immune microenvironment infiltration in gastric cancer. Front Cell Dev Biol. 2021;9:676485.

    Article  PubMed  PubMed Central  Google Scholar 

  15. Peng Y, Liu C, Li M, Li W, Zhang M, Jiang X, Chang Y, Liu L, Wang F, Zhao Q. Identification of a prognostic and therapeutic immune signature associated with hepatocellular carcinoma. Cancer Cell Int. 2021;21:98.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  16. Khalaf K, Hana D, Chou JT, Singh C, Mackiewicz A, Kaczmarek M. Aspects of the tumor microenvironment involved in immune resistance and drug resistance. Front Immunol. 2021;12: 656364.

    Article  PubMed  PubMed Central  Google Scholar 

  17. Petitprez F, Meylan M, de Reyniès A, Sautès-Fridman C, Fridman WH. The tumor microenvironment in the response to immune checkpoint blockade therapies. Front Immunol. 2020;11:784.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Samstein RM, Lee CH, Shoushtari AN, Hellmann MD, Shen R, Janjigian YY, Barron DA, Zehir A, Jordan EJ, Omuro A, et al. Tumor mutational load predicts survival after immunotherapy across multiple cancer types. Nat Genet. 2019;51:202–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Goldman MJ, Craft B, Hastie M, Repečka K, McDade F, Kamath A, Banerjee A, Luo Y, Rogers D, Brooks AN, et al. Visualizing and interpreting cancer genomics data via the Xena platform. Nat Biotechnol. 2020;38:675–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Rath S, Sharma R, Gupta R, Ast T, Chan C, Durham TJ, Goodman RP, Grabarek Z, Haas ME, Hung W, et al. MitoCarta3.0: an updated mitochondrial proteome now with sub-organelle localization and pathway annotations. Nucleic Acids Res. 2021;49:D1541-7.

    Article  CAS  PubMed  Google Scholar 

  21. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, Paulovich A, Pomeroy SL, Golub TR, Lander ES, et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc Natl Acad Sci U S A. 2005;102:15545–50.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Mootha VK, Lindgren CM, Eriksson KF, Subramanian A, Sihag S, Lehar J, Puigserver P, Carlsson E, Ridderstrale M, Laurila E, et al. PGC-1alpha-responsive genes involved in oxidative phosphorylation are coordinately downregulated in human diabetes. Nat Genet. 2003;34:267–73.

    Article  CAS  PubMed  Google Scholar 

  23. Zhao HB, Zeng YR, Han ZD, Zhuo YJ, Liang YK, Hon CT, Wan S, Wu S, Dahl D, Zhong WD, et al. Novel immune-related signature for risk stratification and prognosis in prostatic adenocarcinoma. Cancer Sci. 2021;112:4365–76.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Wu J, Li L, Zhang H, Zhao Y, Zhang H, Wu S, Xu B. A risk model developed based on tumor microenvironment predicts overall survival and associates with tumor immunity of patients with lung adenocarcinoma. Oncogene. 2021;40:4413–24.

    Article  CAS  PubMed  Google Scholar 

  25. Ru B, Wong CN, Tong Y, Zhong JY, Zhong S, Wu WC, Chu KC, Wong CY, Lau CY, Chen I, et al. TISIDB: an integrated repository portal for tumor-immune system interactions. Bioinformatics. 2019;35:4200–2.

    Article  CAS  PubMed  Google Scholar 

  26. Cerami E, Gao J, Dogrusoz U, Gross BE, Sumer SO, Aksoy BA, Jacobsen A, Byrne CJ, Heuer ML, Larsson E, et al. The cBio cancer genomics portal: an open platform for exploring multidimensional cancer genomics data. Cancer Discov. 2012;2:401–4.

    Article  PubMed  Google Scholar 

  27. Gao J, Aksoy BA, Dogrusoz U, Dresdner G, Gross B, Sumer SO, Sun Y, Jacobsen A, Sinha R, Larsson E, et al. Integrative analysis of complex cancer genomics and clinical profiles using the cBioPortal. Sci Signal. 2013;6: l1.

    Article  Google Scholar 

  28. Gong Z, Zhang J, Guo W. Tumor purity as a prognosis and immunotherapy relevant feature in gastric cancer. Cancer Med. 2020;9:9052–63.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Lin Y, Pan X, Zhao L, Yang C, Zhang Z, Wang B, Gao Z, Jiang K, Ye Y, Wang S, et al. Immune cell infiltration signatures identified molecular subtypes and underlying mechanisms in gastric cancer. Npj Genom Med. 2021;6:83.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  30. Chao J, Fuchs CS, Shitara K, Tabernero J, Muro K, Van Cutsem E, Bang YJ, De Vita F, Landers G, Yen CJ, et al. Assessment of pembrolizumab therapy for the treatment of microsatellite instability-high gastric or gastroesophageal junction cancer among patients in the KEYNOTE-059, KEYNOTE-061, and KEYNOTE-062 clinical trials. Jama Oncol. 2021;7:895–902.

    Article  PubMed  PubMed Central  Google Scholar 

  31. Martincorena I, Campbell PJ. Somatic mutation in cancer and normal cells. Science. 2015;349:1483–9.

    Article  CAS  PubMed  Google Scholar 

  32. Marcus L, Fashoyin-Aje LA, Donoghue M, Yuan M, Rodriguez L, Gallagher PS, Philip R, Ghosh S, Theoret MR, Beaver JA, et al. FDA approval summary: pembrolizumab for the treatment of tumor mutational burden-high solid tumors. Clin Cancer Res. 2021;27:4685–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. FDA approves first drug for cancers with a high tumor mutational burden. American Cancer Society. 2020.

  34. Kim JY, Kronbichler A, Eisenhut M, Hong SH, van der Vliet HJ, Kang J, Shin JI, Gamerith G. Tumor mutational burden and efficacy of immune checkpoint inhibitors: a systematic review and meta-analysis. Cancers. 2019;11:1798.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Owada-Ozaki Y, Muto S, Takagi H, Inoue T, Watanabe Y, Fukuhara M, Yamaura T, Okabe N, Matsumura Y, Hasegawa T, et al. Prognostic impact of tumor mutation burden in patients with completely resected non-small cell lung cancer: brief report. J Thorac Oncol. 2018;13:1217–21.

    Article  CAS  PubMed  Google Scholar 

  36. Riviere P, Goodman AM, Okamura R, Barkauskas DA, Whitchurch TJ, Lee S, Khalid N, Collier R, Mareboina M, Frampton GM, et al. High tumor mutational burden correlates with longer survival in immunotherapy-naïve patients with diverse cancers. Mol Cancer Ther. 2020;19:2139–45.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  37. Choi S, Park S, Kim H, Kang SY, Ahn S, Kim K. Gastric cancer: mechanisms, biomarkers, and therapeutic approaches. Biomedicines. 2022;10:543.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Roth KG, Mambetsariev I, Kulkarni P, Salgia R. The mitochondrion as an emerging therapeutic target in cancer. Trends Mol Med. 2020;26:119–34.

    Article  CAS  PubMed  Google Scholar 

  39. Zong WX, Rabinowitz JD, White E. Mitochondria and cancer. Mol Cell. 2016;61:667–76.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  40. Bi Y, Lei X, Chai N, Linghu E. NOX4: a potential therapeutic target for pancreatic cancer and its mechanism. J Transl Med. 2021;19:1.

    Article  Google Scholar 

  41. Lin X, Yang L, Fu S, Lin W, Gao Y, Chen H, Ge Z. Overexpression of NOX4 predicts poor prognosis and promotes tumor progression in human colorectal cancer. Oncotarget. 2017;8:33586–600.

    Article  PubMed  PubMed Central  Google Scholar 

  42. Zhang J, Li H, Wu Q, Chen Y, Deng Y, Yang Z, Zhang L, Liu B. Tumoral NOX4 recruits M2 tumor-associated macrophages via ROS/PI3K signaling-dependent various cytokine production to promote NSCLC growth. Redox Biol. 2019;22: 101116.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Ramadori G, Ioris RM, Villanyi Z, Firnkes R, Panasenko OO, Allen G, Konstantinidou G, Aras E, Brenachot X, Biscotti T, et al. FKBP10 regulates protein translation to sustain lung cancer growth. Cell Rep. 2020;30:3851–63.

    Article  CAS  PubMed  Google Scholar 

  44. Gong LB, Zhang C, Yu RX, Li C, Fan YB, Liu YP, Qu XJ. FKBP10 acts as a new biomarker for prognosis and lymph node metastasis of gastric cancer by bioinformatics analysis and in vitro experiments. Onco Targets Ther. 2020;13:7399–409.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  45. Liang L, Zhao K, Zhu JH, Chen G, Qin XG, Chen JQ. Comprehensive evaluation of FKBP10 expression and its prognostic potential in gastric cancer. Oncol Rep. 2019;42:615–28.

    CAS  PubMed  PubMed Central  Google Scholar 

  46. Yin Z, Wu D, Shi J, Wei X, Jin N, Lu X, Ren X. Identification of ALDH3A2 as a novel prognostic biomarker in gastric adenocarcinoma using integrated bioinformatics analysis. Bmc Cancer. 2020;20:1.

    Google Scholar 

  47. Zhao Y, Tao Z, Chen X. A three-metabolic-genes risk score model predicts overall survival in clear cell renal cell carcinoma patients. Front Oncol. 2020;10: 570281.

    Article  PubMed  PubMed Central  Google Scholar 

  48. Chen L, Guo L, Sun Z, Yang G, Guo J, Chen K, Xiao R, Yang X, Sheng L. Monoamine oxidase A is a major mediator of mitochondrial homeostasis and glycolysis in gastric cancer progression. Cancer Manag Res. 2020;12:8023–35.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Liao C, Lin T, Li P, Geary LA, Chen K, Vaikari VP, Wu JB, Lin C, Gross ME, Shih JC. Loss of MAOA in epithelia inhibits adenocarcinoma development, cell proliferation and cancer stem cells in prostate. Oncogene. 2018;37:5175–90.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Tang CT, Gao YJ, Ge ZZ. NOX4, a new genetic target for anti-cancer therapy in digestive system cancer. J Dig Dis. 2018;19:578–85.

    Article  CAS  PubMed  Google Scholar 

  51. Nie Y, Liu L, Liu Q, Zhu X. Identification of a metabolic-related gene signature predicting the overall survival for patients with stomach adenocarcinoma. PeerJ. 2021;9: e10908.

    Article  PubMed  PubMed Central  Google Scholar 

  52. Wu M, Xia Y, Wang Y, Fan F, Li X, Song J, Ding J. Development and validation of an immune-related gene prognostic model for stomach adenocarcinoma. Bioscience Rep. 2020;40:BSR20201012.

    Article  CAS  Google Scholar 

  53. Izzi V, Davis MN, Naba A. Pan-cancer analysis of the genomic alterations and mutations of the matrisome. Cancers. 2020;12:2046.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Gok Yavuz B, Gunaydin G, Kosemehmetoglu K, Karakoc D, Ozgur F, Guc D. The effects of cancer-associated fibroblasts obtained from atypical ductal hyperplasia on anti-tumor immune responses. Breast J. 2018;24:1099–101.

    Article  PubMed  Google Scholar 

  55. Labani-Motlagh A, Ashja-Mahdavi M, Loskog A. The tumor microenvironment: a milieu hindering and obstructing antitumor immune responses. Front Immunol. 2020;11:940.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  56. Locy H, de Mey S, de Mey W, De Ridder M, Thielemans K, Maenhout SK. Immunomodulation of the tumor microenvironment: turn foe into friend. Front Immunol. 2018;9:2909.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  57. Boutilier AJ, Elsawa SF. Macrophage polarization states in the tumor microenvironment. Int J Mol Sci. 2021;22:6995.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Farhood B, Najafi M, Mortezaee K. CD8(+) cytotoxic T lymphocytes in cancer immunotherapy: a review. J Cell Physiol. 2019;234:8509–21.

    Article  CAS  PubMed  Google Scholar 

  59. Darvin P, Toor SM, Sasidharan Nair V, Elkord E. Immune checkpoint inhibitors: recent progress and potential biomarkers. Exp Mol Med. 2018;50:1–11.

    Article  PubMed  Google Scholar 

  60. Xu F, Huang X, Li Y, Chen Y, Lin L. m6A-related lncRNAs are potential biomarkers for predicting prognoses and immune responses in patients with LUAD. Mol Ther Nucleic Acids. 2021;24:780–91.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Tong X, Yang X, Tong X, Zhai D, Liu Y. Complement system-related genes in stomach adenocarcinoma: Prognostic signature, immune landscape, and drug resistance. Front Genet. 2022;13: 903421.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  62. Wang Q, Liang J, Hu X, Gu S, Xu Q, Yan J. Early B-cell factors involve in the tumorigenesis and predict the overall survival of gastric cancer. Biosci Rep. 2021;41:20210055.

    Article  Google Scholar 

  63. Sausville E, Lorusso P, Carducci M, Carter J, Quinn MF, Malburg L, Azad N, Cosgrove D, Knight R, Barker P, et al. Phase I dose-escalation study of AZD7762, a checkpoint kinase inhibitor, in combination with gemcitabine in US patients with advanced solid tumors. Cancer Chemother Pharmacol. 2014;73:539–49.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  64. Bar SG, Tsalic M, Gaitini D, Steiner M, Haim N. Etoposide, doxorubicin and cisplatin alternating with 5-fluorouracil, doxorubicin and high-dose methotrexate in patients with advanced adenocarcinoma of the stomach or the gastroesophageal junction. J Chemother. 2002;14:623–6.

    Article  Google Scholar 

  65. Forman HJ, Zhang H. Author correction: targeting oxidative stress in disease: promise and limitations of antioxidant therapy. Nat Rev Drug Discov. 2021;20:652.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  66. Ahmadian E, Eftekhari A, Kavetskyy T, Khosroushahi AY, Turksoy VA, Khalilov R. Effects of quercetin loaded nanostructured lipid carriers on the paraquat-induced toxicity in human lymphocytes. Pestic Biochem Physiol. 2020;167: 104586.

    Article  CAS  PubMed  Google Scholar 

Download references

Acknowledgements

We are grateful to the patients for their contributions to this study.

Funding

This study was supported by the National Natural Science Foundation of China (No. 81902513 and 82100680) and Applied Basic Research Project of Shanxi Province (No. 202103021224228 and 20210302124376), and The Science/Technology Project of Sichuan Province (No. 2021YFS0161).

Author information

Authors and Affiliations

Authors

Contributions

JZ, PY, and HW designed the study, JC, HW, JW, ML, WZ, YH, XZ, JX, and LL conducted the experiments and data analysis. JZ, PY and JC wrote and revised the manuscript. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Pengfei Yu or Jianjun Zhu.

Ethics declarations

Ethics approval and consent to participate

Not available.

Consent for publication

All authors agree to publish.

Competing interests

The authors declare that they have no competing interests for this work.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Figure S1.

Identification of DEGs related to mitochondrion and functional enrichment analysis in STAD. A The GO analysis of 183 mitochondrial-related DEGs, including biological process (BP), cellular component (CC), and molecular function (MF). B The KEGG analysis of 183 mitochondrial-related DEGs. C Protein expressions of the 3 prognosis-related genes in gastric cancer tissues and normal gastric tissues from HPA database. D Kaplan–Meier curves for OS of patients with high or low expression of each prognosis-related genes.

Additional file 2: Figure S2.

Assessing the performance of the prognostic risk model in the validation cohort. A Kaplan–Meier curves of the overall survival (OS) in the validation cohort GSE15459. B ROC curves for 1-, 3-, and 5-year OS of the prognostic risk model in the validation cohort GSE15459. C-D Distribution of risk score, survival status (red dots indicate dead, blue dots indicate alive) and the four genes expression heat map in the validation cohort GSE15459. E Gene expression of 4 prognosis-related genes in GSE66229. P values were showed as: ns not significant; *p < 0.05; **p < 0.01; ***p < 0.001

Additional file 3: Figure S3.

The relationships between the risk score and clinical characteristics of STAD patients. Age, Gender, T stage, N stage, M stage, Tumor stage, Family history of stomach cancer, H pylori infection, Grade, Reflux, Disease type. NX, lymph nodal status could not be determined; MX, metastatic status could not be determined; GX, tumor grade could not be determined; N/A, not available. P values were showed as: ns not significant; *p < 0.05; **p < 0.01; ***p < 0.001

Additional file 4: Figure S4.

Construction of nomogram using TCGA-STAD cohort. A A nomogram was constructed based on risk score and related clinical characteristics. B ROC curves and AUC for 1-, 3-, and 5-year OS of the nomogram. C Calibration curves of 1-, 3-, and 5-year OS in the nomogram and ideal model. D DCA results of risk score and clinical characteristics in the TCGA cohort.

Additional file 5: Figure S5.

Results of GSEA analysis in low-risk and high-risk groups in STAD. A The GSEA findings of the c2 reference gene sets for high-risk groups. B The GSEA findings of the c2 reference gene sets for low-risk groups.

Additional file 6: Figure S6.

The condition of Immune infiltration in low-risk and high-risk group. A The proportion of 22 immune cells quantified by CIBERSORT algorithm in low-risk group. B The proportion of 22 immune cells quantified by CIBERSORT algorithm in high-risk group. C Correlation analysis between the activated CD8 + T cell signatures and CAF signatures, as well as that between the activated CD8 + T cell signatures and matrisome signatures. P values were showed as: ns not significant; *p < 0.05; **p < 0.01; ***p < 0.001

Additional file 7: Figure S7.

Correlation analysis for risk score and TIDE score in STAD.

Additional file 8: Figure S8.

Correlation analysis for risk score and 14 immune cells signatures in STAD.

Additional file 9: Figure S9.

Risk score predicts drug therapeutic benefits in STAD. Proportion of normalized IC50 value of the 138 drugs between the low-risk and high-risk groups. P values were showed as: ns not significant; *p < 0.05; **p < 0.01; ***p < 0.001

Additional file 10: Tables S1–S14.

Table S1. The list of mitochondrial-related genes. Table S2. The list of carcinoma associated fibroblast up signatures (n=24). Table S3. The list of carcinoma associated fibroblast down signatures (n=24). Table S4. The list of ECM and Collagen signatures (n=225). Table S5. The list of matrisome signatures (1026). Table S6. The primer sequences. Table S7. The list of all DEGs in tumor and normal groups (n=2381). Table S8. The list of protein-coding DEGs in tumor and normal groups (n=2145). Table S9. The list of mitochondrial-related DEGs (n=183). Table S10. GO terms in tumor and normal groups. Table S11. KEGG terms in tumor and normal groups. Table 12. The list of DEGs in high and low risk groups (n=298). Table S13. GO terms in high and low risk groups. Table S14. KEGG terms in high and low risk groups .

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Chang, J., Wu, H., Wu, J. et al. Constructing a novel mitochondrial-related gene signature for evaluating the tumor immune microenvironment and predicting survival in stomach adenocarcinoma. J Transl Med 21, 191 (2023). https://doi.org/10.1186/s12967-023-04033-6

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12967-023-04033-6

Keywords