Effect of macrolide courses on the intestinal microbiota

To determine if a single antibiotic course introduced early in life is sufficient to lead to durable changes in both the microbiota and in the host, we compared the effect of exposing mice to a single pulsed antibiotic course (PAT1) vs a 3-course (PAT3) regimen (Fig. 1a ). The first course was administered to both PAT1 and PAT3 pups while nursing at postnatal day 5 (P5) for 5 days. The PAT3 group received two additional courses at P27 and P36 for 3 days (Fig. 1a). Since the first antibiotic course was given to the dams during nursing, we sampled both pups and their mothers (dams) to compare the effects of the course on developing (pups; P5–10) and mature (dams; ~12 weeks old) microbiota.

Fig. 1 Effect of number and timing of antibiotic doses on intestinal microbial communities. a Study design: 5-day-old C57BL/6 pups were treated with one course of tylosin at P5 for 5 days (PAT1 group) through their mother’s milk, or with two additional doses at P27 and P36 for 3 days each (PAT3 group). Twelve-week-old dams were treated with one course at pup P5. Sample sizes for dams were n = 7 (PAT), n = 3 (control). Sample sizes for their pups were n = 12 (control), n = 17 (PAT1), n = 19 (PAT3). b Mean (±SEM) unrarified α-diversity using the phylogenetic diversity (PD) metric in fecal samples of dams and offspring (pups) over the course of the experiment. Solid lines represent female pups n = 3 (control), n = 8 (PAT1), n = 9 (PAT3) and dams; dotted lines are male pups; n = 9 (control), n = 9 (PAT1), n = 10 (PAT3). Statistical analysis performed using two-sample t-test with Monte Carlo permutations, for statistical significance, see Supplementary Table 1. c Unweighted UniFrac analysis of fecal specimens of pups and dams visualized by principal coordinate analysis (PCoA). The three components explain 42.3% of the total variance. d Intergroup unweighted UniFrac distances averaged over independently drawn sample pairs (subsampled without replacement and replicated 999 times) for the time points and groups shown in c, shading indicates initial period of antibiotic exposure (for statistical significance, see Supplementary Table 2). e Mean relative abundance of taxa in control, PAT1 and PAT3 groups over the course of the experiment in fecal (over time), cecal and ileal samples Full size image

Exposure to PAT1 and PAT3 regimens early in life altered intestinal microbial community composition and dynamics. After PAT1, both phylogenetic diversity and community evenness were decreased, persisting for > 6 weeks, with more extensive changes in the PAT3 mice (Fig. 1b, Supplementary Fig. 1a and Supplementary Table 1 ). After weaning, the microbial communities in the PAT1 and PAT3 groups are similar until the additional antibiotic courses in the PAT3 group (Fig. 1c). However, 2 weeks after the final antibiotic pulse, the PAT3 microbial communities reverted to the PAT1 community state (Fig. 1c). Inter-group community variation was higher in both the PAT1 and PAT3 groups compared to control, indicating that the greatest effect was from the first antibiotic course (Fig. 1d, Supplementary Fig. 1b and Supplementary Table 2). At P52, 13 days after the final PAT3 exposure, bacterial load was decreased compared to control, but the PAT1 group had recovered (Supplementary Fig. 1c). Using either one or three antibiotic courses decreased relative abundances in the Bacteroides family S24-7, genus Bifidobacterium, and segmented filamentous bacteria (SFB), and increased relative abundances in Enterobacteriaceae and Akkermansia (Fig. 1e and Supplementary Fig. 1d–e). In contrast to pups, adult females, ~12 weeks of age, exposed to one antibiotic pulse showed intermediate effects on microbial phylogenetic diversity and community composition (Fig. 1b, c). The adult mice showed improved microbial community recovery compared to the antibiotic-exposed pups, shown by a decrease in inter-group UniFrac distance in PAT vs control (Fig. 1d). As in the pups, by sacrifice there was no effect on the total bacterial load 50 days after the PAT1 exposure in the dams (Supplementary Fig. 1c). Taken together, these findings indicate that a single antibiotic course, early in the mouse life, is sufficient to lead to long-term alterations on intestinal microbial communities, greater than in adult mice.

