This study was approved by the local Research Ethics Board of The Wenzhou Medical University Second Affiliated Hospital, and informed patient consent was waived. All procedures performed in studies involving human participants were in accordance with the 1964 Helsinki declaration and its later amendments.

Gene sequencing data and the corresponding clinical data for stomach adenocarcinoma (STAD) and adjacent normal tissue samples were downloaded from The Cancer Genome Atlas database (TCGA; https://portal.gdc.cancer.gov). In addition, the raw GSE62254/ACRG (Asian Cancer Research Group) data, containing 300 AGC patients, including survival information, were downloaded from the Gene Expression Omnibus database (GEO: https://www.ncbi.nlm.nih.gov/geo/). Z-score normalization was performed for all gene sequencing data.

Differentially expressed genes (DEGs) between tumor and adjacent normal tissues were determined with the filtering condition of log2 |fold change | > 1 and a false discovery rate (FDR) < 0.05, using the limma package in R software 3.6.3. The IRG list, which has been validated to be involved in immunity, was downloaded from ImmPort, a platform that provides accurate and timely immunological data. Differentially expressed IRGs were thus extracted from the DEGs.

Patients who have an overall survival (OS) of <90 days as well as missing survival information were excluded from further study. Thus, a total of 364 patients were finally enrolled. Univariate Cox and further LASSO regression analysis was used to establish the prognostic model. The survival receiver operating characteristic (ROC) curve was performed to evaluate the performance of the ASIRGPI and to determine the cutoff for classifying AGC patients as low- or high-risk. Survival analysis associated with the prognostic model was carried out via Kaplan–Meier analysis.

Furthermore, multivariate Cox regression analyses were used to evaluate whether the prognostic model could independently predict the prognosis of AGCs. A nomogram was thus formulated using the coefficients of the multivariable Cox regression model via the rms package in R. Calibration curves were assessed graphically by plotting the observed rates against the nomogram-predicted probabilities.

ESTIMATE, which uses gene expression signatures to infer stromal and immune cell fractions to determine stromal and immune scores, was performed to analyze the TME subcomponent, while xCELL, which also uses gene expression to infer the proportions of 64 tumor-infiltrating immune cell (TIIC) and stromal cell types, was further used to analyze the TME cell type. Additionally, gene set variation analysis (GSVA) was performed to estimate the enrichment scores of 110 immunoregulation-related pathways in AGC while the immunomodulatory gene was estimated via analyzing 78 immunomodulatory genes summarized by the TCGA immune response working group. A detailed flowchart is shown in Figure 1.

Tissue microarrays from 154 local AGC patients (35 tissue pairs and 119 GC tissue) were constructed according to our previous research (18). IHC was performed in the Pathology Laboratory of The Second Affiliated Hospital of Wenzhou Medical University, using rabbit anti-TGFB1 monoclonal antibody (Abcam, cat. no. ab27969) at a dilution of 1:250. The results were initially determined by two pathology experts and were accepted if a third expert also confirmed the result. Otherwise, the data were reviewed by all three experts and discussed until a consensus was reached. Finally, a computer-automated method was conducted (Image-Pro Plus 6.0; Media Cybernetics, Inc.) and the expression level was represented as numbers of positive cells per square millimeter (positive cells number/total area).

All statistical analyses are carried out using R (version 3.6.3) and SPSS 22.0. Pearson correlation analysis was used to determine the correlation, and survival analysis was performed using the log-rank test. p < 0.05 was considered statistically different.

To identify Asian-specific IRGs, we first selected all Asian tissues (74 tumor and 7 adjacent normal tissues) in the TCGA-STAD cohort and the differential expression patterns of 2,498 IRGs in AGC tissues were further analyzed to form an expression profile to analyze differentially expressed IRGs. We determined 685 DEGs, including 515 upregulated and 130 downregulated genes between GC and adjacent normal samples, using the limma package with cutoffs |log fold change| > 1 and FDR < 0.05. After removing genes not detected in the ACRG dataset, 183 upregulated and 67 downregulated IRGs were finally identified (Supplementary Figures S1A, B). The expression profiles of the 250 IRGs in the TCGA-STAD Asian cohort and TCGA-STAD white cohort are shown in Figure 2A. Further, gene set enrichment analysis (GSEA) revealed that TGFβ signaling, coagulation, myogenesis, TNFA signaling, and epithelial–mesenchymal transition (EMT) were enriched in the TCGA-STAD white cohort, while MYC targets V2, DNA repair, G2M checkpoint, MYC targets V1, and E2F targets were enriched in the TCGA-STAD Asian cohort (Figure 2B and Supplementary Figure S2).

To further explore the survival differences between the different ethnicities, a Kaplan–Meier survival curve was performed. After excluding patients with an OS of <90 days and those without survival information, we found that white GC patients (TCGA-STAD white cohort) correlated with poor prognosis compared with AGC patients (ACRG cohort and TCGA-STAD Asian cohort; p = 0.004, Figure 2C).

Meanwhile, statistical differences of infiltrating TME cells between the Asian and white GC patients were investigated; we found that the fractions of B cells naïve, plasma cells, T cells CD8, T cells CD4 memory resting, T cells regulatory, NK cells resting, NK cells activated, monocytes, macrophages M0, and macrophages M2 were significantly increased in AGC, while T cells follicular helper, T cells gamma delta, macrophages M1, dendritic cells resting, dendritic cells activated, mast cells activated, eosinophils, and neutrophils were significantly decreased (Figure 2D).

To further explore the clinical significance of differentially expressed IRGs, univariate Cox regression analysis was performed and 35 IRGs were found to be survival-associated using the 64 AGC samples and clinical information for the TCGA-STAD cohort (Figure 3A). LASSO regression analysis further reduced 35 potential IRGs to 7, which had non-zero coefficients in the regression model (Figures 3B, C). An Asian-specific IRG-based prognostic index was calculated for each patient using the following formula: ASIRGPI = (0.0146 * expression of CRABP1) + (0.2543 * expression of F2R) + (0.0413 * expression of LTB) + (0.1680 * expression of PLSCR1) + (0.2501 * expression of S100B) + (0.0979 * expression of SEMG1) + (0.1692 * expression of TYROBP).

ROC analysis was carried out to assess the ASIRGPI. The high area under the curve (AUC of 0.903) confirmed the high prognostic performance of the ASIRGPI in survival surveillance. Additionally, the maximal Youden index value for an ROC curve was determined as the optimal cutoff point (0.3086; Supplementary Figure S3). AGC patients in the TCGA were thus assigned into low- and high-risk cohorts. High-risk patients (16, 25.0%) had shorter OS (HR = 56.89, p < 0.001) than low-risk patients (48, 75.0%) among the 64 Asian patients with GC (Figure 3E). The expression of the seven prognostic genes and the relationship between ASIRGPI distribution and the survival status are shown in Figures 3D, E.

To further validate our ASIRGPI, its performance was assessed in the ACRG cohort, which consisted of 300 AGC patients. Consistently, high-risk patients (92, 30.7%) had a worse prognosis than low-risk patients (208, 69.3%) in the ACRG validation cohort (HR =1.42, p = 0.037; Figures 4A, B). Interestingly, for the white GCs of the TCGA-STAD cohort, ASIRGPI categorized 143 (75.7%) and 46 (24.3%) patients into the low-risk and high-risk groups, respectively. However, no OS difference was found (HR =1.14, p = 0.630; Figures 4C, D) between the two groups.

The clinical characteristics of the patients in the TCGA and ACRG cohorts are depicted in Table 1. To explore the prognostic value of ASIRGPI, multivariate Cox regression analyses were conducted for all 364 AGC patients. We found that the ASIRGPI, age, and TNM stage serve as independent predictors of AGC patient survival outcomes (Figure 5A). Similarly, high-risk patients (88, 29.7%) had shorter OS (HR = 1.79, p < 0.001) than low-risk patients (256, 70.3%; Figure 5B).

To further elaborate the clinical significance of the ASIRGPI, the relationship between the ASIRGPI and clinical and demographic characteristics, including age and TNM stage (T stage, lymphatic invasion, and distant metastasis), according to the International Union against Cancer was analyzed. We found that the ASIRGPI was positively correlated with both age and TNM stage (Figures 5C, D).

Finally, to provide a more intuitive clinical application tool, a nomogram was integrated with the ASIRGPI; age and TNM stage were constructed based on multivariate Cox analysis (Figure 5E). Calibration plots showed that the nomogram could predict the probability of 3- and 5-year OS well (Figure 5F).

To explore potential ASIRGPI-related survival mechanisms, ESTIMATE was performed and stromal/immune scores were inferred via ssGSEA for the entire cohort. Pearson’s correlation analysis showed that both the stromal and immune scores were positively correlated with ASIRGPI [r = 0.74, p < 0.001 (Figure 6A) and r = 0.70, p < 0.001 (Figure 6B), respectively]. Additionally, using the median score as the cutoff values, survival analysis showed that high-stromal-scored AGC patients had a worse OS than low-stromal-scored patients (HR = 1.59, p = 0.003; Figure 6A), but no OS difference was observed between low- and high-immune-scored patients (HR = 1.16, p = 0.330; Figure 6B). Furthermore, TME cell-type analysis was performed and 64 TME cell types were inferred. Among these, nine cell types (Figure 6C) were both significantly related to overall survival (log-rank test, p < 0.05) and the ASIRGPI (|r| ≥ 0.40, p < 0.05). Remarkably, four stromal cell types with the largest proportions as well as five other cell types, including lymphoid, hematopoietic stem, and epithelial cells, were all positively associated with poor prognosis and ASIRGPI.

Next, to explore potential ASIRGPI-associated therapeutic targets, GSVA was performed to estimate the enrichment scores of the 110 immunoregulation related pathways. We found 13 panimmune gene sets to be significantly associated with overall survival (P < 0.05) and with ASIRGPI (|r| ≥ 0.40, P < 0.05); 11 were positively associated with worse outcomes and ASIRGPI risk scores (Figure 6D). The remaining two genes were positively associated with worse outcomes and negatively associated with ASIRGPI risk scores, while the Module11_Prolif_score was the opposite. Interestingly, TGFβ-related gene sets (TGFB_PCA_17349583 and TGFB_score_21050467), which are closely associated with the remodeling of stromal components in the TME, were highlighted during this screening.

To further analyze ASIRGPI-related molecular targets, 29 immunomodulatory genes that could be detected in the entire cohort of 364 patients were explored ; only two (EDNRB and TGFB1) were found to be positively correlated with the ASIRGPI (r ≥ 0.40, P < 0.05) and significantly associated with poor outcomes (Figure 7A). OS analysis (Figures 7B, C) further proved that expression of the therapeutic target TGFβ1 was only significantly associated with poor outcome of the AGC patients (TCGA-STAD Asian cohort and ACRG cohorts), while no such phenomenon was found in the TCGA-STAD White cohort (Figure 7D).

To verify the effect of TGFB1 on Asian gastric cancer, IHC detection of TGFB1 was performed in the TMA. The TGFB1 expression in 35 pairs of matched tissues was analyzed (Supplementary Figure S4A). The TGFB1 expression was significantly higher in the tumor tissue (p = 0.031, Supplementary Figure S4B). Among all samples, 83.8% (129 of 154) were TGFB1-positive overall. A cutoff point of 8.68, which was optimized using a ROC, was chosen to categorize patients into TGFB1 high and low expressed subgroups. Further survival analysis found that the high-TGFB1-expressed patients had a shorter OS than the low-expressed patients (p < 0.001; Figures 7E, F).