Mucosa-Colonizing Microbiota Correlate With Host Autophagy Signaling in Patients With Inflammatory Bowel Disease

Abstract

Both bacteria and autophagy are implicated in inflammatory bowel disease (IBD) pathogenesis. However, how bacteria crosstalk with autophagy signaling remains largely known, especially in intestinal mucosa. This study aimed to profile the internal complex autophagy signaling cascade and their external correlation with these bacteria, and consequently provide a systematic and precise target for future IBD diagnosis and therapy. We found the Ulcerative colitis (UC) patients exhibited more severe dysbiosis than the Crohn’s disease (CD) patients, as represented by alpha diversity, community phenotypes, and functional annotation compared with the control population. Meanwhile, CD patients showed greater transcriptional signaling activities of autophagy, endoplasmic reticulum (ER) stress, and bile acid production. Dominant bacteria (e.g., Rhodococcus, Escherichia, Shigella, and Enterococcus) were positively correlated and low-abundance bacteria (e.g., Bacillus, Acidovorax, Acinetobacter, and Stenotrophomonas) were negatively correlated with the autophagy signaling cascade (184 autophagy genes, 52 ER stress genes, and 22 bile acid production genes). Our observations suggested UC patients showed temporary and widespread microbiota turbulence and CD patients showed processive and local autophagy activity during IBD progression. Intestinal mucosa-colonizing bacteria were correlated with the bile/ER stress/autophagy signaling axis in IBD pathogenesis.

Introduction

Intestinal microbes and their metabolic products play key roles in human diseases, especially inflammatory bowel disease (IBD). Meanwhile, both adherent-invasive and diffusely adherent bacteria stick to the intestinal mucosa and induce host immune cell activity and are linked to Crohn’s disease (CD) and ulcerative colitis (UC), respectively (Mirsepasi-Lauridsen et al., 2019). Notably, mucosa samples can directly reflect bacteria-host interactions during IBD progression. However, few multi-omics studies have investigated these interactions using mucosal biopsies. Furthermore, both the classification and mechanistic dissection of pathways involved in IBD remain challenging due to the complex correlative and interactive networks between host genetics and microbes (Plichta et al., 2019). Therefore, mucosal biopsy specimens are required to further explore these networks between host genetic factors and microbiota.

Genome-wide association studies have indicated that autophagy is an important mechanistic dissection of IBD pathogenesis, especially that of CD (Nguyen et al., 2013). Autophagy not only decreases intestinal epithelial permeability by inducing lysosomal degradation of tight junction proteins (Nighot et al., 2015), but also modulates programmed cell death in the intestinal epithelium (Matsuzawa-Ishimoto et al., 2017). In addition, autophagy is responsible for the elimination of intercellular bacteria from endoplasmic reticulum (ER) function defect-induced intestinal barrier leakage. For instance, autophagy can be induced by Autophagy related 16 like 1 (ATG16L1)-Transmembrane protein 59 (TMEM59) interactions in response to bacterial infection (Boada-Romero et al., 2016). Furthermore, activation of the bile acid receptor can strongly suppress the induction of autophagy (Lee et al., 2014). Thus, the profiles of mucosal bacteria-correlated autophagy gene networks and their interactions with ER stress and bile production signaling should be clarified in IBD patients.

Here, we collected intestinal mucosal biopsies from control population and IBD patients for dual-omics analysis of mucosal transcriptome and mucosa-colonizing bacterial diversity. Analysis focused on the autophagy signaling cascade and their complex correlations with intestinal mucosa-colonizing bacteria in IBD patients and control population. Combined with previous studies, our analysis suggested that mucosal bacteria regulated both ER stress and bile acid production and consequent autophagic activity, finally acting on host IBD progression. We systematically profiled the mucosal bacteria-autophagy correlation network in IBD patients, especially in those with CD, and provided a reliable analysis model of microbe-correlated intestinal diseases.

Materials and Methods

Study Design

The control populations and IBD patients were recruited from the First People’s Hospital of Yunnan Province, China, from October 2017 to December 2019. Each participant provided information on age, occupation, and smoking and alcohol drinking history. The control population had no history of digestive tract disease or serious medical illnesses. Each IBD patient was diagnosed with UC or CD and had received no IBD therapy or antibiotics (within 3 months). All participants provided signed informed consents and completed a questionnaire regarding their age, sex, occupation, and antibiotic use, with family assistance if necessary.

Diagnosis of IBD

We used clinical feature inquiry, laboratory examination, and endoscopy and biopsy histopathological analysis to diagnose UC or CD according to clinical practice guidelines (Gomollón et al., 2017; Magro et al., 2017; Inflammatory Bowel Disease Group, Chinese Society of Gastroenterology, Chinese Medical Association, 2021). Firstly, outpatients with enteric symptoms of abdominal pain, diarrhea, blood in stools, and loose stools for more than 6 weeks, and/or parenteral manifestations, such as fatty liver and cholecystitis, were suggested for laboratory examination, which included a full blood count, electrolyte, liver function, and inflammatory marker [C-reactive protein (CRP)] tests. Secondly, when clinical features and laboratory tests suggested further investigation, an endoscopy (colon and/or intestinal endoscopy) was performed to obtain evidence and collect mucosal biopsies for histopathological analysis. Thirdly, histopathological investigations were performed to provide reliable information on mucosal architecture, lamina propria cellularity, neutrophil granulocyte infiltration, and epithelial abnormality. Finally, UC or CD was diagnosed when comprehensive analysis of clinical features, laboratory examinations, endoscopy, and histopathology gave positive responses.

Exclusion Criteria and Control Population Selection

To eliminate the possible effects of other intestinal diseases on UC and CD diagnosis, we carefully discriminated infectious enteritis, Clostridium difficile infection, intestinal tuberculosis, Behçet syndrome, amebiasis, schistosomiasis, and other intestinal diseases from IBD. Those participants without IBD or other digestive tract diseases and without serious diseases of other tissues and organs were chosen as the control population group. A total of 138 participants were initially recruited, 38 of which were excluded according to exclusion criteria (i.e., UC and/or CD and antibiotic use within the last 3 months). Finally, 25 control population, 26 CD patients, and 51 UC patients were enrolled in the study.

Sample Collection

The IBD patients and control population were first given an intravenous injection of propofol and etomidate, and then an enteroscopy was performed for mucosal sampling. All samples were collected from October 2017 to December 2019. Upon collection, the mucosal samples were immediately placed on ice and frozen at −80°C within 1 h for microbiome and transcriptome analyses.

Bacterial DNA Isolation and 16S rRNA Gene Sequencing

