Immune signature profiling identified prognostic factors for gastric cancer.

OBJECTIVE: Tumor microenvironment, especially the host immune system, plays a pivotal role in tumor initiation and progression. Profiling of immune signature within tumor might uncover biomarkers for targeted therapies and clinical outcomes. However, systematic analysis of immune-related genes in gastric cancer (GC) has not been reported.
METHODS: Expressions of a total of 718 immune-related genes were generated in 372 stomach adenocarcinoma (STAD) patients from The Cancer Genome Atlas (TCGA) database using RNA-sequencing data. Integrated bioinformatics analyses were performed to identify prognostic factors as well.
RESULTS: Survival analyses revealed 73 genes, which were significantly associated with patient's overall survival (OS). Taken together with clinicopathological parameters, we established a predictive model, containing 10 immune-related genes, which were NRP1, C6, CXCR4, LBP, PNMA1, TLR5, ITGA6, MICB, PBK and TNFRSF18, with powerful efficiency in distinguishing satisfactory or poor survival of STAD patients. Moreover, the top 3 ranked prognostic genes, NRP1, TGFβ2 and MFGE8, were also significantly associated with patient's OS by an independent validation achieved from Kaplan-Meier plotter database.
CONCLUSIONS: We profiled prognostic immune signature and established prognostic predictive model for GC, which could reflect immune disorders within tumor microenvironment, and also may provide novel predictive and therapeutic targets for GC patients in the near future.

Introduction

Gastric cancer (GC) is the fourth most frequently diagnosed malignancy and the second leading cause of cancer death worldwide (1,2). Although the incidence of GC has declined for decades, the prognosis of GC remains very poor, especially in China (3). At present, the pathogenesis of GC is unclear, thereby necessitating effective biomarkers and targeted therapeutics. Traditionally, clinicopathological parameters were used in risk stratification of GC outcomes. However, a number of advanced GC patients remained stable for a couple of years, whereas some early-stage patients progressed rapidly (4). Therefore, reliable biomarkers or stratification systems that can be used for more accurate prediction are highly essential.

The critical role of the immune microenvironment in tumor biology has been increasingly appreciated in recent years. Evading immune surveillance is one of the hallmarks of cancer (5,6) and recently, immunotherapeutic agents targeting immunosuppressive proteins [e.g., cytotoxic T lymphocyte-associated antigen-4 (CTLA-4), programmed cell death 1 (PD-1), and programmed cell death 1 ligand 1 (PD-L1)] have shown promising clinical efficacy towards melanoma as well as several other types of cancer (7-10). Considering the poor outcomes after standard treatment and few targeted therapeutics in GC, immunotherapy, as a promising additional approach, is currently under intensive investigation. The most prospective response rates obtained by this class of immunotherapeutic drugs were induced by anti-PD-1/PD-L1 antibodies, which showed efficacy and tolerability in GC patients. The phase 3 trial of nivolumab (anti-PD-1 checkpoint inhibitor) versus placebo in Asian patients with ≥2 prior lines of treatment confirmed that nivolumab can prolong overall survival (OS) (11), leading to the approval of nivolumab in Japan. The results achieved using pembrolizumab (anti-PD-1) administered as 3 lines or later treatment for patients with advanced PD-L1+ GC in KEYNOTE-059 led to accelerated approval in the United States (12). The accumulated data indicated that anti-CTLA-4 antibodies yielded partially satisfactory results; besides, patients with a post-treatment carcinoembryonic antigen (CEA) antigen proliferative response had a median survival time of 17.1 months compared with 4.7 months for non-responders to the anti-CTLA-4, tremelimumab, suggesting a rationale for combination of CTLA-4 blockade with vaccines, targeting GC antigens in the future (13). In addition to CTLA-4 and PD-1/PD-L1, other suppressive immune checkpoint proteins e.g., LAG-3 (a lymphocyte surface protein) and IDO (an enzyme that catabolizes tryptophan, creating an immunosuppressive tumor microenvironment) and inhibitors of these proteins may act in a synergic-manner with anti-PD-1/PD-L1 antibodies to produce a more robust antitumor immune response (14). Meanwhile, tumor-infiltrating lymphocytes have been reported to be useful in predicting the prognosis of GC patients (15-17), suggesting that distinct immune status has a profound influence on outcome of GC patients.

