Mining the Human Gut Microbiota for Immunomodulatory Organisms

Abstract

Within the human gut reside diverse microbes coexisting with the host in a mutually advantageous relationship. Evidence has revealed the pivotal role of the gut microbiota in shaping the immune system. To date, only a few of these microbes have been shown to modulate specific immune parameters. Herein, we broadly identify the immunomodulatory effects of phylogenetically diverse human gut microbes. We monocolonized mice with each of 53 individual bacterial species and systematically analyzed host immunologic adaptation to colonization. Most microbes exerted several specialized, complementary, and redundant transcriptional and immunomodulatory effects. Surprisingly, these were independent of microbial phylogeny. Microbial diversity in the gut ensures robustness of the microbiota's ability to generate a consistent immunomodulatory impact, serving as a highly important epigenetic system. This study provides a foundation for investigation of gut microbiota-host mutualism, highlighting key players that could identify important therapeutics.

Keywords: gnotobiotic; gut bacteria; immunomodulation; immunoprofiling; innate and adaptive immunity; microbiome.

Copyright © 2017 Elsevier Inc. All rights reserved.

Figures

Figure 1:. Experiment design and bacterial colonization

(A) Four week-old GF mice were monocolonized with human gut bacteria and analyzed after 2 weeks for colonization, impact on immune system and genomic activity. (B) Innate and adaptive immune responses were analyzed by flow cytometry of cells extracted from SI, PPs, colons, mLNs, and SLOs. (C) Cladogram of the human gut microbiota. Microbes were identified in the Human Microbiome Project (HMP) database except for SFB. Blue diamonds denote the genera included; red stars mark the species. Species where more than one strain was analyzed are in bold type. The outer ring represents a bar graph of the prevalence of each genus. (D) Average CFU per gram of feces, and their standard deviations. (E) Bar graphs of CFUs in mLNs (per organ, top) and SLO (bottom). See Table S1 and Fig. S1.

Figure 2:. Immunomodulation by gut microbes

(A) Average frequencies of each immunocyte population for every microbe. For cell type frequency determination (y-axis) and microbe identification (x-axis) see Tables S1A,B, S2A,B, and Fig. S1. (B) Heatmap of average fold changes (relative to GF) for cells in the colon and SI following monocolonization, and fecal IgA. Gray- no data. (C) Proportion of colonic immune cell types (compared to GF) with a z-score ≥ 2. (D) Example of colonization influencing the gating configuration but not frequency of cell populations. Flow cytometry plots shown are for CD11b+CD11c+ MNPs and DCs. (E) Cytokine responses in SI and colon. See Fig. S2, Table S2

Figure 3:. Local and systemic immunologic correlations.

(A) Clustered heatmap of Pearson correlation coefficients (r) for immunophenotypes after moncolonization. (B) Average cell frequency correlations: SLO versus colon, for MFs(upper), Tregs (middle) and MNPs (lower). (C) Average cell frequency correlations: SLO versus colon, for monocytes. (D) Hierarchical clustering dendrogram of bacteria based on the Pearson correlation of their overall immunologic impact on the SI and colon. Values were normalized to the mean across all microbes. See also Figure S3.

Figure 4:. Transcriptional responses to colonization

(A) Mean coefficient of variation (CV) in transcripts from the colons of monocolonized mice and GF mice. (B-C). Heatmap of fold changes of transcripts differentially expressed in (B) the colon and (C) SI of monocolonized and SPF mice compared to GF mice. Bacteria (columns) are clustered by hierarchical clustering; Genes (rows) are clustered by K-means clustering. Association of these transcripts with particular immune and non-immune cell types was verified in gene expression databases such as ImmGen and GNF. Enriched pathways were identified using GO. See also Fig. S4, Table S3

Figure 5:. Colonic plasmacytoid dendritic cells are most prolific myeloid responders to the gut microbiota