Microbial community genomic DNA was extracted from intestinal mucosa using a QIAamp DNA Mini Kit (Qiagen, Germany) according to the manufacturer’s instructions. The DNA extract was checked on 1% agarose gel, and DNA concentration and purity were determined using a NanoDrop 2000 UV–vis spectrophotometer (Thermo Scientific, Wilmington, United States). The hypervariable V3–V4 region of the bacterial 16S rRNA gene was amplified with primer pairs 338F (5′-ACTCCTACGGGAGGCAGCAG-3′) and 806R (5′-GGACTACHVGGGTWTCTAAT-3′) using an ABI GeneAmp® 9700 PCR thermocycler (ABI, CA, United States). PCR amplification of the 16S rRNA gene was performed as follows: initial denaturation at 95°C for 3 min, followed by 27 cycles of denaturing at 95°C for 30 s, annealing at 55°C for 30 s, extension at 72°C for 45 s, single extension at 72°C for 10 min, and end at 4°C. The PCR mixture contained 5× TransStart FastPfu buffer 4 μl, 2.5 mM dNTPs 2 μl, forward primer (5 μM) 0.8 μl, reverse primer (5 μM) 0.8 μl, TransStart FastPfu DNA Polymerase 0.4 μl, template DNA 10 ng, and finally ddH₂O up to 20 μl. PCR was performed in triplicate. The PCR products were extracted from 2% agarose gel and purified using an AxyPrep DNA Gel Extraction Kit (Axygen Biosciences, Union City, CA, United States) according to the manufacturer’s instructions and quantified using a Quantus™ Fluorometer (Promega, United States).

Purified amplicons were pooled in equimolar concentrations and paired-end sequenced on the Illumina MiSeq PE300/NovaSeq PE250 platforms (Illumina, San Diego, United States) using standard protocols.

The raw 16S rRNA gene sequencing reads were demultiplexed, quality-filtered using fastp v0.20.0 (Chen et al., 2018), and merged using FLASH v1.2.7 (Magoč and Salzberg, 2011) with the following criteria: (i) 300-bp reads were truncated at any site receiving an average quality score of <20 over a 50-bp sliding window, with truncated reads shorter than 50 bp and reads containing ambiguous characters discarded; (ii) only overlapping sequences longer than 10 bp were assembled according to their overlapping sequence. The maximum mismatch ratio of the overlapping region was 0.2. Reads that could not be assembled were discarded; and (iii) Samples were distinguished according to the barcode and primers, and the sequence direction was adjusted, exact barcode matching, two nucleotide mismatches in primer matching.

Operational taxonomic units (OTUs) with 97% similarity cutoff (Stackebrandt and Goebel, 1994; Edgar, 2013) were clustered using UPARSE v7.1, and chimeric sequences were identified and removed. The taxonomy of each OTU representative sequence was analyzed using RDP Classifier v2.2 (Wang et al., 2007) against the 16S rRNA database (Silva v132) with a confidence threshold of 0.7. To remove possible contamination, we sequenced three samples of pure water as a negative control with the same procedures, including DNA extraction, PCR amplification, cDNA library construction, and final sequencing. Based on the negative control sequencing results, the same DNA sequence count detected in the negative control was extracted from all control and IBD samples. To confirm the reliability of 16S rRNA sequencing data, we also performed real-time PCR and found the data are replicable.

Bioinformatics Analysis of 16S rRNA Gene Sequencing Data

Alpha Diversity Analysis

Community richness (Sobs) and diversity (Shannon) indices were used to estimate α diversity and examined using the Welch’s t-test. Relative abundance of the top 47 genera was displayed using a histogram.

Community Hierarchical Clustering

Community clustering heatmap was conducted to show community variation at the genus level (top 47). Multiple sample comparisons were performed using Tukey–Kramer one-way ANOVA with false-discovery rate (FDR) correction and 95% confidence. Significance was assumed for adjusted values of p ≤ 0.05 (R software, vegan package).

Analysis of Abundance Differences

Significant differences between control population and IBD patients were determined using the Kruskal-Wallis H test at genus level. The FDR and Tukey–Kramer methods (CI = 0.95) were used for multiple testing correction and post hoc tests, respectively.

Functional Microbial Composition Analyses

The Functional Annotation of Prokaryotic Taxa (FAPROTAX) database extrapolates functions of cultured prokaryotes to estimate metabolic and other ecologically relevant functions (Louca et al., 2016). We used FAPROTAX to assess how the metabolic activities of intestinal mucosa-colonizing bacteria affect the host signaling network using the Tukey–Kramer method (CI = 0.95).

Phenotypic Prediction of Intestinal Mucosa-Colonizing Bacteria

Intestinal mucosa-colonizing bacteria phenotypes were predicted and compared using BugBase (Ward et al., 2017) and the Kruskal-Wallis H test. Briefly, BugBase uses an OTU table as an input file, which is normalized by the predicted 16S copy number. The preprocessed database and BugBase tool then automatically select thresholds to predict bacterial phenotypes.

RNA Extraction and Sequencing Data Processing

RNA Extraction

Total RNA was extracted from tissue using TRIzol® Reagent according to the manufacturer’s instructions (Invitrogen), and genomic DNA was removed using DNase I (Takara). RNA quality was then determined using a 2100 Bioanalyzer (Agilent) and quantified using a ND-2000 spectrophotometer (NanoDrop Technologies). Only high-quality RNA samples (OD260/280 = 1.8–2.2, OD260/230 ≥ 2.0, RIN ≥ 6.5, and 28S:18S ≥ 1.0, >1 μg) were used to construct the sequencing library.

Library Preparation and Illumina Hiseq X Ten/NovaSeq 6000 Sequencing

The RNA-seq transcriptome library was prepared using a TruSeq™ RNA Sample Preparation Kit (Illumina, San Diego, CA, United States) with 1 μg of total RNA. In brief, messenger RNA (mRNA) was isolated according to the polyA selection method by oligo(dT) beads and then fragmented using fragmentation buffer. The double-stranded cDNA was synthesized using a SuperScript Double-Stranded cDNA Synthesis Kit (Invitrogen, CA, United States) with random hexamer primers (Illumina). The synthesized cDNA was then subjected to end repair, phosphorylation, and “A” base addition according to the library construction protocols of Illumina. Libraries were size-selected for 300-bp cDNA target fragments on 2% low-range ultra-agarose gel, followed by PCR amplification using Phusion DNA polymerase (NEB) for 15 PCR cycles. After quantification by TBS380, the paired-end RNA-seq library was sequenced using an Illumina HiSeq X Ten/NovaSeq 6000 sequencer (2 × 150-bp read lengths).

Read Mapping

The raw paired-end reads were trimmed and quality controlled using SeqPrep¹ and Sickle² with default parameters. The clean reads were then separately aligned to the reference genome in orientation mode using HISAT2³ software (Kim et al., 2015). The mapped reads of each sample were assembled using StringTie⁴ with a reference-based approach (Pertea et al., 2015).

Differential Expression Analysis and Functional Enrichment