Molecular profiles of tumor cells and cancer-related cells within their microenvironments represent promising candidates for predictive and prognostic biomarkers (9,18-21). Immune signature profiling identified predictive and prognostic factors for esophageal squamous cell carcinoma, pancreatic ductal adenocarcinoma, ovarian cancer, and non-small cell lung cancer (NSCLC) (22-25). For instance, in esophageal squamous cell carcinoma, 9 immune-related genes (ABL1, ATF2, ATG5, C6, CD38, HMGB1, ICOSLG, IL12RB2 and PLAU) were identified associating with patients’ OS, of which a prognostic model [including three independent factors (ABL1, CD38, and ICOSLG)] was established. Understanding of interactions between tumor and the host immune system therefore holds a great promise to uncover biomarkers for targeted therapies and clinical outcomes (22). However, systematic analysis of immune signatures in GC patients should be conducted.

Therefore, systematically investigation of immune signature within GC microenvironment is of great importance to provide predictive and prognostic biomarkers, and also guide effective immunotherapies in GC.

Materials and methods

Data acquisition

RNA-seq data of stomach adenocarcinoma (STAD) cohorts were downloaded from The Cancer Genome Atlas (TCGA) database (https://tcga-data.nci.nih.gov/tcga/). A total of 718 immunology-related human genes, curated from nCounter® PanCancer Immune Profiling Panel (NanoString Technologies Inc., Seattle, WA, USA), were then implemented as candidate genes in this studies. For validation, independent survival analysis was performed with the help of Kaplan-Meier-plotter (http://kmplot.com/). A total of 876 GC patients were included from five Gene Expression Omnibus (GEO) datasets (GSE14210, GSE15459, GSE22377, GSE29272, and GSE51105).

Bioinformatics analysis

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment for the dysregulated genes were analyzed by GeneMANIA (26). Univariate and multivariate Cox regression analyses were performed to search for survival-associated genes using survival package (Version. 2.41). In the multivariate Cox regression, 10 immune-related genes, which were significant in the univariate Cox regression, together with age, T-stage, and M-stage were potential independent factors in the regression. To compare the ability of the prognostic predictors, survival receiver operating characteristic (ROC) package (Version. 1.0.3) in R, allowing time-dependent illustration of ROC curve with censored data (27), was used to calculate the value of the area under the curve (AUC) for each parameter.

Statistical analysis

The survival curves were compared using Kaplan-Meier method and log-rank test. All tests were two-sided, and a P-value of less than 0.05 was considered statistically significant. Data were analyzed using R software (Version 3.3.1; R Foundation for Statistical Computing, Vienna, Austria).

Results

Prognostic immune signatures in STAD

Before studying the prognostic values of immune signatures, univariate survival tests were conducted to assess the relationship between clinical parameters and outcomes in this STAD cohort. As shown in , tumor, node, metastasis (TNM) stages were significantly associated with OS (P<0.05). Meanwhile, age was also significantly associated with OS (P=0.045,). The results of this preliminary assessment indicated that the survival data for the TCGA-STAD cohort were informative and appropriate to be used in the further analysis.

Next, we separated the cohort into high expression and low expression groups based on expression of the 9 genes, respectively. Survival analysis was then conducted using survival package (Version 2.41.3) in R (Version 3.3.1). Consequently, a total of 73 genes were statistically significant in the univariate analysis (). Among these genes, 22 genes were favorable prognostic factors, while 51 genes were adverse prognostic factors (). The top 3 ranked genes were NRP1, TGFβ2 and MFGE8, and Kaplan-Meier plots of these three genes are shown in −. GO and KEGG analysis highlighted that these genes were significantly enriched in important immune responses, including positive regulation of defense responses, leukocyte migration, activation of innate immune responses, pattern-recognition receptor signaling pathway, toll-like receptor signaling pathway, inflammatory response, and positive regulation of mitogen-activated protein kinase (MAPK) cascade ().

Prognostic immune signatures in stomach adenocarcinoma (STAD) and enriched Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG). (A) Hazard ratio and P values of 73 immune related genes which were significantly associated with patients’ overall survival; (B) Significantly enriched GO and KEGG which were analyzed using GeneMANIA.

