Association between Gut Microbiota and Digestive System Cancers: A Bidirectional Two-Sample Mendelian Randomization Study

Accumulating evidence indicates that gut microbiota closely correlates with the tumorigenesis of digestive system cancers (DSCs). However, whether the causality between gut microbiota and DSCs exists is unknown. Genome-wide association study (GWAS) summary statistics for gut microbiota and DSCs and the bidirectional two-sample Mendelian randomization (MR) analysis were utilized to assess the causality between gut microbiota and DSCs. Sensitivity analyses were performed to evaluate the robustness of our results. We found that the genus Eggerthella (OR = 0.464, 95%CI: 0.27 to 0.796, p = 0.005) was negatively associated with the risk of gastric cancer. The genetically predicted genus Lachnospiraceae FCS020 group (OR = 0.607, 95%CI: 0.439 to 0.84, p = 0.003) correlated with a lower risk of colorectal cancer, and genus Turicibacter (OR = 0.271, 95%CI: 0.109 to 0.676, p = 0.005) was a protective factor for liver cancer. In the reverse MR, DSCs regulated the relative abundance of specific strains of gut microbiota. We comprehensively screened the association between gut microbiota and DSCs using a bidirectional two-sample MR analysis and identified the causality between several microbial taxa and DSCs. Our discoveries are beneficial for the development of novel microbial markers and microbiota-modifying therapeutics for DSC patients.


Introduction
Digestive system cancers (DSCs) comprise gastrointestinal cancers (esophageal cancer, gastric cancer, small intestinal cancer, colorectal cancer) and hepatobiliary-pancreatic cancers (liver cancer, gallbladder cancer, biliary duct cancer, pancreatic cancer). In general, DSCs are the most commonly diagnosed cancers and the leading cause of cancer deaths throughout the world, as they account for >26.4% of new cancer cases and >36.3% of new cancer deaths, according to the GLOBOCAN estimates of cancer incidence and mortality in 2020 by the International Agency for Research on Cancer [1]. Hence, to effectively increase life expectancy, it is a priority for us to further elucidate the specific mechanisms underlying the occurrence and progression of DSCs. To date, growing evidence has demonstrated that gut microbiota could regulate the tumorigenesis of DSCs [2,3].
The human gastrointestinal tract is a natural habitat for trillions of microorganisms, i.e., gut microbiota, which comprises bacteria, fungi, archaea, viruses, and yeast. With the rapid advances in multi-omics and sequencing technologies in recent years, our understanding of the complex relationships between gut microbiota and various disorders (e.g., metabolic diseases, gut diseases, and multiple cancers) has been greatly facilitated [4]. Notably, cumulative findings demonstrate that gut microbes and their metabolites not only regulate the development and progression of several DSCs but also modulate the DSC patients' responsiveness to chemotherapy, radiotherapy, molecular targeted therapy, and immunotherapy largely by affecting innate and adaptive immunity [2,5,6]. For instance, specific bacteria strains, e.g., Fusobacterium nucleatum, Bacteroides fragilis, and Escherichia coli, are demonstrated to be closely correlated with colorectal cancer (CRC) [7]. Furthermore, Porphyromonas gingivalis and Porphyromonas asaccharolytica could promote the onset of CRC by inducing butyrate-related cellular senescence [8]. In addition, a two-sample Mendelian randomization (MR) study reported that Ruminococcaceae, Porphyromonadaceae, and Bacteroidetes might be protective factors for liver cancer, and the researchers further validated the findings using a case-control study [9]. Moreover, Clostridium species can recruit NKT immune cells to reduce liver cancer growth and metastasis by mediating bile acid metabolism [10]. Gut microbiota and their metabolites have "double-edged sword" effects on DSC patients' responsiveness to pharmacological treatments. On the one hand, gut microbiota contributes to the development of drug resistance in DSC patients. For example, F. nucleatum promotes colorectal cancer chemoresistance by activating cell autophagy [11]. On the other hand, gut microbiota could boost the therapeutic efficacy and reduce the therapy-associated toxicity or side effects for DSC patients. In mouse models, Bifidobacterium pseudolongum was found to enhance immune checkpoint inhibitor efficacy by producing inosine to promote Th1 transcriptional differentiation [12]. Nonetheless, the current studies did not systematically examine the relationships between gut bacteria and different DSCs. More importantly, these observational studies on the association between gut microbiota and DSCs are not capable of identifying the causal link.
As a widely accepted statistical method, MR can exploit genetic variants, i.e., single nucleotide polymorphisms (SNPs), as instrumental variables (IVs) to explore the causality between interested exposures and specific outcomes [13]. As genetic variants are subject to random allocation at conception and almost remain static after birth, an MR study is capable of largely avoiding confounding factors and reverse causality in traditional observational studies [14]. To obtain a better understanding of the role of gut bacteria in DSCs, we conducted a bidirectional two-sample MR study to systematically investigate the possible causal associations between gut bacteria and diverse DSCs.

