Identification of molecular features correlating with tumor immunity in gastric cancer by multi-omics data analysis.

BACKGROUND: Although immunotherapy has achieved success in treating various refractory malignancies including gastric cancers (GCs) with DNA mismatch repair deficiency, only a subset of cancer patients are responsive to immunotherapy. Therefore, the identification of useful biomarkers or interventional targets for improving cancer immunotherapy response is urgently needed.
METHODS: We investigated the associations between various molecular features and immune signatures using three multi-omics GC datasets. These molecular features included genes, microRNAs (miRNAs), long non-coding RNAs (lncRNAs), proteins, and pathways, and the immune signatures included CD8+ T cell infiltration, immune cytolytic activity (ICA), and PD-L1 expression. Moreover, we investigated the association between gene mutations and survival prognosis in a gastrointestinal (GI) cancer cohort receiving immunotherapy and two GC cohorts not receiving such a therapy.
RESULTS: The mutations of some important oncogenes and tumor suppressor genes were appreciably associated with immune signatures in GC, including PIK3CA, MTOR, RNF213, TP53, ARID1A, PTEN, ATM, and CDH1. Moreover, a number of genes exhibited a significant expression correlation with immune signatures in GC, including CXCL9, CXCL13, CXCR6, CCL5, GUCY2C, MAP3K9, NEK3, PAK6, STK35, and WNK2. We identified several proteins whose expression had a significant positive correlation with immune signatures in GC. These proteins included caspase-7, PI3K-p85, PREX1, Lck, Bcl-2, and transglutaminase. In contrast, acetyl-CoA carboxylase (ACC) had a significant inverse expression correlation with immune signatures in GC, suggesting that inhibiting ACC could enhance antitumor immunity in GC. Furthermore, we identified numerous miRNAs and lncRNAs with a significant expression correlation with GC immunity, including hsa-miR-150, 155, 142, 342, 146, 101, 511, 29, AC022706.1, LINC01871, and AC006033.2. We also identified numerous cancer-associated pathways whose activity was associated with GC immunity, including mTOR, PI3K-AKT, MAPK, HIF-1, and VEGF signaling pathways. Interestingly, we found seven genes (ARID1A, BCOR, MTOR, CREBBP, SPEN, NOTCH4, and TET1) whose mutations were associated with better OS in GI cancer patients receiving anti-PD-1/PD-L1 immunotherapy but were not associated with OS in GC patients without immunotherapy.
CONCLUSIONS: The molecular features significantly associated with GC immunity could be useful biomarkers for stratifying GC patients responsive to immunotherapy or intervention targets for promoting antitumor immunity and immunotherapy response in GC. 2020 Annals of Translational Medicine. All rights reserved.

Introduction

Gastric cancer (GC) is the second leading cause of cancer deaths in the world and particularly prevails in Asian countries (1). Traditional therapeutic strategies (surgery, chemotherapy, radiotherapy, and targeted therapy) often have limited efficacy against advanced GCs (2). Recently, with the advance of immunotherapy in treating various advanced malignancies including the cancers with mismatch repair deficiency (3-7), immunotherapy is becoming a potentially promising treatment strategy for GC, as evidenced by that pembrolizumab and nivolumab (immune checkpoint inhibitors) are being clinically used to treat microsatellite-instable gastrointestinal (GI) cancers (7). However, currently, immunotherapeutic strategies are beneficial to only a subset of cancer patients, with a large proportion of patients failing to respond to such therapies (8). To this end, certain genetic or genomic features have been identified that are associated with cancer immunotherapy response, such as PD-L1 expression (9), tumor mutation burden (TMB) (10), deficient mismatch repair (dMMR) (11), and TP53 mutations (12). In addition, it has been shown that the tumor immune microenvironment (TIME) plays a crucial role in determining the cancer immunotherapy response (13). The “hot” tumors with dense T cell infiltration in TIME are more likely to respond to immunotherapy than the “cold” tumors lack of immune cell infiltration (14). To improve the cancer therapy response, the combination of different immunotherapeutic strategies and the combination of immunotherapy with traditional therapeutic methods have been widely investigated (15-19).