Kaplan-Meier plots and receiver operating characteristic (ROC) curves of survival associated genes and prognostic predictor. (A−L) Kaplan-Meier plots of survival associated genes. Patients were divided into high expression (red line) and low expression (blue line) based on their gene expression by median cut; (M) Kaplan-Meier plot of prognostic predictor built with 10 independent factors (NRP1, C6, CXCR4, LBP, PNMA1, TLR5, ITGA6, MICB, PBK and TNFRSF18). Patients were divided into high risk (red line) and low risk (blue line) by the predictor; (N, O) ROC curves of each parameters with area under the curve (AUC) scores.

Prognostic predictor was established with survival-associated immune genes

To move any genes that might not be independent factors for STAD patients, multivariate Cox regression analysis with an established prognostic model was applied to these candidate immune-related genes together with clinicopathological parameters. As shown in , 10 immune-related genes together with age, T-stage, and M-stage were potential independent factors. Consistently, Kaplan-Meier plots revealed that all these 10 genes were notably associated with patients’ outcomes (). Among these 10 genes, higher expression of NRP1, C6, CXCR4, LBP, PNMA1 and TLR5 indicated poor OS, while ITGA6, MICB, PBK and TNFRSF18 were favorable prognostic genes in STAD (,). Therefore, we established the final prognostic predictors with the help of these 10 genes, age, T-stage, and M-stage, and ROC curves were also applied to compare the efficiency of these predictive models and genes (−). The AUC of the ROC was 0.853, which was remarkably higher than other individual genes (,), indicating that the final predictor indeed showed a significantly superior performance in distinguishing satisfactory or poor survival in STAD patients ().

Validation of top 3 ranked prognostic factors

To validate survival-associated immune signatures in TCGA-STAD cohort, the top 3 ranked genes, which were NRP1, TGFβ2 and MFGE8, were chosen as candidates. Kaplan-Meier plotter established a robust database, realizing swift validation of previous and future gastric cancer survival biomarker candidates (28). Therefore, we chose this database, consisting of 876 GC patients, to perform an independent validation study. Consistently, NRP1, TGFβ2, and MFGE8 were significantly associated with patient’s OS (). All these three genes were adverse prognostic genes, in which higher expression indicated worse OS of GC patients ().

Independent validation of top 3-ranked prognostic genes by Kaplan-Meier plotter. (A−C) Kaplan-Meier plots of top 3-ranked prognostic genes,NRP1, TGFβ2, and MFGE8, in an independent cohort of 876 gastric cancer (GC) patients. HR, hazard ratio.

Taken together, our results demonstrated that our analyses were largely consistent across databases, suggesting that our predictive model was reliable, with a great potential in distinguishing satisfactory or poor outcome of STAD patients.

Discussion

Growing evidence has suggested that tumor microenvironment plays a pivotal role in tumor initiation and progression (29). Dysfunction of the immune system in cancer patients allows cancer cells to avoid immune surveillance. The occurrence, development, and prognosis of GC are closely associated with cross-talks between different immune cells and GC cells (30). Immunotherapeutic agents have emerged as a novel approach for cancer immunotherapy in several malignancies. However, patient prognosis and disease progression relevant to immune-related genes in GC microenvironment have not been elucidated. High throughput technology provides objective data for systematic profiling of immune signatures in clinical samples (31). In the present study, we curated immune signatures, including 718 immune-related genes from TCGA RNA-sequencing data from 372 STAD patients. Survival analysis revealed 10 independent prognostic factors, which were NRP1, C6, CXCR4, LBP, PNMA1, TLR5, ITGA6, MICB, PBK, and TNFRSF18, as well as a predictive model established with high performance in distinguishing satisfactory or poor outcomes for STAD patients. Moreover, the top 3 ranked prognostic genes, NRP1, TGFβ2 and MFGE8, were further validated in independent data from Kaplan-Meier plotter database, in which the results were largely consistent with our results in an independent cohort of 876 GC patients. These findings suggest that our findings might provide reliable novel insights in characterization of tumor microenvironment of STAD.

More importantly, we provided a list of novel candidates in targeted therapy for GC patients. Among these survival-associated genes, some have been reported with preliminary result. Additionally, NRP1, a top ranked prognostic factor in our results, was correlated with clinical staging, tumor differentiation, and pathological types of GC, and was also involved in the growth and metastasis of GC cells (32). Overexpression of PBK in GC plays a crucial role in tumor malignant potential and demonstrates to be a prognostic factor and a potential therapeutic target. Genetic variants in TLR5 may modify the role of Helicobacter pylori infection in the process of causing GC.