Effects of macrolide courses on host immune features

Since our prior studies of PAT showed substantial changes in gene expression in the ileum16, consistent with predictions based on the role of the microbiota in conventional vs germ-free animals17, we first assessed the effect of the number of antibiotic courses on ileal gene expression. Unsupervised hierarchical clustering of differential ileal genes revealed strong antibiotic exposure signatures, with all of the PAT mice segregating from the controls (Fig. 2a and Supplementary Data 1). Compared with controls at P52, after FDR correction, 148 genes were differentially expressed in either of the two PAT groups; 137 (93%) were in the PAT3 mice, with 63 shared with PAT1 and 10 unique to the PAT1 mice (Fig. 2a). Thus, although the greatest effects occurred in the PAT3 mice, even a single PAT course affected gene expression, continuing > 40 days after the last exposure.

Fig. 2 Administration of one or three macrolide courses alters ileal gene signatures and T-cell populations. a Heatmap with unsupervised clustering of FDR-corrected, differential ileal gene expression profiles using Nanostring technology (Control (n = 7), PAT1 (n = 8) and PAT3 (n = 9)). After one-way ANOVA, Tukey HSD multiple comparison testing and FDR-correction, 148 genes remained significantly different between the groups. b Mean (±SD) frequency of small intestine lamina propria lymphocytes after PAT1, PAT3 or control exposure. Cells were isolated at sacrifice from pups (mean P52) and their respective dams (P61). Populations were gated on live CD45+CD4+ cells, and representative frequencies of CD4+IL17A+ cells shown. c Frequency of splenic lymphocytes after PAT1, PAT3 or control exposure. Cells were isolated from pups at sacrifice, gated on CD4+ and CD8+ populations, and representative frequencies shown. Data compiled from three experiments (n = 9–19 mice/group). d Secretory IgA (μg/ml) (mean ± SD) in fecal samples of control and PAT-exposed dams at P27 and their pups at P27 and P40. e Serum IgA (μg/ml) (mean ± SD) in control and PAT-exposed mice in pup and dam samples at mean P52 or P61, respectively. For panels (b–e) significance testing performed using the Mann–Whitney non-parametric test or the Kruskal–Wallis non-parametric test used with Dunnett’s multiple comparisons test. f Volcano plots showing significant global ileal transcriptomic alterations in pups and dams, determined by RNAseq; n > 16,000 genes examined, control (n = 3) and PAT (n = 6) dams, respectively; control (n = 3) or PAT (n = 4) pups. g Differentially expressed genes in relation to PAT exposure in P52 pups, and their dams at P61, determined by RNAseq, and depicted by direction of differences. For all panels: *p < 0.05, **p < 0.01, ***p < 0.001, ****p < 0.0001 Full size image

Since a single PAT treatment showed significant decreases in genes involved in immunity, we asked whether the antibiotic exposure altered local and peripheral T-helper populations. Total numbers of small intestinal lamina propria (SI-LP) and splenic leukocytes were not significantly different between the control, PAT1 and PAT3 pups (p > 0.05, Kruskal–Wallis non-parametric test used with Dunnett’s multiple comparisons test) (Supplementary Fig. 1g). However, intestinal CD4+ IL17A+ lymphocytes were decreased in both pup groups, with no differences in the dams (Fig. 2b), Th1 (CD4+ IFNγ+)-expressing cells were not affected but in the PAT3 group, Foxp3+ regulatory T cells were increased (Supplementary Fig. 1h). Splenic CD8+ lymphocytes were reduced in both pup groups, while CD4+ cells were decreased after PAT3 exposure (Fig. 2c).

