A retrospective analysis was conducted on 7,192 patients who underwent gastrectomy at Zhejiang Cancer Hospital from January 2010 to December 2019. Among them, 5,567 patients met the following inclusion criteria: 1.Preoperative pathological biopsy confirmed primary gastric cancer; 2.Underwent radical or palliative gastrectomy; 3.No concomitant severe diseases such as acute cardiovascular and cerebrovascular diseases, liver cirrhosis, and chronic renal failure. Exclusion criteria: 1.Received neoadjuvant treatments such as preoperative radiotherapy, chemotherapy, or immunotherapy; 2.Number of dissected lymph nodes<16; 3.Presence of other heterogeneous tumors; 4.Other types of gastric cancer (e.g.,neuroendocrine carcinoma, squamous cell carcinoma, adenosquamous carcinoma); 5.Patients with missing critical clinical data. The median follow-up time was 85 months (interquartile range:71 months). All eligible patients underwent radical gastrectomy according to the Japanese gastric cancer treatment guidelines (Association JGC, 2020). Surgical methods included proximal, total, and distal gastrectomy. Postoperatively, specimens were reviewed by pathology experts at the Cancer Hospital of the Chinese Academy of Sciences. Pathological tumor-lymph node metastasis (pTNM) staging was based on the 8th edition of the American Joint Committee on Cancer (AJCC) TNM staging system (Amin et al., 2017). Potential curative resection was defined as R0 resection. Survival time was calculated from the date of surgery to the date of GC-related death or the most recent follow-up. The follow-up cut-off date was August 1, 2023. Perioperative management followed routine procedures, with no differences between groups. Patients meeting the above criteria were divided into two groups according to the Preoperative BMI standards set by WHO: the low BMI group (BMI<18.5 kg/m²) and the non-low BMI group (BMI≥18.5 kg/m²). Various clinicopathological characteristics, surgery-related indicators, and postoperative outcome factors were collected for analysis, including gender, height, and weight. BMI was calculated based on the patients’ height and weight. Tumor location was classified according to the center of the lesion as Upper 1/3 (cardia, fundus), Middle 1/3 (body), Lower 1/3 (antrum, including the angular incisure and pylorus), or involving the entire stomach (Total) (tumor involving more than 2/3 of the stomach wall). Tumor size was determined by the maximum diameter of the tumor. The positive levels of tumor markers were defined as CA199 ≥ 37 U/ml and CEA ≥ 5 ng/ml.

Samples were collected from 335 patients between January 2013 and December 2018 from Zhejiang Cancer Hospital. After screening according to the above clinical cohort criteria and ensuring that no antibiotics or intestinal microecological agents had been used in the previous month, 198 eligible GC patients were included in the final analysis. All patients were followed up by telephone and outpatient clinics with a follow-up cut-off date of 1 August 2023.The study was conducted by the Zhejiang Cancer Hospital (ZCH). The study was approved by the Ethics Committee of Zhejiang Cancer Hospital (approval number: IRB-2023-791) and written informed consent was obtained from all participants. Gastric samples were collected from patients who underwent gastrectomy, with peritumoral tissue 2-5 cm from the tumor margin. Notably, for metabolomics analysis, tissue specimens were subjected to cold ischemia for less than 30 minutes before freezing at -80 degrees Celsius. For 16S rRNA sequencing and transcriptome analyses, tissue specimens were immersed in an RNA-protecting solution at 4°C overnight, and then frozen at -80°C. Specimens for each histology were collected simultaneously. All tissue samples were collected at the time of surgical specimen removal. Histological sections at the top and bottom of each specimen were reviewed by a senior board-certified pathologist to confirm whether the tissue was tumor tissue or adjacent non-tumor tissue. For the purposes of this study, tumor samples had to have an average of 60% tumor cell nuclei and less than 20% necrosis to qualify.