Study Design
A flowchart summarizing the bidirectional two-sample MR analysis that we used in this study is shown in Figure 1. The bidirectional two-sample MR analysis was utilized to explore the causal relationships between gut microbiota and cancers of the digestive system using the genetic variants as the IVs. Three assumptions need to be fulfilled to verify the validity of the MR analysis: (1) IVs correlate strongly with the exposure; (2) there is no correlation between IVs and any confounding variables; and (3) IVs could only affect the outcome via the exposure [15].

Data Sources and Study Population
The genome-wide association study (GWAS) statistics on gut microbiota were obtained from the IEU open GWAS database (https://gwas.mrcieu.ac.uk/, accessed on 3 October 2022), a database of 245,341,232,597 genetic associations from 42,335 GWAS summary datasets. A total of 211 bacterial traits (classified into specific phylum, class, order, family, and genus) were obtained, and the sample size was 14,306. Since 15 bacterial traits did not have specific species names, we excluded them and selected 196 bacterial traits for analysis. We initially selected seven cancers of the digestive system including esophagus cancer, gastric cancer, small intestine cancer, colorectal cancer, pancreatic cancer, gallbladder cancer, and liver cancer. The GWAS summary data on the cancers of the digestive system noted above were accessed from the IEU open GWAS project (https://gwas.mrcieu.ac.uk/, accessed on 3 October 2022). The details of the GWAS datasets we selected are shown in Table 1. Esophageal cancer was excluded from this study since the IVs in datasets did not meet the requirements of MR analysis. All participants in the above datasets except gallbladder cancer were of European ethnicity. Gallbladder cancer was excluded from this study due to the incompatibility of ethnicity in the dataset with other datasets. Ethical approval and consent to participate were obtained in all original studies.

Selection of IVs
Firstly, the SNPs in the GWAS summary data for exposures that were genetically associated with the traits with genome-wide significance (p < 5 × 10 −8 ) were sifted out as instrumental variables. In such cases where the number of IVs was fewer than 4, the significance threshold was relaxed to 5 × 10 −6 to avoid inaccurate results due to insufficient SNPs. Secondly, to exclude some undesirable SNPs (r 2 > 0.001, window size < 10,000 kb), linkage disequilibrium clumping was utilized. Thirdly, the IVs strongly associated with outcomes (p < 5 × 10 −8 ) were removed according to the third assumption of MR. Finally, harmonization of the exposure and outcome datasets was accomplished, and the palindromic SNPs with intermediate allele frequencies were removed.
The strength of genetic instruments for exposures was ensured by calculating the F statistic. The F statistics were computed using the admittedly reliable formula: , where R 2 and N refer to the cumulative explained variance in the selected SNPs and sample size, respectively. R 2 was calculated using the formula: where MAF refers to the minor allele frequency [13]. If F > 10, there is sufficient strength to avoid the problem of weak instrument bias in the two-sample model [14].

Statistical Analysis
The statistical analyses were performed in R 4.1.3 with the R packages "TwoSam-pleMR" (version 0.5.6) and "MRPRESSO" (version 1.0). p < 0.05 was considered statistical significance for evidence of a potential causal effect. The multiple methods for the main analyses (inverse variance weighted, weighted median, and MR Egger) and several sensitivity analyses (heterogeneity test, pleiotropy test, and leave-one-out sensitivity test) were utilized in our study [16,17]. We utilized the inverse-variance-weighted (IVW) method as the primary approach to estimate the causal effect. Weighted median and MR Egger were used as auxiliary methods. The sensitivity analyses were performed to evaluate the robustness of the findings. The heterogeneity test was performed using Cochran's Q test. For the pleiotropy test, the MR-Egger intercept and MR-pleiotropy residual sum and outlier (MR-PRESSO) methods were used to evaluate pleiotropy. The MR-Egger intercept was used to evaluate the potential pleiotropy in the IVW model, and MR-PRESSO was used for testing horizontal pleiotropy by identifying and removing outlying instrumental variables (NbDistribution = 3000, SignifThreshold = 0.05) [18]. The leave-one-out sensitivity test was also used to evaluate the robustness of this study's findings. The estimates were presented as odds ratio (OR) or β with their 95% confidence intervals (CIs) per one standard deviation (SD) increase in the exposure.

Instrumental Variables
To analyze the effects of cancers on gut microbiota, we selected five SNPs for GC, five SNPs for small intestinal cancer, nine SNPs for CRC, three SNPs for pancreatic cancer, and two SNPs for liver cancer as instrumental variables with the genome-wide significance of p < 5 × 10 −6 . To analyze the effects of gut microbiota on cancers, we selected at least two and up to twelve SNPs for gut microbiota species as instrumental variables. Some analyses failed because no SNPs were left after harmonization. The F statistics for IVs indicated that the estimates were less likely to suffer weak instrumental bias (F > 10, Supplemental Table S1).

The Bidirectional Causal Associations between Gut Microbiota and DSCs
We utilized bidirectional MR analysis to investigate the causality between gut microbiota and five types of cancers (i.e., GC, small intestinal cancer, CRC, pancreatic cancer, and liver cancer). All the significant results are shown in Tables 2 and 3.

Gastric Cancer
The estimates calculated using the IVW test indicated that the genetically predicted relative abundance of family Bacteroaceae

Liver Cancer
The estimates calculated using the IVW test indicated that the genetically predicted relative abundance of family Rhodospirillaceae The estimates calculated using the IVW test indicated that the genetic predisposition to liver cancer was negatively associated with the relative abundance of genus Butyricimonas

Sensitivity Analysis
We performed sensitivity analyses to evaluate the robustness of the findings for the exposures with over four SNPs as IVs. We found ten significant robust results after all the sensitivity analyses were performed (Figure 2 and Supplemental Table S2). The significant causal associations between the genus Eggerthella and the lower risk of GC, between the genus Lachnospiraceae FCS020 group and the lower risk of CRC, and between the genus Turicibacter and the lower risk of liver cancer were robust. Moreover, the significant causal associations between GC and the higher relative abundance of genus Dialister and genus Eubacterium ventriosum group, between small intestinal cancer and the higher relative abundance of genus Collinsella, genus Howardella, and genus Lachnospiraceae UCG008, between small intestinal cancer and the lower relative abundance of genus Peptococcus, and between CRC and the lower relative abundance of genus Lactococcus were robust. For these results, pleiotropic effects among the selected SNPs were not found using the MR-Egger and MR-PROSSO global tests (p > 0.05, Supplemental Table S2). No evidence of heterogeneity was observed between the selected IVs and cancers using Cochran's Q test (p > 0.05, Supplemental Table S2). The robust results verified using the leave-one-out sensitivity test showed that removing any SNP had no effect on these results (Figure 3). heterogeneity was observed between the selected IVs and cancers using Cochran's Q test (p > 0.05, Supplemental Table S2). The robust results verified using the leave-one-out sensitivity test showed that removing any SNP had no effect on these results (Figure 3).

Discussion
Our study is the first attempt to systematically analyze the causal associations between gut microbiota and DSCs using a bidirectional two-sample Mendelian Randomization study. We found that the genetically increased relative abundance of specific gut microbial genera was related to the lower risk of several DSCs, i.e., genus Eggerthella for GC, genus Lachnospiraceae FCS020 group for CRC, and genus Turicibacter for liver cancer, respectively. Furthermore, the reverse MR analysis identified significant associations between genetic predisposition to DSCs and the relative abundance of several bacterial genera, i.e., a positive association between GC and genus Dialister and genus Eubacterium ventriosum

Discussion
Our study is the first attempt to systematically analyze the causal associations between gut microbiota and DSCs using a bidirectional two-sample Mendelian Randomization study. We found that the genetically increased relative abundance of specific gut microbial genera was related to the lower risk of several DSCs, i.e., genus Eggerthella for GC, genus Lachnospiraceae FCS020 group for CRC, and genus Turicibacter for liver cancer, respectively. Furthermore, the reverse MR analysis identified significant associations between genetic predisposition to DSCs and the relative abundance of several bacterial genera, i.e., a positive association between GC and genus Dialister and genus Eubacterium ventriosum group, a positive association between small intestinal cancer and genus Collinsella, genus Howardella and genus Lachnospiraceae UCG008, a negative association between small intestinal cancer and genus Peptococcus, and a negative association between CRC and genus Lactococcus.
Gut dysbiosis is widely recognized to be implicated in DSCs, and current studies have partially discovered how the dysregulation of gut microbiota promotes the occurrence and development of DSCs [19]. On the one hand, microbiota dysbiosis can induce the oncogenic transformation of host cells with direct effects on them. For instance, campylobacter jejuni can directly evoke DNA double-strand breaks by producing a genotoxin with DNase activity, a cytolethal-distending toxin, which aids the onset of CRC [20]. In addition, microbiota taxa such as Helicobacter pylori might facilitate genomic instability by inhibiting DNA mismatch repair, which contributes to gastric carcinogenesis [21]. On the other hand, gut microbiota dysregulation could shape anti-tumor immunity by exerting effects on the innate and adaptive immune responses implicated in DSCs. Fusobacterium nucleatum could suppress anti-tumor immunity by directly interacting with the inhibitory T-cell receptor TIGIT via FAP2 and inhibiting natural killer (NK) cell-mediated tumor killing [22,23]. Additionally, metabolites synthesized or processed by gut microbiota, e.g., short-chain fatty acids (SCFAs) and secondary bile acids (BAs), can play an immunomodulatory role. As the representative products of dietary fiber fermentation by gut microorganisms, SCFAs (e.g., butyrate, acetate, and propionate) can not only promote the generation of regulatory T cells by impeding histone deacetylase (HDAC) activity but also enhance the expansion of macrophage precursors and the function of CD8 + T cells by regulating cellular metabolism [24,25], which results in the paradoxical role of SCFAs in tumors [7].
Our results, for the first time, suggested that the genus Eggerthella had protective effects against GC. Recently, mounting evidence has jointly underlined that Eggerthella was associated with various diseases, such as asthma [26], multiple sclerosis [27], systemic lupus erythematosus [28], rheumatoid arthritis [29], and CRC [30]. For instance, in cigarette smoke-exposed mice, Eggerthella lenta might increase the bile acid metabolite, i.e., taurodeoxycholic acid (TDCA), to activate oncogenic MAPK/ERK signaling and induce gut barrier dysfunction, thus playing a protumorigenic role in CRC [30]. Since bile acid and its metabolites could influence gastric carcinogenesis [31], the genus Eggerthella is probably implicated in GC, but its exact role needs further exploration.
Additionally, we found that the genus Lachnospiraceae FCS020 group was a protective factor against CRC. Likewise, researchers found that the abundance of genus Lachnospira, a member of the Lachnospiraceae family, was reduced in CRC patients [32], and the Lachnospiraceae FCS020 group was negatively related to lymph node metastasis in CRC [33]. Furthermore, the Lachnospiraceae FCS020 group was reported to be linked with lower plasma trimethylamine N-oxide levels [34], which was a risk factor for CRC [35] and might facilitate CRC by inducing chronic inflammation, oxidative stress, and DNA damage [36]. The Lachnospiraceae family also may exert protective effects against CRC by enabling the immune surveillance of CD8 + T cells [37] and inhibiting the colonization of CRC-associated oral bacteria [38]. Although most studies report the anti-cancer effects of Lachnospiraceae, some studies suggest the opposite results. Secondary bile acids (BAs), including taurodeoxycholic acid and deoxycholic acid, are identified as canonical carcinogenic bile acids [30,31,39,40]. Hence, it is plausible that Lachnospiraceae plays a tumorigenic role in CRC as it could convert primary bile acids into secondary bile acids to increase the abundance of secondary BAs [39]. Moreover, as discussed above, SCFAs could play a dual role (pro-carcinogenic and anti-cancer effects) in the development of tumors. Considering the SCFA-producing capability of Lachnospiraceae and the aforementioned research results, its specific effects on CRC remain to be further elucidated.
Moreover, we found that the genus Turicibacter was negatively associated with the risk of liver cancer. However, the association between the genus Turicibacter and liver cancer has not been reported or fully defined in the published literature. Turicibacter has been proposed as a beneficial bacterium with anti-inflammation effects [41]. In HBV-CLD (chronic hepatitis B virus infection-associated liver diseases) patients, Turicibacter was inversely cor-related with serum aspartate aminotransferase, total bilirubin, and direct bilirubin, which suggested its negative association with disease progression [42]. Researchers demonstrated that the genus Turicibacter markedly increased in NAFLD (nonalcoholic fatty liver disease) patients with significant liver fibrosis compared to those without liver fibrosis [43], but they did not prove cause and effect. Since HBV-CLD and NAFLD are significant risk factors and common causes of liver cancer, we presume that the genus Turicibacter may contribute to the development of liver cancer with its effects on these two liver diseases, which warrants further mechanistic studies.
Using the reverse MR analysis, we discovered that GC, small intestine cancer, and CRC were associated with the abundance of several specific gut bacterial genera. For instance, we found that GC was markedly correlated with the increased abundance of genus Dialister in the gut. This was supported by the published findings that the relative abundance of Dialister was significantly higher in GC than in other benign gastric lesions [44,45]. Additionally, Eubacterium ventriosum group was first identified as a GC-related microbe in our study, which was reported to be associated with a lower risk of IBD [46]. Notably, we revealed several previously unreported small intestine cancer-related bacteria. Specifically, for the first time, we identified that the genus Collinsella, genus Howardella, and genus Lachnospiraceae UCG008 increased and that the genus Peptococcus decreased in small intestine cancer patients, in spite of the reported associations between Collinsella and anti-PD-1 efficacy [47], between Collinsella and oesophageal adenocarcinoma [48], between Lachnospiraceae and anti-PD-1 efficacy [49], between Lachnospiraceae and CRC [37], and between Peptococcus and CRC [50]. Moreover, our study showed that CRC was associated with a lower abundance of genus Lactococcus. However, contradictory findings about the genus Lactococcus in CRC exist. Researchers found that Lactococcus lactis, a bacterium belonging to the genus Lactococcus, was more abundant in normal colorectal tissues than in CRC tissues [51], while the others demonstrated that genus Lactococcus was enriched in CRC patients compared with healthy individuals [52].
Accumulating evidence suggests that gut microbiota are crucial drivers and potential therapeutic targets for DSCs. A series of clinical studies demonstrate that antibiotics (e.g., vancomycin), fecal microbiota transplantation, probiotics, and dietary interventions may be conducive to improving the efficacy of DSC therapies and reducing treatment complications [53]. For instance, fecal microbiota transplantation could reconstitute the gut microbiota and induce a relative increase in regulatory T cells within the colonic mucosa to relieve immune checkpoint inhibitor-associated colitis [54]. Lactobacillus rhamnosus GG is able to inhibit the proliferation of colon cancer cells and exerts antitumor activity by sensitizing cancer cells to 5-Fluorouracil and Irinotecan [55]. Oral administration of live Lactobacillus rhamnosus GG has also been found to increase tumor-infiltrating DCs and T cells and further augment the antitumor activity of anti-programmed cell death 1 immunotherapy [56]. Noteworthily, some bacterial taxa could enhance the antitumor effects of cancer treatments, while others could impair the efficacy of antitumor agents. Hence, it remains a research priority for us to achieve a better understanding of the specific gut microbes involved in DSCs and of the mechanisms underlying their action, which may be translated into possible clinical benefits for DSC patients in the future. Nevertheless, some limitations still existed in our study. Firstly, as our GWAS data on gut microbiota are from phylum to genus levels, we may have missed some specific bacteria that causally correlate to DSCs at a more specialized level such as the species or strain level. Secondly, all participants in the datasets we used were of European ethnicity, which limits the extrapolation of our findings to other ethnicities. Thirdly, the number of SNPs available for some microbial taxa is limited, which may have led to biased results.

Conclusions
We accomplished a comprehensive screening of gut microbiota implicated in DSCs and identified the causal relationship between several microbial taxa and DSCs using a bidirectional two-sample MR analysis. Our discoveries can provide more theoretical support for modulating the composition of gut microbiota to improve the treatment efficacy and promote the development of microbiota-based therapeutic strategies as well as microbial biomarkers for DSCs.