Due to known roles of the microbiome and Th17 lymphocytes on IgA expression, we next asked whether the single antibiotic course affected mucosal and systemic IgA levels. In the pups, at P27, after both PAT groups had received only a single course, fecal sIgA levels were lower than in the control mice (p < 0.001; Kruskal–Wallis non-parametric test with Dunnett’s multiple comparisons test) (Fig. 2d). By P40, sIgA levels rose in the controls as expected, reflecting increased PIGR-mediated active transport18, but was lower in the PAT1 group, and more so in the PAT3 mice; however, this phenomenon was absent in the PAT-dams (Fig. 2d). We also assessed fecal sIgA levels in the control and PAT dams before treatment to determine if there was variability in sIgA expression that could have mitigated the PAT-induced sIgA effects in the dams and her subsequent offspring. Fecal samples collected before breeding and during gestation showed no differences in fecal sIgA in dams that later were randomized into control and PAT groups (Supplementary Table 3). At sacrifice, IgA levels in serum were decreased in both the PAT1 and PAT3 pups, but not in their PAT1-exposed dams (Fig. 2e).

Since the dams of the PAT1 pups were exposed at the same time to the single antibiotic pulse but were adults rather than neonates (Fig. 1a), we examined whether host age at antibiotic exposure affected gene expression. Using RNAseq to compare global ileal gene expression between control and PAT-exposed animals, we found fewer differentially expressed genes in the dams than in the parallel comparisons of the pups (Fig. 2f, g). In total, only 47 of the 103 genes that were different between PAT and control mice were shared by the PAT-exposed dams and pups. Using Ingenuity Pathway Analysis, we examined the differentiating pathways that were shared between pups and dams. Of the 12 most significant gene expression pathways differentiating the PAT and control pups, most affected immune processes and were down-regulated. However, in the PAT dams, essentially all significant pathways were up-regulated genes and related to metabolic processes (Supplementary Table 4).

In total, a single antibiotic regimen was sufficient to alter mucosal and systemic immunological markers in the pups; comparisons with their similarly exposed dams provide evidence that the host effects depended on the age at exposure.

Microbiota presence is necessary for PAT effects

In addition to their anti-bacterial activity, macrolides are also known for immunomodulatory properties19, 20. To determine if the observed immunologic effects were due to the direct effects of the antibiotic on the host, or due to antibiotic-induced changes in the microbiome, we exposed germ-free (GF) mice with the same P5-10 tylosin pulse (PAT) or not (control), and followed mice until sacrifice at P50 (Fig. 3a); in parallel, we exposed SPF mice to PAT or not. As we now expected in the PAT-exposed SPF mice, microbial diversity, community structure, and composition were altered throughout the experimental course (Fig. 3b, c ).