Microbial DNA was extracted using an E.Z.N.A. Tissue DNA Kit (D3396-01; Omega, Norcross, Georgia, USA) following the manufacturer’s instructions as described previously. The DNAs were quantified using a Qubit 2.0 Fluorometer (Invitrogen, Carlsbad, CA, USA), and molecular size was estimated using agarose gel electrophoresis. Primers targeting the hypervariable V3-V4 region of the 16S rRNA gene were used to amplify the extracted DNA samples. The forward primer was 5’-CCTACGGGNGGCWGCAG-3’ and the reverse primer was 5’-GACTACHVGGGTATCTAATCC-3’. AxyPrep PCR Clean-up Kit (AP-PCR-500G; Corning, NY, USA) was used to separate, extract and purify the PCR products, and the products were quantified using a Quant-iT PicoGreen dsDNA Reagent (P7581, Thermo Scientific, Waltham, MA, USA). After quality determination, libraries passing quality control were sequenced with Novaseq sequencer for 2 x Two terminal sequencing of 250 bp at LC-Bio Co., Ltd.

Species annotation of the colonies was performed using the Greengene database v13.8, and then the ASV/OUT data of the colonies were extracted using the phyloseq package v1.26.1. We used the α-diversity index to characterize the diversity of the flora, where Shannon and Simpson indices were used to characterize species richness, homogeneity, and concentration reflecting species diversity, respectively. Beta diversity was calculated based on weighted Unifrac distances, and principal coordinate analysis (PCoA) was used in order to assess differences in microbial community composition. Linear discriminant analysis (LDA) was performed using the Mann-Whitney U test, and linear discriminant analysis effect size (LEfSe) analysis was performed using lefse software v1.0.0 to screen for species most likely to explain differences between groups, while LDA scores were used to assess effect sizes for species with significant differences between groups, with |LDA |> 2 and P < 0.05 as the thresholds of difference to screen for differences between species, and ggplot 2 software was used to assess differences in the composition of microbial communities. The results were also analyzed as bar graphs using the ggplot 2 software package v3.4.0. The results were presented as bar graphs. The α-diversity, β-diversity indices between the two groups were compared using Mann-Whitney U rank sum test through vegan software package v2.5.6. All the above analyses were carried out in R software v4.3.1, and the above P-values were two-tailed tests, and differences were considered statistically significant when P<0.05.

Paired tumor tissues from 108 GCs were subjected to mRNA sequencing (RNA-seq).In the end, 64 samples met the screening criteria. Total RNA was isolated from tumor tissues and NATs using TRIzol reagent (Invitrogen, Carlsbad, CA, USA) in an RNA protection solution. the amount and purity of RNA from each sample was quantified using a NanoDrop ND-1000 (NanoDrop, Wilmington, DE, USA). RNA integrity was assessed using an Agilent 2100 with a RIN>7.0. For mRNA sequencing, libraries were prepared on 1 μg of DNase I-treated total RNA using the TruSeq kit (Illumina) and processed for 150 bpb on the Illumina HiSeq X Ten instrument at LC-Bio Technology Co. (Hangzhou, China) on an Illumina HiSeq X Ten instrument with 150-bp paired-end sequencing. (Hangzhou, China) performed 150-bp paired-end sequencing on an Illumina HiSeq X Ten machine according to the protocol recommended by the vendor.