Benefit from the rapid advances in genomics technologies, many cancer genomics data are publicly available, such as TCGA (https://portal.gdc.cancer.gov/) and ICGC (https://dcc.icgc.org/). Several studies have explored tumor immunity in GC using publicly available cancer genomic datasets (20-22). Zeng et al. identified three GC subtypes based on the tumor microenvironment (TME) cell infiltration patterns (20). Our previous study revealed that the TP53 mutation contributed to an immune-suppressive TME in GC (21). Park et al. classified GC patients into three subtypes with significant prognostic differences and found that the GC subtype with the best clinical outcome was enriched in the PD-L1 signaling pathway (22).

Despite these prior explorations (20-22), the identification of genetic and genomic features significantly correlating with GC immunity and immunotherapy response is worth further investigation. For example, although the Epstein-Barr virus positive (EBV+) and the microsatellite instability (MSI) are associated with a high GC immunity and an active immunotherapy response, many EBV- and microsatellite-stable (MSS) GCs equally exhibit highly enriched immune signatures and favorable immunotherapy response (21,22). This suggests that there are other factors correlating with GC immunity and immunotherapy response. In this study, we aimed to explore the associations of various molecular features with GC immunity based on multi-omics datasets from the TCGA GC cohort (23) and the Asian Cancer Research Group (ACRG) GC cohort (24). These molecular features included genes, microRNAs (miRNAs), long non-coding RNAs (lncRNAs), proteins, and pathways. We identified the molecular features that were significantly associated with GC immunity and immunotherapy response. These molecular features may have clinical implications for GC immunotherapy. We present the following article in accordance with the MDAR reporting checklist (available at http://dx.doi.org/10.21037/atm-20-922).

Methods

Datasets

We obtained multi-omics (gene somatic mutation, gene expression, protein expression, miRNA expression, and lncRNA expression) and clinical data for the TCGA GC cohort from the GDC data portal (https://portal.gdc.cancer.gov/). The multi-omics (gene somatic mutation and gene expression) and clinical data for the ACRG GC cohort were from the associated publication by ACRG (24). The gene somatic mutation and clinical data for a GI cancer cohort [Samstein cohort (25)] with immunotherapy were from the associated publication (25). In addition, we obtained cancer-associated pathways and their gene sets from KEGG (26).

Quantification of immune signature enrichment levels in GC samples

We analyzed three immune signatures, including CD8+ T cells, immune cytolytic activity (ICA), and PD-L1 expression. The enrichment level of an immune signature in a GC sample was quantified as the mean expression level of immune signature marker genes. The marker genes of three immune signatures included CD8A for CD8+ T cells, GZMA and PRF1 for ICA, and PD-L1 for PD-L1.

Associations of immune signature enrichment levels with various molecular features in GC

We identified the genes whose mutations were significantly associated with immune signature enrichment levels (ISELs) using the Mann-Whitney U test and the genes, proteins, miRNAs, and lncRNAs whose expression levels were significantly correlated with ISELs using the Pearson’s correlation test. We quantified the activity of a pathway in a GC sample by the single-sample gene-set enrichment analysis (ssGSEA) (27) score of the gene set in the pathway and identified the cancer-associated pathways whose activity was significantly associated with ISELs using the Spearman’s correlation test. The false discovery rate (FDR) (28) was used to adjust for multiple tests.

For each GC sample, we defined its TMB as the total somatic mutation count in the sample and its tumor aneuploidy level (TAL) as the tumor ploidy score calculated by ABSOLUTE (29). We divided all GC samples into the high-TMB/TAL (> median) group and the low-TMB/TAL (< median) group and compared their ISELs using the Mann-Whitney U test.

Gene-set enrichment analysis

We used GSEA (30) to identify the KEGG (26) pathways that were significantly associated with the genes having an important expression correlation with ISELs using a threshold of FDR <0.05.

Survival analysis

We compared the overall survival (OS) between gene-mutated and gene-wildtype cancers in the three cohorts (TCGA-GC, ACRG-GC, and Samstein). We used Kaplan-Meier survival curves to show the survival time differences and the log-rank test to assess the significance of survival time differences.

Results

Identification of genes whose mutations are associated with antitumor immune response in GC

Based on TCGA-GC cohort, we identified eight genes whose mutations were significantly associated with heightened CD8+ T cell infiltration levels in GC (Mann-Whitney U test, FDR <0.05) (). The eight genes included PIK3CA, ABCB4, SOGA2, ARID1B, CDH1, ABCA12, TENM2, and YLPM1. In contrast, TP53 mutations were associated with reduced CD8+ T cell infiltration levels in GC (Mann-Whitney U test, FDR =0.042). This is consistent with our previous study showing that TP53 mutations inhibited GC immunity (21). Moreover, we identified 106 genes whose mutations were associated with elevated immune cytolytic activity (ICA) in GC and one gene (TP53) whose mutations were associated with reduced ICA (Mann-Whitney U test, FDR <0.05) (https://cdn.amegroups.cn/static/application/3378c85657fe63fe5b598b068b3f0bbb/10.21037atm-20-922-1.pdf). The 106 genes included the aforementioned eight genes whose mutations were positively associated with CD8+ T cell infiltration levels in GC. Notably, some well-recognized cancer-associated genes were presented in the list of 106 genes, including tumor suppressor genes ARID1A, ATM, CDH1, FBXW7, and PTEN, oncogenes CIC, CREBBP, MLL, NUP98, PIK3CA, and RNF213, and protein kinase encoding genes ERBB3, MTOR, and WNK1. These results indicate that the mutations in a wide range of genes are associated with the enhanced antitumor immunity in GC. However, the mutation of TP53, one of the most important tumor suppressor genes, appears to suppress antitumor immunity in GC. The associations between the mutations of genes (ARID1A, PIK3CA, and TP53) and GC immunity were confirmed to be significant in ACRG-GC cohort (Mann-Whitney U test, FDR <0.05).

Figure 1

Associations of gene mutations with antitumor immune response and immunotherapy response in GC. (A) Eight genes whose mutations are associated with the elevated CD8+ T cell infiltration levels in TCGA-GC cohort (Mann-Whitney U test, FDR <0.05). (B) Kaplan-Meier survival curves show that the mutations of seven genes are associated with better overall survival (OS) in Samstein cohort with anti-PD-1/PD-L1 immunotherapy while they showed no a significant correlation with OS in either of TCGA- and ACRG-GC cohorts not receiving such a therapy (log-rank test, P<0.05). (C) The mutations of seven genes are consistently associated with the elevated expression of PD-L1 in TCGA-GC cohort (Student’s t test, FDR <0.05). (D) TP53-mutated gastrointestinal (GI) cancer patients have significantly worse OS than TP53-wildtype patients in Samstein cohort (log-rank test, P<0.05) while the TP53 mutation has no a significant correlation with OS in either of TCGA- and ACRG-GC cohorts, and the TP53 mutation is associated with reduced PD-L1 expression in GC. *, FDR <0.05; **, FDR <0.01; ***, FDR <0.001. FDR, false discovery rate. They also apply to following figures.

We identified 180 genes whose mutations were associated with elevated PD-L1 expression in GC (https://cdn.amegroups.cn/static/application/3378c85657fe63fe5b598b068b3f0bbb/10.21037atm-20-922-1.pdf). Interestingly, a large proportion (86%) of the genes whose mutations were associated with elevated ICA were contained in the set of 180 genes. This indicates that the mutations in many genes may result in increased antitumor immune response and enhanced protumor immunosuppression simultaneously.

We further explored the associations between the mutations of these genes and OS in TCGA-GC, ACRG-GC, and Samstein cohorts. We found that the mutations of MLL3, FAT2, FAT3, ARID2, and CSMD1 were consistently associated with better OS in ACRG-GC cohort (log-rank test, P<0.05) (Figure S1A). In TCGA-GC cohort, the mutations of 18 genes were associated with better OS (Figure S1B). The positive associations between the mutations of these genes and OS could be attributed to the elevated antitumor immune activity in the GCs harboring these mutations. Interestingly, we found seven genes (ARID1A, BCOR, MTOR, CREBBP, SPEN, NOTCH4, and TET1) whose mutations were associated with better OS in Samstein cohort but did not show a significant correlation with OS in either of TCGA- and ACRG-GC cohorts ( and Figure S1C). A possible explanation of this disparity could be that the GI cancer patients with the mutations of these genes were more responsive to anti-PD-1/PD-L1 immunotherapy than those without such mutations. This explanation was bolstered by the fact that the mutations of these genes were consistently associated with the elevated expression of PD-L1 in TCGA-GC cohort (). Interestingly, the TP53-mutated GI cancer patients had a significantly worse OS than the TP53-wildtype patients in Samstein cohort while the TP53 mutation had no a significant correlation with OS in either TCGA- or ACRG-GC cohort (). Again, this could be attributed to the more active response to anti-PD-1/PD-L1 immunotherapy in TP53-wildtype GI cancer patients since the TP53 mutation was associated with reduced PD-L1 expression in GC ().

Therefore, the identification of genes whose mutations were significantly associated with tumor immune signatures may provide predictive biomarkers for immunotherapy response in GC.

Identification of genes whose expression levels are associated with antitumor immune response in GC

We found 136 and 123 genes whose expression levels exhibited a strong positive correlation with CD8+ T cell enrichment levels in TCGA-GC cohort and ACRG-GC cohort, respectively (Pearson correlation coefficient r>0.7) (https://cdn.amegroups.cn/static/application/4b0111ff02b66c45bc28328104e88994/10.21037atm-20-922-2.pdf). A total of 76 genes showed a strong positive expression correlation with CD8+ T cell enrichment levels in both cohorts (https://cdn.amegroups.cn/static/application/4b0111ff02b66c45bc28328104e88994/10.21037atm-20-922-2.pdf), most of which were immune relevant, such as three granzyme genes (GZMA, GZMH, and GZMK), and numerous cytokine and cytokine receptor genes, including CXCL9, CXCL13, CXCR6, CCL5, IL2RB, IL10RA, and IL21R. GSEA (30) identified 12 KEGG (26) pathways significantly associated with the 76 genes, including cytokine-cytokine receptor interaction, natural killer cell mediated cytotoxicity, chemokine signaling, T cell receptor signaling, graft-versus-host disease, cell adhesion molecules, Jak-STAT signaling, hematopoietic cell lineage, primary immunodeficiency, prion diseases, allograft rejection, and type I diabetes mellitus (). Evidently, all these pathways are immune relevant.

Figure 2

Genes whose expression levels are significantly associated with CD8+ T cell enrichment levels or PD-L1 expression levels and their associated pathways. (A) GSEA (30) identified 12 KEGG (26) pathways significantly associated with the 76 genes which have a strong positive expression correlation with CD8+ T cell enrichment levels in both TCGA- and ACRG-GC cohorts [Pearson correlation coefficient (r) >0.7]. The FDR for each pathway is shown in its bar. (B) Four genes (GBP1, GBP5, IDO1, and WARS) whose expression levels exhibit a strong positive correlation with PD-L1 expression levels in both TCGA- and ACRG-GC cohorts (r>0.7). The Pearson's correlation test P values and correlation coefficients are shown. r: Pearson correlation coefficient. It also applies to following figures.

We identified 146, 734, and 83 genes whose expression levels had a significant inverse correlation with CD8+ T cell enrichment levels in TCGA-GC cohort, ACRG-GC cohort, and both cohorts, respectively (r<−0.3) (https://cdn.amegroups.cn/static/application/4b0111ff02b66c45bc28328104e88994/10.21037atm-20-922-2.pdf). Among the 83 genes, six (GUCY2C, MAP3K9, NEK3, PAK6, STK35, and WNK2) encodes protein kinases. Considering that kinase inhibitors are a class of intensively investigated anticancer drugs that have been widely used in clinics or clinical trials (31), these protein kinases could be potentially useful targets for cancer therapy in that inhibiting them may promote antitumor immune response. GSEA identified one KEGG pathway (steroid biosynthesis) significantly associated with the 83 genes, suggesting that the steroid biosynthesis pathway activity may have an inverse correlation with antitumor immune response in GC.

Likewise, we found 59 and 57 genes having a strong positive and a significant inverse expression correlation with ICA in both cohorts, respectively, and identified 16 immune-related pathways associated with the 59 genes (https://cdn.amegroups.cn/static/application/4b0111ff02b66c45bc28328104e88994/10.21037atm-20-922-2.pdf). Furthermore, we found four genes (GBP1, GBP5, IDO1, and WARS) whose expression levels exhibited a strong positive correlation with PD-L1 expression levels (). The associations of these genes with PD-L1 and GC immunity have been investigated (32). In addition, we identified 80 genes whose expression levels had a significant negative correlation with PD-L1 expression levels (https://cdn.amegroups.cn/static/application/4b0111ff02b66c45bc28328104e88994/10.21037atm-20-922-2.pdf). GSEA revealed that the 80 genes were significantly associated with three KEGG pathways, including drug metabolism-cytochrome P450, peroxisome, and terpenoid backbone biosynthesis.

Collectively, the identification of genes with a significant expression correlation with antitumor immune response or protumor immunosuppression may provide useful intervention targets for promoting antitumor immunity and immunotherapy response in GC.

Identification of proteins whose expression levels are associated with antitumor immune response in GC

Based on TCGA-GC cohort, we identified six proteins (Caspase-7, PI3K-p85, PREX1, Lck, Bcl-2, and transglutaminase) whose expression levels exhibited a significant positive correlation with the enrichment levels of CD8+ T cells in GC (r>0.3) (). Notably, both Caspase-7 (33) and Bcl-2 (34) are involved in the apoptosis pathway which has been positively associated with antitumor immune signatures in cancer (35). Transglutaminase has been associated with tumor development and tumor immunity (36). The associations of PI3K-p85 (37), PREX1 (38), and Lck (39) with tumor immunity have been explored. In addition, we found that the expression levels of acetyl-CoA carboxylase (ACC) showed a significant inverse correlation with CD8+ T cell enrichment levels (r=−0.3). Previous studies have shown that ACC was essential for cancer cell survival and was an important target for cancer therapy (40,41). Our result suggests that ACC may promote cancer development by inhibiting antitumor immune response. Moreover, we found that ACC showed a significant inverse expression correlation with ICA in GC (r=−0.33), bolstering the notion that ACC suppresses antitumor immunity.

Figure 3

Correlation of protein expression levels with CD8+ T cell enrichment levels, immune cytolytic activity, and PD-L1 expression levels in TCGA-GC cohort. The Pearson correlation coefficients are shown.

Furthermore, we found six proteins (Caspase-7, PREX1, Lck, transglutaminase, Bax, and MEK1) whose expression levels had a significant positive correlation with PD-L1 expression levels in GC (r>0.3) (). Among them, the elevated expression of Caspase-7, PREX1, Lck, and transglutaminase was also significantly associated with highly enriched CD8+ T cells and ICA in GC. The positive correlation between the expression of these proteins and PD-L1 expression indicates that their upregulation could be a predictor for the active response to anti-PD-1/PD-L1 immunotherapy in GC since PD-L1 expression is a biomarker for predicting the anti-PD-1/PD-L1 immunotherapy response in cancer (42).

All together, the identification of proteins significantly correlating with antitumor immune signatures may provide predictive biomarkers or intervention targets for improving GC immunotherapy.

Identification of miRNAs whose expression levels are associated with antitumor immune response in GC

Based on TCGA-GC cohort, we identified 19 and 8 miRNAs whose expression levels had a significant positive and negative correlation with the enrichment levels of CD8+ T cells in GC, respectively (|r|≥0.3) (). 14 out of the 19 miRNAs also had a significant positive expression correlation with ICA in GC (r≥0.3) (https://cdn.amegroups.cn/static/application/a7d3d9d3c316a9cc2278b8b6881156f4/10.21037atm-20-922-3.pdf). A literature review showed that many of these miRNAs were associated with tumor immune response, such as hsa-miR-150, 155, 142, 342, 146, 101, 511, 29 (43). Notably, hsa-miR-155-5p expression levels were significantly associated with CD8+ T cell infiltration levels and ICA in GC (r>0.5). The positive association between hsa-miR-155 expression and tumor immunity has been confirmed in a recent study (44). The hsa-miR-150 genes hsa-miR-150-5p and hsa-miR-150-3p exhibited a significant positive expression correlation with CD8+ T cell infiltration levels and ICA in GC (r>0.5). Previous studies have shown that hsa-miR-150 played a key role in tumor pathogenesis and tumor immunity (45). hsa-miR-342 was importantly associated with tumor immunity (46). Our data showed that the expression levels of hsa-miR-342 had a significant positive correlation with CD8+ T cell enrichment levels and ICA in GC (r>0.5).

Figure 4

19 and 8 miRNAs whose expression levels show a significant positive and negative correlation with CD8+ T cell enrichment levels in TCGA-GC cohort, respectively (|r|≥0.3).

We found 15 miRNAs whose expression levels were significantly and positively associated with PD-L1 expression levels in GC (r≥0.3). The 15 miRNAs included hsa-let-7i-3p, hsa-let-7i-5p, hsa-miR-142-3p, hsa-miR-142-5p, hsa-miR-146b-3p, hsa-miR-146b-5p, hsa-miR-150-3p, hsa-miR-155-3p, hsa-miR-155-5p, hsa-miR-342-3p, hsa-miR-342-5p, hsa-miR-4772-3p, hsa-miR-4772-5p, hsa-miR-511-5p, and hsa-miR-7702 (https://cdn.amegroups.cn/static/application/a7d3d9d3c316a9cc2278b8b6881156f4/10.21037atm-20-922-3.pdf). Previous studies have shown that many of these miRNAs regulated PD-L1 expression on tumor cells (47).

Collectively, the identification of miRNAs correlating with antitumor immune signatures or tumor immunosuppressive signatures may provide potential miRNA biomarkers for GC immunotherapy.

Identification of LncRNAs whose expression levels are associated with antitumor immune response in GC

LncRNAs play important roles in immune regulation and tumor immunity (48-50). Based on TCGA-GC cohort, we identified 107 and 19 lncRNAs with a significant positive and negative expression correlation with the enrichment levels of CD8+ T cells in GC, respectively (|r|≥0.3) (https://cdn.amegroups.cn/static/application/98131eca64c3c7dbb54af3a6732c4f4b/10.21037atm-20-922-4.pdf). Notably, AC022706.1 had the strongest expression correlation with CD8+ T cell enrichment levels in GC (r=0.64) (). It also had a strong expression correlation with ICA in GC (r=0.58) (). LINC01871 had the second strongest expression correlation with CD8+ T cell enrichment levels (r=0.63) and the strongest expression correlation with ICA in GC (r=0.71) (). The expression levels of LINC01871 were also significantly correlated with PD-L1 expression levels in a positive direction in GC (r=0.59) (). AC006033.2 had the third strongest expression correlation with CD8+ T cell enrichment levels (r=0.60), the second strongest expression correlation with ICA in GC (r=0.59), and the third strongest expression correlation with PD-L1 expression in GC (r=0.54) ().

Figure 5

Three long non-coding RNAs (lncRNAs) AC022706.1, LINC01871, and AC006033.2 whose expression levels exhibit a significant positive correlation with CD8+ T cell enrichment levels, immune cytolytic activity, and PD-L1 expression levels in TCGA-GC cohort. Pearson’s correlation test P values and correlation coefficients are shown.

We identified 58 and 79 lncRNAs with a significant positive and negative expression correlation with ICA in GC, respectively (|r|≥0.3) (https://cdn.amegroups.cn/static/application/98131eca64c3c7dbb54af3a6732c4f4b/10.21037atm-20-922-4.pdf). Moreover, we found 58 and 17 lncRNAs whose expression levels were positively and negatively associated with PD-L1 expression levels in GC, respectively (|r|≥0.3) (https://cdn.amegroups.cn/static/application/98131eca64c3c7dbb54af3a6732c4f4b/10.21037atm-20-922-4.pdf).

Because the associations between lncRNAs and tumor immunity have not been sufficiently explored, the identification of lncRNAs associated with antitumor immune signatures may provide useful lncRNA candidates for furthering such exploration.

Identification of cancer-associated pathways whose activity is associated with antitumor immune response in GC

Based on TCGA-GC cohort, we identified nine cancer-associated pathways whose activity was positively correlated with CD8+ T cell infiltration levels in GC [Spearman correlation coefficient (ρ) ≥0.3] (). The nine pathways included focal adhesion, mTOR signaling, PI3K-AKT signaling, MAPK signaling, HIF-1 signaling, VEGF signaling, apoptosis, Jak-STAT signaling, and cytokine-cytokine receptor interaction. The positive association between a majority of these pathways (such as focal adhesion, MAPK signaling, VEGF signaling, apoptosis, and Jak-STAT signaling) and GC immunity has been revealed in our previous study (21). The HIF-1 signaling pathway is hyperactivated in various cancers and is associated with the activation of tumor‐associated immune signatures (51). This is consistent with our result that there was a significant positive correlation between this pathway and CD8+ T cell infiltration in GC. Moreover, our previous study (35) showed that tumor glycolysis had a significant positive association with immune signatures in a wide variety of cancers. This is consistent with that the glycolysis-stimulating HIF-1 signaling was positively associated with immune signatures in GC. In addition, our results showed that the mTOR signaling pathway was associated with CD8+ T cell infiltration levels in GC. This is in line with previous studies showing that mTOR signaling played a significant role in regulating tumor microenvironment, thereby affecting tumor immunity (52). We further identified a number of cancer-associated pathways which had a positive correlation with ICA and PD-L1 expression levels in TCGA- and ACRG-GC cohorts, including p53 signaling, HIF-1 signaling, apoptosis, Jak-STAT signaling, and cytokine-cytokine receptor interaction (https://cdn.amegroups.cn/static/application/a41ae079640a0bc9bd35cf88aec82e8d/10.21037atm-20-922-5.pdf).

Figure 6

Nine cancer-associated pathways whose activity is positively associated with CD8+ T cell infiltration levels in TCGA-GC cohort [Spearman correlation coefficient (ρ) ≥0.3].

Collectively, the identification of antitumor immunity-associated pathways provided potential intervention targets for promoting antitumor immunity and immunotherapy response in GC.

Discussion

Our bioinformatics analysis identified a number of genes, miRNAs, lncRNAs, proteins, and pathways whose mutations, expression, or activity were significantly associated with antitumor immunity in GC. Of note, we found some important oncogenes and tumor suppressor genes whose mutations were appreciably associated with antitumor immunity in GC, including PIK3CA, MTOR, RNF213, TP53, ARID1A, PTEN, ATM, and CDH1. These significant associations may partially explain why the mutations of some of these genes were associated with OS in GC. Interestingly, we found seven genes (ARID1A, BCOR, MTOR, CREBBP, SPEN, NOTCH4, and TET1) whose mutations were associated with better OS in GI cancer patients receiving anti-PD-1/PD-L1 immunotherapy but were not associated with OS in GC patients not receiving such a therapy (). A potential mechanism underlying this difference is the elevated expression of PD-L1 in the cancer patients with the mutation of these genes. To date, many studies have revealed that PD-L1 expression is a predictive biomarker for anti-PD-1/PD-L1 immunotherapy in GC patients (53,54). A clinical trial study (KEYNOTE-059) revealed that pembrolizumab, an anti-PD-1 monoclonal antibody, had encouraging antitumor activity in PD-L1-positive advanced gastric/gastroesophageal cancer patients (53). Jin et al. reported that a PD-L1-positive unresectable locally advanced GC with MSI exhibited pathological complete response to a single dose of anti-PD-1 immunotherapy in combination with chemotherapy (54). Hence, the identification of gene mutations associated with tumor immunity may furnish biomarkers for the optimal stratification of GC patients responsive to immunotherapy. Likewise, we identified some important genes and proteins whose expression levels were significantly associated with antitumor immunity in GC and thus may serve as biomarkers or intervention targets for promoting GC immunotherapy. For example, ACC exhibited a significant inverse expression correlation with CD8+ T cell infiltration levels and ICA in GC, suggesting that inhibiting ACC could enhance antitumor immunity in GC. Furthermore, we identified a number of miRNAs and lncRNAs having a significant expression correlation with GC immunity, such as the miRNAs hsa-miR-150, 155, 142, 342, 146, 101, 511, 29 and the lncRNAs AC022706.1, LINC01871, and AC006033.2. These non-coding RNAs could be potentially useful biomarkers for GC prognosis and immunotherapy. Finally, we identified numerous cancer-associated pathways whose activity was significantly associated with antitumor immunity in GC, including the mTOR, PI3K-AKT, MAPK, HIF-1, and VEGF signaling pathways. These pathways were mainly involved in cell proliferation, growth and migration, metabolism, and immune response.

TMB often has a positive association with tumor immunity and immunotherapy response (25,55) while TALs often inversely correlate with them (56). As expected, we observed a significant positive correlation of TMB with PD-L1 expression levels in TCGA-GC cohort (Student’s t test, P=0.03) () and a significant negative correlation of TALs with CD8+ T cell enrichment levels, ICA, and PD-L1 expression levels in TCGA- and ACRG-GC cohorts (Student’s t-test, P<0.001) (). This indicates that TMB and TALs are useful biomarkers for predicting GC patients responsive to immunotherapy.

Figure 7

Associations of tumor mutation burden (TMB) and tumor aneuploidy levels (TALs) with CD8+ T cell enrichment levels, immune cytolytic activity (ICA), and PD-L1 expression levels. (A) TMB has a significant positive correlation with PD-L1 expression levels in TCGA-GC cohort (Student’s t test, P=0.03). (B) TALs have a significant negative correlation with CD8+ T cell enrichment levels, immune cytolytic activity, and PD-L1 expression levels in both TCGA- and ACRG-GC cohorts (Student’s t test, P<0.001). For each tumor sample, its TMB is the total somatic mutation count and TAL is the tumor ploidy score calculated by ABSOLUTE (29). ***, P<0.001.

In conclusion, the molecular features significantly associated with GC immunity could be useful biomarkers for stratifying GC patients responsive to immunotherapy or potential intervention targets for promoting antitumor immunity and immunotherapy response in GC.

The article’s supplementary files as