Fig. 3 Physiologic effects of PAT in germ-free and conventional specific pathogen-free mice. a Study design: Germ-free (GF) and Specific pathogen-free (SPF) C57BL/6 mice were treated with one course of tylosin on P5 for 5 days; SPF mice (control (n = 15), PAT (n = 13)), GF mice (control (n = 7)), (PAT (n = 6)). All mice were sacrificed at P50. b Mean (±SD) unrarified α-diversity using the phylogenetic diversity (PD) metric of fecal samples at P21 and cecal and luminal contents at P50 in control and PAT SPF mice. Statistical analysis performed using two-sample t-test with Monte Carlo permutations, for statistical significance, see Supplementary Table 1. c Unweighted UniFrac analysis of fecal specimens in male (darker circles) and female (lighter circles) PAT or control SPF mice, depicted by principal coordinate analysis (PCoA). The three components explain 58% in total variance. Statistical analysis performed using an Adonis test, p = 0.001. d Aggregate cecal size in SPF and GF mice at P50 (mean (±SD)). e Secretory IgA (mean (±SD)) in luminal contents of offspring at P50, SPF mice (control (n = 5), PAT (n = 6), GF mice (control (n = 3), PAT (n = 4)). f PCoA of differential gene expression profiles in GF and SPF Control and PAT-exposed groups, with 76.2 and 8.1 % of the variance represented in PC1 and PC2, respectively. g Heat map of unsupervised hierarchal clustering of differential ileal gene expression profiles. After one-way ANOVA, Tukey HSD multiple comparison testing and FDR-correction, 145 genes were differentially expressed either between control vs PAT SPF, or GF vs SPF control groups (n = 3/group). h Frequency of small intestine lamina propria lymphocyte (LPL) populations in PAT-treated GF and SPF mice (mean (±SD)). Intestinal LPLs were isolated from offspring at sacrifice. Representative flow cytometry plots were gated on live TCRβ+ cells, with frequencies of TCRβ+, IL17A+CD4+, and Rorγt+CD4+ cells shown. Data are compiled from two experiments. Statistical analysis performed using the Kruskal–Wallis non-parametric test with Dunnett’s multiple comparisons test. For all panels: *p < 0.05, **p < 0.01, ****p < 0.0001 Full size image

Cecal size was increased in GF mice, as expected, and trended larger in the PAT-exposed SPF mice (Fig. 3d). At sacrifice, cecal sIgA levels in the GF-control mice were lower than the SPF-control group, as expected, confirming the important role of the microbiota in priming the sIgA response21; sIgA levels in the PAT-exposed and unexposed GF mice were not different (Fig. 3e ), indicating lack of direct antibiotic effect.

To assess potential PAT effects on immunologic development, we examined ileal gene expression in the four study groups at P50. Between all groups, 145 genes examined were different. Between the (antibiotic-free) control SPF and control GF mice, 141 (26%) of the 547 genes were different, reflecting the expected profound influence of the GF state on immunological development4. In aggregate, the gene expression profile of the antibiotic-unexposed SPF mice (normal) clustered far apart from the SPF-PAT-exposed, and the germ-free mice (Fig. 3f). By unsupervised hierarchical clustering, all GF mice, whether or not PAT-exposed, formed a single cluster (Fig. 3g and Supplementary Data 2 ), indicating lack of macrolide effect in mice without microbiota. In contrast, SPF mice showed the expected strong PAT effect compared to the unexposed controls; importantly, the profiles for PAT-exposed SPF mice clustered with the GF mice, and not with the SPF controls (Fig. 3g), even 40 days after the antibiotic exposure, indicating a persistent effect. Between the control and PAT-exposed SPF mice, multiple (n = 55) genes had differential expression, but in the GF mice, there were no genes that were differentially expressed between PAT and control. Thus, gene expression profiling supported the substantial broad immunological differences between the GF and SPF states, the similarity of PAT-exposed SPF to GF even > 40 days after the antibiotic exposure, and the lack of PAT effect in the GF state.

We examined small intestinal lamina propria T-cell populations in the four groups of mice. As we now expected in the PAT-exposed compared to the (unexposed) SPF controls, TCRβ+, TCRβ+ CD4+ Rorγt+, and TCRβ+ CD4+ IL17A+ SI-LP lymphocytes were decreased (Fig. 3h). In the GF mice, all of these cell populations were depleted, without significant differences between PAT and untreated controls. Moreover, ileal immune gene expression, sIgA expression, and intestinal immune populations were substantially low in GF mice and were not further affected by PAT treatment. This work provides evidence that the immunological effects caused by early-life PAT in mice are not a direct antibiotic effect but an effect of an antibiotic-altered microbiota.

A PAT-altered microbiota is sufficient to impair immune features

Since a PAT-perturbed microbiota is necessary to impair host immunity, we next asked whether or not it is sufficient for the immunologic effects observed. We focused on changes in fecal sIgA, with the strong evidence from the prior experiments, and because we could obtain serial data from individual mice. To accomplish this, first we exposed mice to PAT or not (control) to serve as donors for microbiota transfer (Fig. 4a). On P12, we harvested cecal contents from the PAT and control mice as two separate pools, and gavaged these donor materials into P35 GF-recipient mice, conventionalizing them. For the first 6 days after transfer, all recipients had low sIgA levels (Fig. 4b ). Fecal sIgA levels started rising from post-transfer day 12 and leveled off after day 25. However, the area under the curve was significantly smaller in the recipients of the PAT-perturbed microbiota (random effects ANOVA, F 1,10 = 43.36, p < 0.0001) (Fig. 4b). At sacrifice (day 77 after transfer), the recipients of the PAT-perturbed microbiota had lower splenic TCRβ+ populations than control recipients (Fig. 4c). The differential immunologic phenotypes were independent of microbial density (Supplementary Fig. 2a ), since 16S levels were similar for both recipient groups. Importantly, SFB were negligible in all samples until day 77 after transfer, so the early differential effects between the PAT and control microbiota recipients were SFB-independent (Supplementary Fig. 2b ). This study provides evidence that a PAT-perturbed microbiota (regardless of SFB status) is sufficient to impair host immunological development in recipient mice, involving both humoral and cell-mediated immune characteristics.

Fig. 4 Dynamics of PAT-perturbed microbiota after transfer and effects on host immune phenotypes. a Study design: donor mice received non-acidified water alone (controls) or drinking water with tylosin (PAT) from P5 to P10 and were sacrificed at P12, then cecal contents were harvested and pooled. Germ-free (GF) C57BL/6 mice at P35 were gavaged with cecal contents from the P12 control or PAT-exposed donors (n = 7 per group). b Fecal secretory IgA from day 0 to 76 post-gavage in the now-conventionalized mice that had received control or PAT cecal contents. A random effects repeated measures analysis of variance (ANOVA) was used to test for differences in IgA values between the PAT and the control recipients, taking into account the correlated structure of the measurements within subjects, (F 1,10 = 43.36, p > 0.0001). c Frequency of splenic TCRβ+ lymphocytes at sacrifice, 77 days post-gavage, (mean (±SD)). Statistical analysis was performed using the Mann–Whitney U non-parametric test; *p < 0.05. d Relative abundance of taxa in inoculum (in), fecal, cecal (ce) and colon (co) samples over the course of the experiment. e Median (±IQR) unrarified microbiota α-diversity over the course of the experiment in recipients of the control or PAT-perturbed inoculum, in fecal samples, or ileal (il), cecal (ce) and colonic (co) contents at sacrifice. f Unweighted Unifrac analysis of inoculum and fecal specimens from groups conventionalized with control (blue) or PAT-perturbed (red) microbiota, visualized in principal coordinate analysis (PCoA); the three components explain 43% of the total variance for each panel. Statistical analysis of intergroup UniFrac distances performed by Adonis test, with p-values shown Full size image

Dynamics of antibiotic-perturbed microbiota after transfer

We next assessed how microbial communities re-assemble in naive (GF) recipients to detail differences due to an inoculum that had been antibiotic-perturbed in the prior host (Fig. 4d, e and Supplementary Fig. 2c). Fecal samples were collected serially after gavage to evaluate the characteristics of the microbiota after transfer and to determine whether the PAT and control microbiota remain distinct. The inoculum from the PAT-exposed donors had lower phylogenetic richness and evenness than the inoculum from the controls (Fig. 4e, Supplementary Fig. 2c and Supplementary Table 1), reflecting the strong immediate PAT effects. With conventionalization, both inocula colonized the recipients with similar efficiency, with stable total bacterial levels reached by day 12 post-transfer (Supplementary Fig. 2a ); however, community α-diversity, composition and community structure markedly differed (Fig. 4d, e and Supplementary Fig. 2d). One-day after transfer, α-diversity in fecal samples fell for both groups of recipients, and then progressively increased until the time of sacrifice (Fig. 4e and Supplementary Table 1). However, community evenness was lower in the PAT-conventionalized mice until day 43, and remained lower in the ileum at sacrifice (Supplementary Fig. 2c and Supplementary Table 1). The community structures of the control and PAT inocula were distinct, and fecal samples in the recipients were different from day 12 after transfer through the experiment’s end (Fig. 4f ). PAT-conventionalized mice exhibited sex-specific differential taxonomic profiles beginning 6 days after transfer and remaining until the end of the experiment; nonetheless, community diversity in these mice remained distinct from control. These consistent differences indicate the fidelity of the microbial compositions colonizing naïve hosts over a period of months. The initial microbiota colonizing the control inocula recipients were enriched in Ruminococcus and Erysipelotrichaceae species, S24-7, and Dorea, whereas the recipients of the PAT inocula had higher representation of Akkermansia, Lactobacillus, Rikenellaceae, and particular Clostridia species (Supplementary Fig. 2d). One day after transfer, genus Lactobacillus bloomed in both groups of recipients, showing their role as pioneers even in the absence of milk (Fig. 4d). Over time, in both recipient groups, Enterobacteriaceae, Akkermansia and family S24-7 bloomed, with greater S24-7 abundances in the controls (Fig. 4d).

IgA-coated bacteria after control or PAT microbiota transfer

The transfer experiment allowed exploring the interplay between the PAT-induced sIgA alterations and the microbial populations bound to sIgA. Using IgA-Seq, we identified sIgA-coated bacteria after the conventionalization with control or PAT microbiota. The community composition of the sIgA-recognized bacteria differed between the PAT and control microbiota recipients (Supplementary Fig. 3a). Akkermansia was highly abundant in the sIgA+ fractions from both PAT and control samples for most of the experiment, indicating a strong sIgA affinity to Akkermansia if present in the microbial community, independent of exposure (Supplementary Fig. 3b ). Next, using the IgA-coating index (ICI), we identified taxa highly represented in the sIgA-negative fractions, including family Lachnospiraceae, S24-7, and Clostridium citroniae, whereas genus Lactobacillus and a particular Clostridia taxon were consistently over-represented in the sIgA-positive fractions (Supplementary Fig. 3c, d). In the sIgA-positive fractions in the control inocula-recipients, there was enhanced recognition of S24-7, Dorea, Lactobacillus, Erysipelotrichaceae, Clostridial taxa and Mollicutes, but in the PAT microbiota-recipients, only Bacteroides acidifaciens, genus Sutterella, and certain Clostridia were enriched, and only late in the experiment (Supplementary Fig. 3d). In total, these results indicate that compared to control, the PAT inocula recipients have a distinct IgA repertoire associated with a markedly altered community of host-interactive bacteria.

Associations between sIgA and microbial taxa

To assess the relationships between particular operational taxonomic unit﻿s (OTUs) and sIgA levels we applied multi-level, sparse PLS (sPLS) models. We focused on the fecal microbiota in the dosing experiment (Fig. 1 ), and used out-of-sample prediction of fecal microbes from the transfer experiment (Fig. 4) to drive model selection. Overall, 8 of the 53 OTUs commonly abundant in both experiments were selected by the sPLS model (Supplementary Fig. 4a, b ). Use of these taxa provided high goodness-of-fit (within sample r 2 = 0.77) and statistically significant model coefficients, at α = 0.05, and explained 31% of the sIgA variation from the transfer experiment. Three Lachnospiraceae family OTUs were negatively associated with ΔsIgA (Supplementary Fig. 4c), including Clostridium hathewayi, previously associated with reduced cytokine production in a murine colitis model22. Three OTUs showed PAT-dependent positive associations with ΔsIgA: Rikenellaceae, Enterobacteriaceae, and Lachnospiraceae/Blautia. The control-dependent ΔsIgA signature identified a Bifidobacterium sp. and a Clostridium genus (in the Erysipelotrichaceae family). In total, these analyses identify taxa individually associated with the altered sIgA phenotypes.

PAT induces conserved microbial network topology alterations

Next, to define the effect of PAT-induced perturbations on the structure of microbe–microbe associations, we developed a network model based on combining analysis of the dosing experiment (Fig. 1) and transfer experiment (Fig. 4). We inferred networks using SPIEC-EASI and used the stability metric to rank confidence in individual edges for the transfer experiment networks. Using this ranking to generate an ROC curve, we compared the ranked edges from the transfer experiment against the ‘true’ corresponding edges from the dosing experiment using the set of OTUs common to both experiments. Between the two experiments, we could compare Control-Control, PAT1-all PAT, and PAT3-all PAT relationships. Assessing the maximum harmonic mean of Precision and Recall (the F 1 score), more edges were recovered than expected by at-random prediction (dotted red line), under all conditions (Supplementary Fig. 4d). Overall, edges in the control mice were more likely to be recovered in the transfer experiment group (F 1 = 0.75). PAT vs PAT1 resulted in higher edge recovery (F 1 = 0.72) than PAT vs PAT3 (F 1 = 0.62), the expected outcome of additional antibiotic pulses. These results indicate that together with composition, condition-dependent association networks were largely conserved across the fecal transfer experiments, and lead to the hypothesis that the networks themselves help establish the new communities that drive the observed phenotypes.

To further examine this hypothesis, we considered the global topological properties of the inferred networks. Using graphlet correlation distance as a global measure comparing networks, and using isometric MDS as an analytic tool, we inferred low-dimensional embedding of the microbial association networks (Fig. 5). Notably, as expected, there was linear separation between the two experiments, but we also observed linear separation between the Control and PAT conditions across the experiments. These findings indicate that while transfer of distinct entire microbial communities affects both composition and inferred network topology in the recipients, comparing higher order network structures also can reveal global, transferable patterns of association.

Fig. 5 Network properties recapitulate experimental perturbation. Networks were inferred using the SPIEC-EASI pipeline. To compare graphs, a two-dimensional embedding of graphlet correlation distances (using multidimensional scaling (MDS)) was used with the network positions shown as circles (from the dosing experiment (see Fig. 1)) or triangles (from the transfer experiment (see Fig. 4)). The corresponding networks on the MDS plot were overlaid near their respective positions in the embedding, and shown in a force-directed layout. Within each graph, nodes are colored by OTU family, with size proportional to (geometric) mean abundance, and edge width proportional to confidence score, as determined by stability selection. Each panel shows the number of OTUs (identified as those present in > 20% of samples), corresponding number of edges and calculated edge density (predicted edges/number of pairs of nodes) Full size image

Roles of keystone taxa in the community dynamics

We next hypothesized that antibiotic interventions promote network fragility by altering representation of keystone taxa that mediate multiple intra-community fitness interactions. As such, we asked whether overall network robustness mediated via the keystone taxa was affected by the differential exposures. To test this hypothesis, we used natural connectivity as a stability metric of the inferred networks after simulated network “attacks”. Considering theoretical keystone indices, we assessed network stability after removing OTUs ordered by degree or betweenness centrality (Supplementary Fig. 5a), or at random (as a control) (Supplementary Fig. 5b). In the random attack setting, overall networks inferred from the dosing (Dose) experiment (Fig. 1) were more robust than transfer (Trans) experiment (Fig. 4) networks. In the degree and betweenness-based attack schemes, the Dose-PAT3 networks were relatively fragile, and were at the same level as both groups from the transfer experiment (Supplementary Fig. 5a–d ). Within the groups from the dosing experiment, Control and PAT were equally robust in the random and degree-based removal approaches, but PAT1 was more fragile to betweenness-attack schemes. These results indicate that if experimental interventions (i.e. antibiotic exposure or transfer to new hosts) targeted bottleneck keystone taxa, communities would be destabilized.

Finally, we also considered the particular OTUs identified as potential keystone taxa in the SPIEC-EASI-inferred association networks. For each network, we identified the three highest-scoring OTUs that were both hubs (having high node degree) and bottlenecks (as characterized by the highest betweenness centrality) in the network (Supplementary Fig. 5e ). Despite the differences between the experiments revealed by the network topology and robustness analysis, the taxonomic identities of these species were conserved within treatment type, across both the dosing and transfer experiments; in Controls, these were Lachnospiraceae, Blautia, Blautia producta, and Enterobacteriaceae, but only S24-7 in PAT1. Taken together, these results allow the prediction that perturbation of these taxa impacts the stability of the overall community, especially since their keystone network attributes were conserved even after transfer into GF mice.