(A) Representative flow cytometry dot plots of a pDC ‘low inducer’, Propionibacterium granulosum (Pgran.A042) and a ‘high inducer’ B. vulgatus (Bvulg.ATCC8482). Cells were gated as CD45+CD19-CD11b-. (B) Frequencies of pDCs in the colon induced by monocolonization (C) Pearson correlation between pDC’s in SI vs. colon (p=0.0006). (D) Pearson correlation between colonic pDCs and Tregs (p=0.003). (E-F) Correlation coefficients were calculated between the expression value of each gene from the whole tissue transcriptome (SI, and colon) and the proportions of pDCs for each monocolonizing microbe (SI and colon). (E) Genes related to the interferon signature are marked in red. (F) Genes having similar expression patterns and correlating best in both SI and colon are highlighted in green. (Right) Bar graph of the enriched biological pathways of these highly correlating genes as analyzed by Enrichr. Most significant pathways determined by GO Molecular Function (p

Figure 6:. Antimicrobial peptides exhibit divergent patterns of expression in the SI and colon

(A) Coefficient of variation (CV) vs. mean expression in GF mice for all genes in the SI (left panel) and colon (right panel). Only genes expressed above background level are shown. Antimicrobial peptides (AMPs) are highlighted and color-coded according to the categories listed. (B) The CV of all expressed genes in the colons of GF vs monocolonized mice, as shown in Fig. 4A, here AMP genes highlighted. (C-D) Heatmaps illustrating the differential expression of AMPs in the SI (C) and colon (D) in various monocolonized mice compared to GF. Heatmap colors represent the log2 fold change values relative to GF. Only AMPs expressed above background levels are shown. (E) Gene programs correlated with AMP expression in the colon. For every gene expressed in the colon, its correlation with colonic AMP genes (Reg3 family and α-defensins) is plotted for GF vs. monocolonized mice (left panel). Top correlated genes (Spearman rho>0.6) are highlighted in black and parsed for enrichment of biological pathways using Enrichr. Top pathways from GO Molecular Function, with corresponding adjusted p.values and gene names, are shown (right panel). See also Table S5.

Figure 7:. Host response to Fusobacterium varium

(A) Amplified gene expression preferential to F. varium (Fvari.AO16), based on the conservative gene list established in Fig. 4B–C. Fold change (FC) of Fvari.AO16 over GF (y-axis) was compared to the maximum induced FC by any other microbe over GF (x-axis). Top – SI, bottom- colon. (B) Functional analysis of genes suppressed by F. varium. STRING-db clustering and functional categories of significantly altered genes (FC≤0.5 in SI; FC≤0.67 in colon vs. GF; FDR 0.1). Blue dots - genes from (A) preferentially suppressed by Fvari.AO16; gray dots - all other suppressed genes in the Fvari.AO16 response that formed connected clusters. Functional categories determined by GO and KEGG are shown: “Retinol metabolism” FDR 2.25e-15. “Bile acid metabolism” FDR 2.6e-7. “Immune response” FDR 0.0138. (C) Functional analysis of genes induced by F. varium. STRING-db clustering and functional categories of significantly altered genes (SI FC≥2, colon FC≥1.5 vs. GF; FDR 0.1). Red dots – genes from (A) preferentially induced by Fvari.AO16-; gray dots - all other induced genes in Fvari.AO16 response that formed connected clusters. Functional categories determined by GO and KEGG: “Regulation of TRP channels” FDR 0.00313; “AA metabolism” FDR 0.0241; “Globin” FDR 3.78e-8; “Triglyceride metabolism” FDR 0.0184; “Glycerolipid metabolism” FDR 1.32e-7. (D) Representative flow cytometry plots of CD4 and CD8 expression in GF and Fvari.AO16, gated on CD45+ CD19-TCRβ+ cells. (E) Frequencies of T4, T8, and DN T cells normalized to the mean frequency of all microbes in all monocolonizations. See also Tables S6A–B.