We aligned reads of all samples to the < research species > reference genome using HISAT2 (https://daehwankimlab.github.io/hisat2/, version:hisat2-2.0.4) package, which initially remove a portion of the reads based on quality information accompanying each read and then maps the reads to the reference genome. HISAT2 allows multiple alignments per read (up to 20 by default) and a maximum of two mismatch when mapping the reads to the reference. HISAT2 build a database of potential splice junctions and confirms these by comparing the previously unmapped reads against the database of putative junctions. The mapped reads of each sample were assembled using StringTie (http://ccb.jhu.edu/software/stringtie/, version:stringtie- 1.3.4d) with default parameters. Then, all transcriptomes from all samples were merged to reconstruct a comprehensive transcriptome using gffcompare software (http://ccb.jhu.edu/software/stringtie/gffcompare.shtml, version: gffcompare-0.9.8). After the final transcriptome was generated, StringTie and ballgown (http://www.bioconductor.org/packages/release/bioc/html/ballgown.html) were used to estimate the expression levels of all transcripts and perform expression abundance for mRNAs by calculating FPKM (fragment per kilobase of transcript per million mapped reads) value.

Differentially expressed genes (DEGs) were screened by DESeq261. genes with P< 0.05 and FC ≥ 2 or FC ≤ 0.5 were considered statistically significant DEGs. enriched functional pathways and modules were analyzed by using KEGG and CO databases. The Mann-Whitney U test was used to compare differences between groups.

CIBERSORT is a computational method for analyzing the composition of immune cells from RNA sequencing data based on the expression profiles of immune cell-specific genes and uses machine learning algorithms to analyze and classify the expression patterns of these genes. We use the CIBERSORT R-script v1.03 to construct a support vector regression-based model using the known expression data of the reference genes and to-be-estimated gene expression data of the mixed samples, constructed the optimization problem by the correlation matrix consisting of the cellular composition, and solved it in the form of a sparse solution. Thus, the cellular composition ratio of the mixed samples is estimated. The FPKM matrix obtained by transcriptome sequencing was transformed into a matrix of relative content of 22 different types and functional states of immune cells. The flora matrix was combined with the immune cell abundance matrix and the correlation coefficients between the columns in the combined matrix were calculated by calling the rcorr function. The type of correlation coefficient was Spearman’s correlation coefficient.

The samples were taken out of the -80°C freezer and thawed on ice, and metabolite were extracted with 80% methanol buffer. Briefly, 50 mg of sample was extracted with 0.5 ml of precooled 80% methanol. The extraction mixture was then stored in 30 min at -20°C. After centrifugation at 20,000 g for 15 min, the supernatants were transferred into new tube to and vacuum dried. The samples were redissolved with 100 μL 80% methanol and stored at -80°C prior to the LC-MS analysis. In addition, pooled QC samples were also prepared by combining 10 μL of each extraction mixture. The extracted samples were then sorted for machine analysis with randomization. QC samples were inserted before, in the middle, and after the samples to evaluate experimental technical replicates. The samples underwent mass spectrometry positive and negative ion scans. All samples were acquired by the LC-MS system followed machine orders. Firstly, all chromatographic separations were performed using an UltiMate 3000 UPLC System (Thermo Fisher Scientific, Bremen, Germany). An ACQUITY UPLC T3 column (100mm×2.1mm,1.8μm, Waters, Milford, USA) was used for the reversed phase separation. The column oven was maintained at 40°C. A high-resolution tandem mass spectrometer TripleTOF 6600 (SCIEX, Framingham, MA, USA) was used to detect metabolites eluted form the column. The Q-TOF was operated in both positive and negative ion modes. The curtain gas was set 30 PSl,lon source gas1 was set 60 PSI,lon source gas2 was set 60 PSI, and an interface heater temperature was 500°C. For positive ion mode, the lonspray voltage floating were set at 5000 V, respectively. For negative ion mode, the lonspray voltage floating were set at -4500V, respectively. The mass spectrometry data were acquired in IDA mode. The TOF mass range was from 60 to 1200 Da. The survey scans were acquired in 150 ms and as many as 12 production scans were collected if exceeding a threshold of 100 counts per second (counts/s) and with a 1⁺charge-state. Dynamic exclusion was set for 4s. During the acquisition, the mass accuracy was calibrated every 20 samples. Furthermore, in order to evaluate the stability of the LC-MS during the whole acquisition, a quality control sample (Pool of all samples) was acquired after every 10 samples.

The raw data from mass spectrometry were converted into readable data mzXML format using Proteowizard’s MSConvert software. XCMS software was utilized for peak extraction, and peak extraction quality control was conducted. Subsequently, substances extracted were annotated using CAMERA for adduct and ion annotation, followed by primary identification using the metaX software. Identification was performed separately using the mass spectrometry first-level information and matching the mass spectrometry second-level information with an in-house standard compound database. Differential metabolites were identified by Mann-Whitney U test and partial least squares discriminant analysis (PLS-DA).Metabolites with variable importance in projection (VIP) > 1 and p < 0.05 and FC ≥ 2 or FC ≤ 0.5 were considered differential metabolites. The functions of these metabolites and metabolic pathways were analyzed using the KEGG database.

Continuous variables with normal distribution are expressed as mean ± standard deviation (x ± s) or Mean ± SD and analyzed using t-test or Mann-Whitney U test. Categorical variables are presented as counts (n, %) and analyzed using Chi-square test or Fisher’s exact test. Propensity score matching (PSM) was used to account for differences in patient backgrounds, with a 1:4 ratio set to minimize selection bias between the two groups. Survival rates were calculated using the Kaplan-Meier method and survival curves were compared using the log-rank test. A Cox proportional hazards model with forward stepwise regression was employed to identify independent prognostic factors. Spearman correlation was used for the joint analysis of microbiome with transcriptome and metabolome. All data were analyzed using SPSS software version 26.0 (IBM USA), the Medsta statistical platform (www.medsta.cn/software), OmicStudio tools (https://www.omicstudio.cn/tool), and R version 4.3.1. All statistical tests were two-sided, and a p-value < 0.05 was considered statistically significant.

In this study, data from 5567 patients who met the criteria and had complete follow-up information were collected from 7192 hospitalized patients with GC ( Figure 1A ). There were no statistically significant differences between the two groups of BMI patients in terms of smoking history, alcohol consumption history, extent of resection, type of pathology, pM stage or recurrent metastasis (all P > 0.05). Analysis revealed that relative to NLBMI patients, LBMI patients had a smaller percentage of family history of GC; more tumors were located in the lower 1/3 and the whole stomach and less in the upper 1/3, and there was a greater percentage of larger and more poorly differentiated tumors, and a greater percentage of open surgeries (all P < 0.05); the level of pre-CA199 positivity was significantly greater (P=0.039), and the pre-CEA positivity level was similar; and the percentage of nerve invasion was greater (P=0.011), while there was no significant difference in vascular invasion. Moreover, in the LBMI group, the percentage of female patients aged ≥60 years, incidence of complications, deep tumor infiltration, high number of lymph node metastases, late pathological stage and low percentage of receiving postoperative adjuvant chemotherapy were significantly greater than that those in the NLBMI group (P<0.001). ( Supplementary Table S1 ).

Univariate and multivariate COX analysis revealed that LBMI is an independent poor prognostic factor for overall survival (OS) in GC patients (HR=1.28, 95%CI: 1.13-1.45, P<0.001) ( Table 1 ). Kaplan-Meier survival analysis based on BMI classification showed that LBMI patients had worse prognosis compared to NLBMI patients before PSM (5-yr OS: 50.8% vs. 66.2%, P<0.001) ( Figure 1B ). After adjusting for clinicopathological characteristics that influence prognosis (P<0.05) using PSM (ratio 1:4), the clinical characteristics of the two groups were comparable (P>0.05, Supplementary Table S2 ). Similarly, LBMI patients had worse prognosis (5-yr OS: 50.8% vs. 60.5%, P<0.001) ( Figure 1C ). Stratified analysis by TNM stage show no significant difference in OS between the two BMI groups in stage I and IV disease ( Figures 1D, G ); however, in stage II and III patients, LBMI disease have worse OS compared to NLBMI patients ( Figures 1E, F ). Subgroup analysis based on receiving postoperative chemotherapy set the PSM ratio to 1:4, and included clinicopathological data that influence prognosis. After PSM, the clinicopathological characteristics of the two BMI groups were comparable ( Supplementary Table S3 ), and LBMI patients had worse OS both before and after PSM ( Figures 1H, I ).

To evaluate whether there were differences in microbial diversity, abundance, and composition between LBMI and NLBMI gastric cancer patients, we included 189 eligible gastric cancer patients, including 27 in the LBMI cohort and 162 in the NLBMI cohort ( Supplementary Figure S1 ). As shown in Supplementary Table S4 , the clinicopathologic data were balanced and comparable between the two cohorts. On the basis of the species sparsity curves ( Supplementary Figures S2B, C ), we found that the curves of the four groups in both metrics flattened out. The Venn diagram ( Supplementary Figure S2A ) revealed that there are many overlaps in the microbial environments among the four groups. The diversity within the cancer tissues was significantly greater in both LBMI and NLBMI carcinomas than in the LBMI and NLBMI paracarcinomas, whereas there was no significant difference between the two cohorts of LBMI and NLBMI cancer tissues ( Figures 2A, B ). Principal coordinate analysis (PCoA) revealed significant differences in both BMI carcinomas and paracancers in both groups; however, there was no significant difference between LBMI-CT and NLBMI-CT (P=0.855) ( Figure 2C ). PLS-DA analysis, revealed that the intratumoral microbiome of two groups of BMI carcinomas could be divided into two different clusters ( Figure 2D ). Regarding species composition, the differences were a smaller between the tumor tissues of different BMI groups, while the differences were a greater between tumor and peritumoral tissues of the same BMI group ( Figures 2E, F ).

To determine the differentially dominant flora in GC patients with different BMIs, LEfSe analysis was performed (LDA> 2.0, P<0.05), which revealed 59 ( Supplementary Figure S3 ) and 230 ( Supplementary Figure S4 ) differentially dominant flora in the LBMI group and the NLBMI group, respectively, compared with the paracancerous tissues. There were 32 differentially dominant flora in the LBMI-CT group compared with the NLBMI-CT group. At the phylum level only p_Nitrospinae dominated the flora in LBMI, whereas at the genus level g_Lachnoanaerobaculum,g_Brevundimonas and g_Stomatobaculum dominated the flora in the NLBMI group, whereas g_Acidiphilium,g_Thiobacillus and g_Abiotrophia and 9 other genera were the dominant flora in the LBMI group. At the species level, s_Knoelia_sp_BA2_2011 and 15 other species were dominant flora in the LBMI group ( Figure 3A ; Supplementary Figure S5 ). At the genus level, the abundances of two groups of differentially bacteria, g_Abiotrophia and g_Lachnoanaerobaculum, significantly differed ( Figures 3B-J ). Spearman correlation analysis revealed that g_Abiotrophia was positively correlated with g_Lachnoanaerobaculum and g_Stomatobaculum and that g_Brevundimonas was negatively correlated. These findings suggest a possible complementary relationship between the dominant differential flora between the two BMI groups ( Supplementary Figure S2D ).

RNA sequencing analysis was performed on 64 tumor tissues from both groups, and PCA revealed that there was no significant difference in BMI the between the two groups (P=0.136) ( Figure 4A ). Compared with NLBMI, 343 genes were significantly upregulated and 320 genes in LBMI were significantly downregulated ( Figure 4B ). KEGG and GO analyses were performed on the BMI differential genes of the two groups, and KEGG analysis revealed that the LBMI group was enriched mainly in the Wnt signaling pathway, gastric cancer, and African trypanosomiasis ( Figure 4C ); similarly GO enrichment analysis was performed mainly in the extracellular region, extracellular space, plasma membrane and Wnt signaling pathway ( Figure 4D ). Correlation analysis of the DEGs associated with the dominant flora at the genus level revealed that g_Abiotrophia was significantly positively correlated with 11 genes, such as LGR6,and significantly negatively correlated with 30 genes, such as P2RY12 and SCN4B,in the LBMI group ( Figure 4E ; for details, see Additional File S1 ). The above results revealed that GC patients with different BMIs presented different transcriptomic landscapes and that many of these genes were closely related to differential intratumoral microbiota, suggesting that differential intratumoral microbiota may regulate the progression of GC by influencing the genes of the host.

To explore the associations between BMI-associated intratumoral microbiota and tumor-infiltrating immune cells, we analyzed the composition of immune cells in 64 GC patients via transcriptome sequencing information and BMI data and plotted bar graphs of immune cell abundance ( Figures 5A, B ) to discover the unique features of the TIME of GC patients with different BMIs.

Correlation analysis revealed that in the NLBMI group, g_Lachnoanaerobaculum showed a significant positive correlation with T cell follicular helper, while Mast cells resting exhibited a significant negative correlation. g_Stomatobaculum demonstrated a significant positive correlation with T cell follicular helper, whereas T cells CD4 memory resting showed a significant negative correlation. In the LBMI group, g_Enterobacter displayed a significant positive correlation with B cells naive, while Dendritic cells activated and T cells CD4 memory resting showed significant negative correlations.g_Abiotrophia showed a significant negative correlation with eosinophils ( Figure 5C ). The above results indicated that the BMI-related dominant intratumoral microbiota of GC patients were significantly associated with various tumor-infiltrating immune cells, suggesting that they may play a role in regulating the immune microenvironment of GC.

Untargeted metabolomic analysis was performed on 57 tumor tissues in the transcriptome( Supplementary Table S5 ), and a total of 2688 metabolites were identified, of which 122 metabolites were significantly different between the LBMI and NLBMI groups (P<0.05,FC ≥2 or FC ≤0.5) ( Figure 6A ), and the PLS-DA scoring plot revealed that the different metabolites in the LBMI versus NLBMI tumors could be classified into two different clusters (R2Y = 0.432,Q2Y = 0.368) ( Figure 6B ). Tests of the PLS-DA model revealed that R2 > Q2 and the Q2 regression line had a negative intercept (R2 = [0.0, 0.354],Q2 = [0.0, -0.421]) ( Figure 6C ). The heatmap revealed that compared to the NLBMI group, the LBMI group had higher abundance of intratumoral purine metabolites, such as idp ( Supplementary Figure S6 ).Differentially abundant metabolite KEGG enrichment analysis revealed that the LBMI group was enriched mainly in pathways such as purine metabolism and the caffeine metabolism pathway( Figure 6D ). Genus-level differential dominant bacteria and differentially abundant metabolite correlation analysis, as shown in Figure 6E , revealed that the differential dominant bacteria in the NLBMI group, g_Lachnoanaerobaculum, were significantly negatively correlated with 12 differentially abundant metabolites, such as 8-methoxykynurenate;g_Stomatobaculum was significantly negatively correlated; g_Brevundimonas was significantly positively correlated with eleutheroside b1 and 2-dodecylbenzenesulfonic acid and significantly negatively correlated with mimosine and latamoxef. g_Abiotrophia in the LBMI group presented a significant negative correlation with demethoxyfumitremorgin c, whereas it presented a significant positive correlation with guanine and idp; g_Dubosiella presented a significant positive correlation with caffeine and four others; g_Enterobacter presented a significant positive correlation with cyclic n-acetylserotonin glucuronide and 8-methoxykynurenate presented a significant positive correlation; g_Sphingopyxis showed significant negative correlation with 2-piperidinone;g_Prevotellaceae_UCG-001 and g_Methylocystis demonstrated a significant negative correlation with differentially abundant metabolites that were not significantly correlated(see Additional File S2 ). The above results revealed significant correlations between the two groups of intratumoral microbiota and metabolites, suggesting that they may further affect the biological process of gastric cancer by influencing metabolites.