To identify differentially expressed genes (DEGs) between two samples, the expression level of each transcript was calculated according to the transcripts per million reads (TPM) method. RSEM⁵ was used to quantify gene abundances (Li and Dewey, 2011). Essentially, differential expression analysis was performed using DESeq2 (Love et al., 2014)/EdgeR (Robinson et al., 2010) with Q value ≤ 0.05 and DEGs with |log2FC| > 1 and Q value ≤ 0.05 (DESeq2 or EdgeR)/Q value ≤ 0.001 (DEGseq) deemed significant. In addition, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) functional enrichment analyses were performed to identify DEGs significantly enriched in GO terms and metabolic pathways (Bonferroni-corrected value of p ≤ 0.05) compared to the whole transcriptome background. GO functional enrichment and KEGG pathway analysis were carried out using GOATOOLS⁶ and KOBAS (http://kobas.cbi.pku.edu.cn/home.do; Xie et al., 2011).

Alternative Splicing Event Identification

All alternative splicing events that occurred in the samples were identified using the recently released program rMATS (http://rnaseq-mats.sourceforge.net/index.html; Shen et al., 2014). Only isoforms similar with the reference or comprising novel splice junctions were considered; and splicing differences were detected as exon inclusion and exclusion, alternative 5′ and 3′, and intron retention events.

Transcriptome Analysis

Time Series Expression Trend Analysis

STEM is used for clustering, comparing, and visualizing time series gene expression data (Ernst and Bar-Joseph, 2006). The STEM clustering algorithm is a supervised algorithm, that is, clustering is classified into artificially set trends. First, the software simulates n of the most representative possible trends according to the preset, and then calculates the correlation coefficient between each gene and the preset trends. Finally, each gene is classified into the trend to which it is most similar. We used STEM (v1.3.11) to explore the gene expression pattern of the healthy controls and CD and UC patients. Analysis explored functional enrichment of genes with a certain expression pattern and predicted the genetic regulatory network of intestinal mucosa.

Gene Expression Correlation Analysis

Genes obtained from expression trend analysis and/or related to autophagy, ER stress, and bile acid production were used to calculate Spearman correlation coefficients based on the correlation of gene expression. In the correlation networks, the larger the node, the greater the number of correlations between the expression of the gene and other genes.

Protein–Protein Interaction Analysis

We used the STRING database⁷ to perform Protein–Protein Interaction (PPI) network analysis of genes of interest. The interactions corresponding to genes of interest were directly extracted from the database to construct the network. NetworkX (Python) was used to visualize the network of the genes of interest.

KEGG Functional Enrichment Analysis

We used R script and Fisher’s precision probability test to perform KEGG pathway enrichment analysis. To control the false positives of enrichment analysis, multiple testing corrections (Benjamini-Hochberg, BH) were carried out. KEGG pathways reaching a corrected p value of 0.05 were defined as significantly enriched.

Ethics Approval

All study protocols and procedures were approved by the Medical Ethics Board of the First People’s Hospital of Yunnan Province (GXBSC-2021001, 2021 updated), China, and were carried out in accordance with all relevant provincial, national, and international guidelines, including the Declaration of Helsinki. Written informed consent was obtained from all participants prior to their inclusion in the study.

Statistical Analysis

GraphPad Prism v9.0.0 (San Diego, CA, United States) and the R stats package were used to analyze all data (R Core Team, 2013). The Kruskal-Wallis H test with FDR correction (CI = 0.95) for multiple comparisons was conducted to compare abundant bacterial taxa. Spearman rank-order correlation was used to evaluate associations between genus-level bacterial relative abundances and autophagy-related gene expression levels. In all analyses, p < 0.05 was considered statistically significant.

Results

Characteristics of IBD Patients and Controls

We initially recruited 138 participants, 38 of which were excluded due to antibiotic use within the last 3 months and/or suffering from other serious illnesses. In total, 23 healthy control population, 26 CD patients, and 51 UC patients were enrolled (Table 1). Most IBD patients were male (68 out of 77 subjects) and had no history of smoking. The control population was recruited from outpatients suffering abdominal discomfort, but with no gastrointestinal or other serious illness. For transcriptome analysis, an additional biopsy sample was collected from six control population, five CD patients, and 12 UC patients (Supplementary Table S1) with informed consent.

IBD Patients, Especially UC Patients, Show Intestinal Mucosa-Colonizing Microbial Dysbiosis and Corresponding Dysfunction

The intestinal microbiota play a key role in IBD progression. IBD patients exhibit clear dysbiosis, mainly represented by a decrease in bacterial diversity (Ott et al., 2004; Manichanh et al., 2006). Mucosal biopsies are advantageous as they can be used to explore the real interactions between dysbiosis and host responses. Here, we collected intestinal mucosal biopsies from control population and IBD patients to investigate the compositional movement and functional variation potential of intestinal mucosa-colonizing bacteria. Consistent with previous study (Ott et al., 2004), we found that alpha diversity (Sobs and Shannon indices) of the intestinal mucosal bacterial community was significantly decreased in IBD patients compared with the control population. Furthermore, UC patients showed poorer alpha diversity than CD patients (Figures 1A,B). The decrease in the Sobs index mainly originated from the reduced abundances of Proteobacteria and Bacteroidetes (Supplementary Figure S1). Correspondingly, FAPROTAX analysis indicated that the metabolic activities of nitrite, nitrate, and fumarate were significantly decreased in IBD patients, especially UC patients (Figure 1D; Supplementary Figure S3). Notably, these metabolic activities had brought focuses because of their key roles in IBD incidence, prevention, and therapy (Dykhuizen et al., 1996; Saijo et al., 2010; Jädert et al., 2014; Casili et al., 2016). Redox of these nitrogenous salts is also associated with IBD pathology (Bourgonje et al., 2020). Based on community phenotype analysis, we also found that IBD patients showed a decrease in oxidative stress tolerance, primarily due to the lower abundances of Escherichia, Shigella, and Pseudomonas (Figure 1E; Supplementary Figures S4A,B). Unexpectedly, pathogen causing diarrhea and gastroenteritis potentials also decreased significantly in IBD patients. Potential pathogens, such as Rhodococcus, Streptococcus, Enterococcus, Veillonella, Ruminococcus torques, and Ruminococcus gnavus, were enriched in the intestinal mucosa of IBD patients (Figure 1C; Supplementary Figure S2), which contributed to high Gram-positive performance, another community phenotype (Supplementary Figures S4C,D).

Figure 1

Active Intestinal Mucosa Signaling in CD Patients Exceeds That in UC Patients

Functional genomics network analysis is a powerful tool for identifying the key regulatory networks involved in IBD progression (Peters et al., 2017). Here, gene expression clustering analysis revealed that the intestinal mucosa signaling network was comprised of eight clusters based on differential expression in control population and UC and CD patients (Figure 2A). Three clusters, i.e., cluster 4 (3,339 genes, p = 0.0000012), cluster 6 (6,355 genes, p = 1.2e–156), and cluster 7 (9,773 genes, p = 0), were significant according to the differential transcription levels in the control population and UC and CD patients. Cluster 6 genes were highly expressed in the intestinal mucosa of both UC and CD patients, whereas cluster 7 and 4 genes showed a gradual decrease in expression in UC patients (Figures 2B–D). Thus, these findings indicate that intestinal mucosa signaling activity was higher in CD patients than in UC patients.

Figure 2

Environmental Factor Correlation Analysis Reveals Bacterial Community-Matched Autophagy Signaling in Intestinal Mucosa of IBD Patients

Increasing studies emphasize that crosstalk between the mucosal microbiome and host signaling network greatly affects IBD progression and clinical outcome (Morgan et al., 2015; Ryan et al., 2020). Recent review work provides new insights of the interplay between autophagy and intestinal bacteria and suggested that IBD-associated autophagy alleles and their interactions with environmental triggers, such as resident microbiota, are crucial for developing new therapeutic strategies for IBD treatment (Larabi et al., 2020). Thus, we performed correlation analysis between intestinal mucosa-colonizing bacteria and autophagy gene expression levels in clusters 6 and 7, which covered most KEGG autophagy-related genes (Supplementary Figure S7). In this analysis, 144 autophagy genes (58 in cluster 6 and 86 in cluster 7) were significantly correlated with intestinal mucosa-colonizing bacteria (Figures 3A,D). Gene expression heatmap analysis showed that these bacteria-correlated autophagy genes were more highly expressed in CD patients than in UC patients (Supplementary Figure S8). Furthermore, these autophagy genes not only formed a correlation network, but also a PPI network. For example, bacteria-correlated autophagy genes Autophagy related 4A cysteine peptidase (ATG4A), Autophagy related 5 (ATG5), ATG16L1, Autophagy related 9A (ATG9A), Signal transducer and activator of transcription 3 (STAT3), NRAS proto-oncogene, GTPase (NRAS), and Heat shock protein 90 alpha family class A member 1 (HSP90AA1) in cluster 6 (Figures 3B,C) and Beclin 1 (BECN1), Mitogen-activated protein kinase 1 (MAPK1), Phosphatidylinositol 3-kinase catalytic subunit type 3 (PIK3C3), Phosphatidylinositol-4,5-bisphosphate 3-kinase catalytic subunit alpha (PIK3CA), Charged multivesicular body protein 2A (CHMP2A), Charged multivesicular body protein 3 (CHMP3), Charged multivesicular body protein 4A (CHMP4A), Vacuolar protein sorting 4 homolog A (VPS4A), Vacuolar protein sorting 25 homolog (VPS25), Translocase Of outer mitochondrial membrane 20 (TOMM20), and AKT serine/threonine kinase 2 (AKT2) in cluster 7 (Figures 3E,F) were core members of both the correlation and interaction networks.

Figure 3

The above bacteria could also be clustered into two groups, i.e., those showing positive correlations with autophagy genes in both cluster 6 and cluster 7 (e.g., Haemophilus, Ruminococcus torques, and Rhodococcus) and those showing negative correlations (e.g., Bacillus, Acidovorax, Acinetobacter, and Stenotrophomonas; Figures 3A,D). It is worth noting that most of the intestinal mucosa-colonizing bacteria showing positive correlation with the above autophagy genes have been reported in previous IBD studies (Goyal et al., 2018; Soltys et al., 2020). For those bacteria negatively correlated with the autophagy genes, Bacillus is recognized as a probiotic bacterium in IBD prevention and therapy (Ghavami et al., 2020; Shinde et al., 2020; Zhang et al., 2020; Liu et al., 2021), but the remaining bacteria (e.g., Acidovorax, Acinetobacter, and Stenotrophomonas) were first reported to involve in IBD progression. Interestingly, the co-occurrence network indicated that low-abundance intestinal mucosa-colonizing bacteria were negatively and densely associated with autophagy genes, whereas high-abundance bacteria were positively and sparsely associated with autophagy genes (Supplementary Figure S6).

The autophagy genes in cluster 4, which were highly expressed in CD patients but not in UC patients compared with control population, showed more positive correlations and fewer negative correlations with the intestinal mucosa-colonizing bacteria than found in clusters 6 and 7. The positive and negative correlations were represented by Haemophilus, Ruminococcus torques, and Rhodococcus, and Bacillus, Acidovorax, Acinetobacter, and Stenotrophomonas, respectively. These cluster 4 bacteria-correlated autophagy genes also showed correlation and interaction networks, both of which dominated by Tumor protein P53 (TP53), Stimulator of interferon response CGAMP interactor 1 (STING), and interferon gamma inducible protein 16 (IFI16). These results are highly consistent with previous study suggesting that STING senses bacterial viability to orchestrate autophagy (Moretti et al., 2017). The cluster 4 co-occurrence network also showed more positive correlations between autophagy genes and intestinal mucosa-colonizing bacteria (e.g., Haemophilus, Ruminococcus torques, Enterococcus, Veillonella, Escherichia, and Shigella), but strong negative correlations between Bacillus and Streptococcus and autophagy genes [e.g., Ring finger protein 41 (RNF41), BNIP3P5, and Interleukin 10 (IL10); Supplementary Figure S5]. However, this correlation only existed in CD patients. Although similar correlations between mucosa-colonizing bacteria and autophagy activity have been reported (Nakagawa et al., 2004; Sudhakar et al., 2019; Wu et al., 2019), we used clinical biopsies to explore real crosstalk between mucosal bacteria and autophagy genes during IBD progression and to clarify the complex signaling pathways connecting intestinal microbes and host disease pathogenesis.

Autophagy-Induced ER Stress Is Correlated With Intestinal Mucosa-Colonizing Bacteria in IBD Patients

Moretti et al. (2017) found that STING initiates autophagy after bacterial infection and that this autophagy is mediated by ER stress. ER stress also regulates autophagy processes during IBD progression (Kaser and Blumberg, 2009; Kaser and Blumberg, 2011; Grootjans et al., 2019). According to our analysis, 15 ER stress-related genes in cluster 6 were significantly correlated with intestinal mucosa-colonizing bacteria (Figure 4A; Supplementary Figure S12A), and the main correlations were presented with a co-occurrence network [e.g., negative correlation between Acinetobacter and Heat shock protein family A (Hsp70) member 1A (HSPA1A); Supplementary Figure S10A]. ER stress-related genes in cluster 6, such as HSPA1A, Heat shock protein family A (Hsp70) member 5 (HSPA5), Heat shock protein 90 beta family member 2 (HSP90B2I), Protein disulfide isomerase family A member 3 pseudogene 1 (PDIA3P1), and Phorbol-12-myristate-13-acetate-induced protein 1 (PMAIP1), formed correlation and interaction networks (Supplementary Figures S10B,C). In addition, 31 ER stress-related genes in cluster 7 were significantly correlated with mucosal bacteria (Figure 4D; Supplementary Figure S12B), and formed correlation and interaction networks (Supplementary Figures S10E,F). The main correlations were also presented with a co-occurrence network (e.g., negative correlation of Acidovorax with BCL2 like 11 (BCL2L11) and Wolframin ER transmembrane glycoprotein (WFS1) in cluster 7; Supplementary Figure S10D). The correlation and interaction networks of ER stress-related genes in cluster 7 mainly consisted of Activating transcription factor 6 (ATF6), Homocysteine inducible ER protein with ubiquitin like domain 1 (HERPUD1), BAG cochaperone 6 (BAG6), BCL2 binding component 3 (BBC3), Mesencephalic astrocyte derived neurotrophic factor (MANF), DnaJ heat shock protein family (Hsp40) member B9 (DNAJB9), BCL2L11, Membrane bound transcription factor peptidase, site 1 (MBTPS1), Eukaryotic translation initiation factor 2 subunit alpha (EIF2S1), Valosin containing protein (VCP), Mitogen-activated protein kinase kinase kinase 5 (MAP3K5), and Thioredoxin domain containing 12 (TXNDC12; Supplementary Figures S10E,F). HSPA5-induced ATF6 activity is a typical representative of ER stress signaling in the intestinal epithelium (Stengel et al., 2020) that bridging clusters 6 and 7 profile in this study. Furthermore, combinative analysis of both correlation and interaction network, also exhibited the complexity of ER stress-autophagy signaling cascade in IBD (Figures 4B,C,E,F). In addition, we also found similar ER stress signaling based on KEGG functional pathway analysis (Supplementary Figure S11).

Figure 4

In cluster 4, six ER stress-related genes [i.e., endoplasmic reticulum protein 27 (ERP27), NHL repeat containing E3 ubiquitin protein ligase 1 (NHLRC1), Fc gamma receptor IIb (FCGR2B), Stanniocalcin 2 (STC2), Caspase 17, pseudogene (CASP17P), and ChaC glutathione specific gamma-glutamylcyclotransferase 1 (CHAC1)] were significantly correlated with intestinal mucosa-colonizing bacteria and were highly expressed in CD patients but not in UC patients compared with the control population. The co-occurrence network showed that both STC2 and FCGR2B were negatively and positively correlated with Proteobacteria (Acinetobacter and Klebsiella), respectively. Interestingly, both STC2 and FCGR2B are involved in ER stress activities during IBD progression (Franke et al., 2016; Coope et al., 2019). In addition, co-occurrence network analyses also linked autophagy signaling of STING and IL10 with the ER stress activities of CASP17P and FCGR2B (Supplementary Figure S9).

ER Stress-Related Bile Acid Production Signaling Is Correlated With Mucosa-Colonizing Bacteria in IBD Patients

Increasing evidence indicates that intestinal dysbiosis induces bile acid dysmetabolism and consequent IBD progression (Duboc et al., 2013; Lloyd-Price et al., 2019; Quinn et al., 2020; Sinha et al., 2020). Bacterial-driven bile acid metabolites have also been investigated to describe the pathophysiological basis of bacteria-bile acid associations (Long et al., 2017; Lavelle and Sokol, 2020). Here, we expanded the association between bile acid production-related genes and intestinal mucosa-colonizing bacteria. Eleven genes [e.g., Leptin (LEP), cytochrome P450 family 46 subfamily A member 1 (CYP46A1), KIT proto-oncogene, receptor tyrosine kinase (KIT), Phospholipase A2 group IB (PLA2G1B), and Retinoid X receptor alpha (RXRA)] in cluster 6 and eight genes [e.g., Sulfotransferase family 2A member 1 (SULT2A1), Fibroblast growth factor receptor 4 (FGFR4), and fatty acid binding protein 1 (FABP1)] in cluster 7 were significantly correlated with the intestinal mucosa-colonizing bacteria (Figures 5A,D). Similar to autophagy and ER stress signaling, genes in both cluster 6 and cluster 7 formed correlation and interaction networks, dominated by RXRA (Supplementary Figures S14B,C) and FABP1 (Supplementary Figures S14D,E), respectively. Notably, KIT and PLA2G1B in cluster 6 (Supplementary Figure S14A) and SULT2A1 in cluster 7 (Supplementary Figure S14D) showed strong correlations with the intestinal mucosa-colonizing bacteria. Bile-regulated lipid metabolism can induce ER stress in intestinal epithelial cells and intestinal bacteria modulate this metabolic process (Ko et al., 2020). Here, we identified bile acid-related genes in cluster 6 [RXRA, Aldo-keto reductase family 1 member C2 (AKR1C2), and Oxysterol binding protein like 7 (OSBPL7)] and cluster 7 [FGFR4, Nuclear receptor coactivator 2 (NCOA2), and FABP1], which were all significantly correlated with intestinal mucosa-colonizing bacteria and formed both correlation and interaction networks with ER stress signaling (Figures 5B,C,E,F). KEGG functional analysis also revealed the upregulation of genes related to bile biosynthesis (Supplementary Figure S15) and secretion (Supplementary Figure S16). Furthermore, the expression differentiation of bile acid-producing gene between UC and CD patients increased in cluster 7 than in cluster 6 (Supplementary Figure S17).

Figure 5

In the other hand, the bacteria that correlated with intestinal mucosa bile acid signaling, mainly comprised of Proteobacteria. The correlations between Stenotrophomonas, Acinetobacter, Bacillus, and Corynebacterium and cluster 6 genes (Figure 5A; Supplementary Figure S14A) and Rhodococcus and Acidovorax and cluster 7 genes (Figure 5B; Supplementary Figure S14B), exhibited the most significant correlations with host bile acid signaling. In addition, FGF19 in cluster 4, which is reported to modulate the connection between intestinal microbiota and host inflammation (Gadaleta et al., 2020), was correlated with multiple intestinal bacteria, including Staphylococcus, Pelomonas, Haemophilus, Ruminococcus torques, and Bacteroides (Supplementary Figure S13).

Discussion

Increasing studies focus on the mucosal microbiome and local host immune activities in IBD patients (Howell et al., 2018; Lloyd-Price et al., 2019; Ryan et al., 2020) and the role of autophagy in IBD progression, especially in CD patients (Khor et al., 2011; Larabi et al., 2020). To date, however, no studies have reported on how mucosal bacteria interact with host autophagy activity based on biopsy samples even that obviously hold a great interest of precise diagnosis and therapy of IBD. Here, we used intestinal mucosal biopsies to elucidate the interaction network between mucosal bacteria and host autophagy signaling, as this network can reveal the bacterial infection-related pathogenesis of IBD and highlight potential and precise therapeutic targets for IBD.

Our results showed that active mucosal signaling, including autophagy, could be divided into cluster 4, 6, and 7 based on signaling extensity in UC patients. In general, all genes in the three clusters were highly activated in CD patients, but only genes in cluster 6 were highly expressed in UC patients. In addition, genes in cluster 7 showed increased expression in UC patients compared with that in the control population. The genes in cluster 4 also showed the same activities in UC patients and control population. These results are consistent with previous study showing that autophagy genes are pivotal for intestinal homeostasis and antimicrobial function in IBD, especially CD (Larabi et al., 2020).

In our study, UC patients exhibited greater disturbance in mucosal microbial diversity, community phenotype, and underlying functional spectrum, whereas CD patients exhibited greater autophagy and associated signaling cascades. These differences may have originated from the different pathological durations of UC and CD in patients. The long disease course in CD patients can strongly induce autophagy signaling to fight pathogenic microbial invasion. In another hand, pathogenic microbial invasion can induce CD progression, whereas diffusion behavior tends to increase UC progression (Mirsepasi-Lauridsen et al., 2019). Certainly, mucosal biopsies are ideal samples to determine bacterial diffusion behavior. Previous studies using fecal samples found that CD patients show more extensive microbial community disturbance than UC patients (Halfvarson et al., 2017; Pascal et al., 2017). Based on our findings, we suggest that CD patients continuously repair their gut defenses and establish new gut immune hemostasis during the long-term pathological progression of the disease, which is different from normal states and UC patients.

We also found that most dominant bacteria were positively correlated with autophagy genes and related signaling cascades, whereas less abundant bacteria showed negative correlations. Dominant bacteria may form stronger community structures that can withstand serious environmental shifts or may change into opportunistic pathogens for survival. In contrast, low-abundance bacteria communities may be more sensitive to environmental shift, and therefore tightly related to host metabolism and healthy conditions and finally process a high sensitivity and great value for intestinal disease diagnosis and therapy, including IBD. However, the low-abundance bacteria also brought some uncertainty. For example, Acidovorax (<0.2% of total abundance) was negatively correlated with most autophagy genes and related signaling cascades, even these genes or signaling play an adversarial role in autophagy. In addition, we observed several bacterial species that were not significantly correlated with any host genes (data not shown); it was difficult to explain how they survive under severe immune stress from host and fierce competition from other microbes if they do not clearly show cooperative (or confrontational) properties.

Our study has several limitations. Both the recruited patients (average age > 45) and control population (average age > 35) were older than participants in other studies, which is likely due to the consciousness deficiency of IBD risk in the local young population. In this study, the microbiota-correlated autophagy genes concurrently showed a strong correlation with ER stress. Integrating previous study (Kökten et al., 2018), the autophagy in our study could be either macroautophagy or chaperone-mediated autophagy, but the specific subtype was not clear and needs to be further explored. The relationship between mucosal bacteria and host autophagy activity needs more mucosa biopsies to conform their correlation profiles, but our matched analysis of microbiome and transcriptome from same subjects and location still drew an interactive map of bacteria-autophagy-IBD progression and could provide more valuable diagnostic and therapeutic targets based on autophagy mechanism.

Conclusion

In conclusion, we found that UC patients exhibited more severe dysbiosis and the functional phenotype of intestinal mucosa-colonizing bacteria, whereas CD patients exhibited more active autophagic signaling compared to the controls. Dominant bacteria and low-abundance bacteria showed positive and negative correlations with host autophagy genes, respectively. In addition, correlation and interaction networks were found between bile acid production and ER stress, with the latter showing an interesting interaction with autophagy activity. Thus, the present study elucidated how the intestinal mucosa-colonizing bacterial community interacts with host bile/ER stress/autophagy signaling cascades during IBD progression, which could aid in disease diagnosis and autophagy-targeted therapy.

Funding

This work was supported by the National Natural Science Foundation of China (81860437 and 82160514), the Yunnan Province Innovation Team of Intestinal Microecology-Related Disease Research and Technological Transformation (202005AE160010), Eminent Doctors Program of Yunnan Province (YNWR-MY-2019-072), Yunnan Digestive Endoscopy Clinical Medical Center Foundation [2X2019-01-02]-(2019LCZXKF-XH05, 2020LCZXKF-XH01, 2021LCZXXF-XH01, and 2021LCZXXF-XH15), and Fundamental Research Projects of Yunnan Province (202101AT070275, 202101AY070001-236, and 2018FE001-130).

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

• 1

  Boada-RomeroE.Serramito-GómezI.SacristánM. P.BooneD. L.XavierR. J.Pimentel-MuiñosF. X.et al. (2016). The T300A Crohn’s disease risk polymorphism impairs function of the WD40 domain of ATG16L1. Nat. Commun.7:11821. doi: 10.1038/ncomms11821

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 2

  BourgonjeA. R.FeelischM.FaberK. N.PaschA.DijkstraG.van GoorH. (2020). Oxidative stress and redox-modulating therapeutics in inflammatory bowel disease. Trends Mol. Med.26, 1034–1046. doi: 10.1016/j.molmed.2020.06.006

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 3

  CasiliG.CordaroM.ImpellizzeriD.BruschettaG.PaternitiI.CuzzocreaS.et al. (2016). Dimethyl fumarate reduces inflammatory responses in experimental colitis. J. Crohns Colitis10, 472–483. doi: 10.1093/ecco-jcc/jjv231

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 4

  ChenS.ZhouY.ChenY.GuJ. (2018). fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics34, i884–i890. doi: 10.1093/bioinformatics/bty560

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 5

  CoopeA.PascoalL. B.BotezelliJ. D.da SilvaF. A. R.AyrizonoM. D. L. S.RodriguesB. L.et al. (2019). ER stress activation in the intestinal mucosa but not in mesenteric adipose tissue is associated with inflammation in Crohn’s disease patients. PLoS One14:e0223105. doi: 10.1371/journal.pone.0223105

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 6

  DubocH.RajcaS.RainteauD.BenarousD.MaubertM.-A.QuervainE.et al. (2013). Connecting dysbiosis, bile-acid dysmetabolism and gut inflammation in inflammatory bowel diseases. Gut62, 531–539. doi: 10.1136/gutjnl-2012-302578

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 7

  DykhuizenR. S.MassonJ.McKnightG.MowatA. N.SmithC. C.SmithL. M.et al. (1996). Plasma nitrate concentration in infective gastroenteritis and inflammatory bowel disease. Gut39, 393–395. doi: 10.1136/gut.39.3.393

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 8

  EdgarR. C. (2013). UPARSE: highly accurate OTU sequences from microbial amplicon reads. Nat. Methods10, 996–998. doi: 10.1038/nmeth.2604

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 9

  ErnstJ.Bar-JosephZ. (2006). STEM: a tool for the analysis of short time series gene expression data. BMC Bioinformatics7:191. doi: 10.1186/1471-2105-7-191

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 10

  FrankeL.el BannoudiH.JansenD. T. S. L.KokK.TrynkaG.DiogoD.et al. (2016). Association analysis of copy numbers of FC-gamma receptor genes for rheumatoid arthritis and other immune-mediated phenotypes. Eur. J. Hum. Genet.24, 263–270. doi: 10.1038/ejhg.2015.95

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 11

  GadaletaR. M.Garcia-IrigoyenO.CarielloM.ScialpiN.PeresC.VetranoS.et al. (2020). Fibroblast growth factor 19 modulates intestinal microbiota and inflammation in presence of Farnesoid X receptor. EBioMedicine54:102719. doi: 10.1016/j.ebiom.2020.102719

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 12

  GhavamiS. B.YadegarA.AghdaeiH. A.SorrentinoD.FarmaniM.MirA. S.et al. (2020). Immunomodulation and generation of tolerogenic dendritic cells by probiotic bacteria in patients with inflammatory bowel disease. Int. J. Mol. Sci.21:6266. doi: 10.3390/ijms21176266

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 13

  GomollónF.DignassA.AnneseV.TilgH.Van AsscheG.LindsayJ. O.et al. (2017). 3rd European evidence-based consensus on the diagnosis and management of Crohn’s disease 2016: part 1: diagnosis and medical management. J. Crohns Colitis11, 3–25. doi: 10.1093/ecco-jcc/jjw168

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 14

  GoyalA.YehA.BushB. R.FirekB. A.SieboldL. M.RogersM. B.et al. (2018). Safety, clinical response, and microbiome findings following fecal microbiota transplant in children with inflammatory bowel disease. Inflamm. Bowel Dis.24, 410–421. doi: 10.1093/ibd/izx035

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 15

  GrootjansJ.KrupkaN.HosomiS.MatuteJ. D.HanleyT.SaveljevaS.et al. (2019). Epithelial endoplasmic reticulum stress orchestrates a protective IgA response. Science363, 993–998. doi: 10.1126/science.aat7186

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 16

  HalfvarsonJ.BrislawnC. J.LamendellaR.Vázquez-BaezaY.WaltersW. A.BramerL. M.et al. (2017). Dynamics of the human gut microbiome in inflammatory bowel disease. Nat. Microbiol.2:17004. doi: 10.1038/nmicrobiol.2017.4

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 17

  HowellK. J.KraiczyJ.NayakK. M.GasparettoM.RossA.LeeC.et al. (2018). DNA methylation and transcription patterns in intestinal epithelial cells from pediatric patients with inflammatory bowel diseases differentiate disease subtypes and associate with outcome. Gastroenterology154, 585–598. doi: 10.1053/j.gastro.2017.10.007

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 18

  Inflammatory Bowel Disease Group, Chinese Society of Gastroenterology, Chinese Medical Association (2021). Chinese consensus on diagnosis and treatment in inflammatory bowel disease (2018, Beijing). J. Dig. Dis.22, 298–317. doi: 10.1111/1751-2980.12994

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 19

  JädertC.PhillipsonM.HolmL.LundbergJ. O.BorniquelS. (2014). Preventive and therapeutic effects of nitrite supplementation in experimental inflammatory bowel disease. Redox Biol.2, 73–81. doi: 10.1016/j.redox.2013.12.012

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 20

  KaserA.BlumbergR. S. (2009). Endoplasmic reticulum stress in the intestinal epithelium and inflammatory bowel disease. Semin. Immunol.21, 156–163. doi: 10.1016/j.smim.2009.01.001

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 21

  KaserA.BlumbergR. S. (2011). Autophagy, microbial sensing, endoplasmic reticulum stress, and epithelial function in inflammatory bowel disease. Gastroenterology140, 1738–1747.e2. doi: 10.1053/j.gastro.2011.02.048

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 22

  KhorB.GardetA.XavierR. J. (2011). Genetics and pathogenesis of inflammatory bowel disease. Nature474, 307–317. doi: 10.1038/nature10209

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 23

  KimD.LangmeadB.SalzbergS. L. (2015). HISAT: a fast spliced aligner with low memory requirements. Nat. Methods12, 357–360. doi: 10.1038/nmeth.3317

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 24

  KoC.-W.QuJ.BlackD. D.TsoP. (2020). Regulation of intestinal lipid metabolism: current concepts and relevance to disease. Nat. Rev. Gastroenterol. Hepatol.17, 169–183. doi: 10.1038/s41575-019-0250-7

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 25

  KöktenT.GibotS.LepageP.D’AlessioS.HablotJ.NdiayeN. C.et al. (2018). TREM-1 inhibition restores impaired autophagy activity and reduces colitis in mice. J. Crohns Colitis12, 230–244. doi: 10.1093/ecco-jcc/jjx129

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 26

  LarabiA.BarnichN.NguyenH. T. T. (2020). New insights into the interplay between autophagy, gut microbiota and inflammatory responses in IBD. Autophagy16, 38–51. doi: 10.1080/15548627.2019.1635384

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 27

  LavelleA.SokolH. (2020). Gut microbiota-derived metabolites as key actors in inflammatory bowel disease. Nat. Rev. Gastroenterol. Hepatol.17, 223–237. doi: 10.1038/s41575-019-0258-z

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 28

  LeeJ. M.WagnerM.XiaoR.KimK. H.FengD.LazarM. A.et al. (2014). Nutrient-sensing nuclear receptors coordinate autophagy. Nature516, 112–115. doi: 10.1038/nature13961

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 29

  LiB.DeweyC. N. (2011). RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome. BMC Bioinformatics12:323. doi: 10.1186/1471-2105-12-323

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 30

  LiuY.YinF.HuangL.TengH.ShenT.QinH. (2021). Long-term and continuous administration of during remission effectively maintains the remission of inflammatory bowel disease by protecting intestinal integrity, regulating epithelial proliferation, and reshaping microbial structure and function. Food Funct.12, 2201–2210. doi: 10.1039/d0fo02786c

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 31

  Lloyd-PriceJ.ArzeC.AnanthakrishnanA. N.SchirmerM.Avila-PachecoJ.PoonT. W.et al. (2019). Multi-omics of the gut microbial ecosystem in inflammatory bowel diseases. Nature569, 655–662. doi: 10.1038/s41586-019-1237-9

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 32

  LongS. L.GahanC. G. M.JoyceS. A. (2017). Interactions between gut bacteria and bile in health and disease. Mol. Asp. Med.56, 54–65. doi: 10.1016/j.mam.2017.06.002

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 33

  LoucaS.ParfreyL. W.DoebeliM. (2016). Decoupling function and taxonomy in the global ocean microbiome. Science353, 1272–1277. doi: 10.1126/science.aaf4507

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 34

  LoveM. I.HuberW.AndersS. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol.15:550. doi: 10.1186/s13059-014-0550-8

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 35

  MagočT.SalzbergS. L. (2011). FLASH: fast length adjustment of short reads to improve genome assemblies. Bioinformatics27, 2957–2963. doi: 10.1093/bioinformatics/btr507

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 36

  MagroF.GionchettiP.EliakimR.ArdizzoneS.ArmuzziA.Barreiro-de AcostaM.et al. (2017). Third European evidence-based consensus on diagnosis and management of Ulcerative Colitis. Part 1: definitions, diagnosis, extra-intestinal manifestations, pregnancy, cancer surveillance, surgery, and Ileo-anal pouch disorders. J. Crohns Colitis11, 649–670. doi: 10.1093/ecco-jcc/jjx008

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 37

  ManichanhC.Rigottier-GoisL.BonnaudE.GlouxK.PelletierE.FrangeulL.et al. (2006). Reduced diversity of faecal microbiota in Crohn’s disease revealed by a metagenomic approach. Gut55, 205–211. doi: 10.1136/gut.2005.073817

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 38

  Matsuzawa-IshimotoY.ShonoY.GomezL. E.Hubbard-LuceyV. M.CammerM.NeilJ.et al. (2017). Autophagy protein ATG16L1 prevents necroptosis in the intestinal epithelium. J. Exp. Med.214, 3687–3705. doi: 10.1084/jem.20170558

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 39

  Mirsepasi-LauridsenH. C.VallanceB. A.KrogfeltK. A.PetersenA. M. (2019). Pathobionts associated with inflammatory bowel disease. Clin. Microbiol. Rev.32, e00060–e00018. doi: 10.1128/CMR.00060-18

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 40

  MorettiJ.RoyS.BozecD.MartinezJ.ChapmanJ. R.UeberheideB.et al. (2017). STING senses microbial viability to orchestrate stress-mediated autophagy of the endoplasmic reticulum. Cell171, 809–823.e13. doi: 10.1016/j.cell.2017.09.034

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 41

  MorganX. C.KabakchievB.WaldronL.TylerA. D.TickleT. L.MilgromR.et al. (2015). Associations between host gene expression, the mucosal microbiome, and clinical outcome in the pelvic pouch of patients with inflammatory bowel disease. Genome Biol.16:67. doi: 10.1186/s13059-015-0637-x

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 42

  NakagawaI.AmanoA.MizushimaN.YamamotoA.YamaguchiH.KamimotoT.et al. (2004). Autophagy defends cells against invading group A Streptococcus. Science306, 1037–1040. doi: 10.1126/science.1103966

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 43

  NguyenH. T. T.LapaquetteP.BringerM.-A.Darfeuille-MichaudA. (2013). Autophagy and Crohn’s disease. J. Innate Immun.5, 434–443. doi: 10.1159/000345129

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 44

  NighotP. K.HuC.-A. A.MaT. Y. (2015). Autophagy enhances intestinal epithelial tight junction barrier function by targeting claudin-2 protein degradation. J. Biol. Chem.290, 7234–7246. doi: 10.1074/jbc.M114.597492

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 45

  OttS. J.MusfeldtM.WenderothD. F.HampeJ.BrantO.FölschU. R.et al. (2004). Reduction in diversity of the colonic mucosa associated bacterial microflora in patients with active inflammatory bowel disease. Gut53, 685–693. doi: 10.1136/gut.2003.025403

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 46

  PascalV.PozueloM.BorruelN.CasellasF.CamposD.SantiagoA.et al. (2017). A microbial signature for Crohn’s disease. Gut66, 813–822. doi: 10.1136/gutjnl-2016-313235

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 47

  PerteaM.PerteaG. M.AntonescuC. M.ChangT. C.MendellJ. T.SalzbergS. L. (2015). StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat. Biotechnol.33, 290–295. doi: 10.1038/nbt.3122

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 48

  PetersL. A.PerrigoueJ.MorthaA.IugaA.SongW.-M.NeimanE. M.et al. (2017). A functional genomics predictive network model identifies regulators of inflammatory bowel disease. Nat. Genet.49, 1437–1449. doi: 10.1038/ng.3947

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 49

  PlichtaD. R.GrahamD. B.SubramanianS.XavierR. J. (2019). Therapeutic opportunities in inflammatory bowel disease: mechanistic dissection of host-microbiome relationships. Cell178, 1041–1056. doi: 10.1016/j.cell.2019.07.045

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 50

  QuinnR. A.MelnikA. V.VrbanacA.FuT.PatrasK. A.ChristyM. P.et al. (2020). Global chemical effects of the microbiome include new bile-acid conjugations. Nature579, 123–129. doi: 10.1038/s41586-020-2047-9

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 51

  R Core Team (2013). R: A Language and Environment for Statistical Computing.

  • Google Scholar

• 52

  RobinsonM. D.McCarthyD. J.SmythG. K. (2010). edgeR: a bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics26, 139–140. doi: 10.1093/bioinformatics/btp616

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 53

  RyanF. J.AhernA. M.FitzgeraldR. S.Laserna-MendietaE. J.PowerE. M.ClooneyA. G.et al. (2020). Colonic microbiota is associated with inflammation and host epigenomic alterations in inflammatory bowel disease. Nat. Commun.11:1512. doi: 10.1038/s41467-020-15342-5

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 54

  SaijoF.MilsomA. B.BryanN. S.BauerS. M.VowinkelT.IvanovicM.et al. (2010). On the dynamics of nitrite, nitrate and other biomarkers of nitric oxide production in inflammatory bowel disease. Nitric Oxide22, 155–167. doi: 10.1016/j.niox.2009.11.009

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 55

  ShenS.ParkJ. W.LuZ. X.LinL.HenryM. D.WuY. N.et al. (2014). rMATS: robust and flexible detection of differential alternative splicing from replicate RNA-Seq data. Proc. Natl. Acad. Sci. U. S. A.111, E5593–E5601. doi: 10.1073/pnas.1419161111

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 56

  ShindeT.PereraA. P.VemuriR.GondaliaS. V.BealeD. J.KarpeA. V.et al. (2020). Synbiotic supplementation with prebiotic green banana resistant starch and probiotic Bacillus coagulans spores ameliorates gut inflammation in mouse model of inflammatory bowel diseases. Eur. J. Nutr.59, 3669–3689. doi: 10.1007/s00394-020-02200-9

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 57

  SinhaS. R.HaileselassieY.NguyenL. P.TropiniC.WangM.BeckerL. S.et al. (2020). Dysbiosis-induced secondary bile acid deficiency promotes intestinal inflammation. Cell Host Microbe27, 659–670.e5. doi: 10.1016/j.chom.2020.01.021

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 58

  SoltysK.StuchlikovaM.HlavatyT.GaalovaB.BudisJ.GazdaricaJ.et al. (2020). Seasonal changes of circulating 25-hydroxyvitamin D correlate with the lower gut microbiome composition in inflammatory bowel disease patients. Sci. Rep.10:6024. doi: 10.1038/s41598-020-62811-4

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 59

  StackebrandtE.GoebelB. M. (1994). Taxonomic note: a place for DNA-DNA reassociation and 16S rRNA sequence analysis in the present species definition in bacteriology. Int. J. Syst. Evol. Microbiol.44, 846–849. doi: 10.1099/00207713-44-4-846

  • CrossRef
  • Google Scholar

• 60

  StengelS. T.FazioA.LipinskiS.JahnM. T.AdenK.ItoG.et al. (2020). Activating transcription factor 6 mediates inflammatory signals in intestinal epithelial cells upon endoplasmic reticulum stress. Gastroenterology159, 1357–1374.e10. doi: 10.1053/j.gastro.2020.06.088

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 61

  SudhakarP.JacominA.-C.HautefortI.SamavedamS.FatemianK.AriE.et al. (2019). Targeted interplay between bacterial pathogens and host autophagy. Autophagy15, 1620–1633. doi: 10.1080/15548627.2019.1590519

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 62

  WangQ.GarrityG. M.TiedjeJ. M.ColeJ. R. (2007). Naive Bayesian classifier for rapid assignment of rRNA sequences into the new bacterial taxonomy. Appl. Environ. Microbiol.73, 5261–5267. doi: 10.1128/aem.00062-07

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 63

  WardT.LarsonJ.MeulemansJ.HillmannB.LynchJ.SidiropoulosD.et al. (2017). BugBase predicts organism-level microbiome phenotypes. biorxiv [Preprint]. 133462. doi: 10.1101/133462

  • CrossRef
  • Google Scholar

• 64

  WuY.WangB.XuH.TangL.LiY.GongL.et al. (2019). Probiotic attenuates oxidative stress- induced intestinal injury via p38-mediated autophagy. Front. Microbiol.10:2185. doi: 10.3389/fmicb.2019.02185

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 65

  XieC.MaoX.HuangJ.DingY.WuJ.DongS.et al. (2011). KOBAS 2.0: a web server for annotation and identification of enriched pathways and diseases. Nucleic Acids Res.39, W316–W322. doi: 10.1093/nar/gkr483

  • Pubmed Abstract
  • CrossRef
  • Google Scholar

• 66

  ZhangX.TongY.LyuX.WangJ.WangY.YangR. (2020). Prevention and alleviation of dextran sulfate sodium salt-induced inflammatory bowel disease in mice with -fermented milk inhibition of the inflammatory responses and regulation of the intestinal flora. Front. Microbiol.11:622354. doi: 10.3389/fmicb.2020.622354

  • Pubmed Abstract
  • CrossRef
  • Google Scholar