Interactions between GC cells and immune cells were also highlighted in our analyses. For instance, macrophages within tumor microenvironment might promote invasion of GC cells though enhancing transforming growth factor β (TGF-β)/bone morphogenic proteins (BMPs) signaling pathway, playing a pivotal role in regulating GC progression (33). TGF-β is correlated with poorer prognosis in GC. Direct interaction between GC cells and peripheral blood mononuclear cells may promote production of TGF-β1, that might facilitate cancer development (34). Moreover, MFGE8, an anti-fibrotic protein which is secreted by mesenchymal stem cells, could strongly inhibit TGF-β signaling pathway and reduce deposition of extracellular matrix (35). Moreover, primary gastric adenocarcinoma tissue could release soluble MICB into the extracellular milieu and induce a significant decrease in the levels of NKG2D receptor on effector natural killer cells and CD8+ T cells, correlating with an impaired cytotoxic function (36,37). Furthermore, CXCL12/CXCR4 derived from cancer-associated fibroblasts (CAFs) can promote invasion of GC cells by enhancing the clustering of integrin β1 in GC cells, resulting in GC progression (38), and cross-talk between CXCR4 and CXCR2 might contribute to epithelial-mesenchymal transition, migration and invasion of GC (39). Taken together, these results suggested that immune cells and other mesenchymal cells were dynamically involved in GC development, and the survival-associated genes might serve as new biomarkers and therapeutic targets.

Conclusions

In this study, we established our prognostic predictor by profiling of immune signatures in GC, which could be regarded as immune-related protective and risky patterns in GC microenvironment. Our study is broadening current observations by exploring multicomponent panel of immune-related genes to ascertain their prognostic effects in GC. These significantly prognostic genes provide alternative targets alone or may be better on combination, attributing to their immune nature and prognostic significance, thereby improving the clinical outcomes.

Acknowledgements

This study was financially supported by the Scientific Research Foundation of Shanxi Province Healthy Commission (No. 2017068) the Doctor Scientific Research Foundation of Shanxi cancer hospital (No. 2017A03), and the Natural Science Foundation of Shanxi Province (No. 201801D221259).

Footnote

Conflicts of Interest: The authors have no conflicts of interest to declare.

1

Univariate and multivariate Cox regression analyses results in TCGA STAD cohort

Variables	Univariate analyses			Multivariate analyses
	HR (95% CI)	P		HR (95% CI)	P
TCGA, The Cancer Genome Atlas; STAD, stomach adenocarcinoma; HR, hazard ratio; 95% CI, 95% confidence interval.
Age	1.387 (1.008−1.910)	0.045		3.079 (1.963−4.829)	<0.001
Sex	0.795 (0.563−1.122)	0.186		0.820 (0.528−1.276)	0.380
T stage	1.713 (1.151−2.548)	0.005		2.753 (1.240−4.091)	0.008
N stage	2.086 (1.384−3.145)	<0.001		1.465 (0.743−2.889)	0.271
M stage	2.040 (1.260−3.303)	0.008		2.151 (1.089−4.248)	0.027
TNM stage	1.823 (1.307−2.541)	<0.001		1.735 (0.961−3.133)	0.068
C6	1.552 (1.123−2.146)	0.008		2.004 (1.240−3.240)	0.005
CXCR4	1.648 (1.189−2.283)	0.003		2.675 (1.465−4.882)	0.001
ITGA6	0.635 (0.460−0.877)	0.006		0.513 (0.310−0.849)	0.009
LBP	1.754 (1.267−2.427)	<0.001		1.691 (1.072−2.669)	0.024
MICB	0.670 (0.486−0.924)	0.015		0.432 (0.257−0.727)	0.002
NRP1	2.025 (1.454−2.820)	<0.001		3.212 (1.700−6.071)	<0.001
PBK	0.710 (0.515−0.978)	0.036		1.795 (1.039−3.099)	0.036
PNMA1	1.612 (1.166−2.229)	0.004		1.819 (1.127−2.935)	0.014
TLR5	1.637 (1.183−2.267)	0.003		1.926 (1.126−3.296)	0.017
TNFRSF18	0.705 (0.512−0.917)	0.032		0.575 (0.345−0.959)	0.034