High-Fiber Whole-Food Diet Changes Gut Bacteria but Not Short-Chain Fatty Acids in Stool

Twenty-six individuals participated in a CURE course at UC Irvine, designed to tandemly investigate pedagogical methods (25) and the role of fiber in the microbiome. We collected nutritional data over 3 weeks from all 26 individuals who initially began the intervention (1 week before and 2 weeks during the dietary intervention) (Fig. 1A). We collated the total amount of macronutrients consumed per day, including fiber, protein, carbohydrates, fats, as well as overall calories (Fig. 1B to F). Additionally, we informally surveyed food items the study participants frequently used to supplement their meal plans, beyond the meals supplied from Thistle, and found that items such as fiber-fortified cereals, lentils or beans, and berries were common (25). For the intervention, subjects increased their average fiber consumption from 21.0 g/day (± 14.2 g/day) before the intervention to 46.4 g/day (± 12.5 g/day) during the intervention (Wilcoxon rank sum test, P < 0.0001) (Fig. 1B). While these dietary shifts increased carbohydrate intake by an average of 84% (36 g) during the intervention (P = 0.013), other macronutrients measured, such as calories, fat, and proteins, did not significantly change (P > 0.05) (Fig. 1C to F).

To evaluate whether increased fiber consumption contributed to shifts in the gut microbiome, we characterized the microbial communities from 20 individuals using 86 shotgun metagenomic libraries collected before and after the fiber intervention (Fig. 2A). The alpha diversity of microbial taxa decreased during the high-fiber diet intervention, as measured by the Shannon diversity index (Fig. 2B) (Wilcoxon rank sum test, P < 0.05). Using alternative approaches to assess taxonomy and diversity (see Materials and Methods) showed either no change or supported the decreasing trend of diversity during the intervention period (Fig. S1).

Despite little difference in alpha diversity, beta diversity changed significantly in response to a high-fiber diet. Multivariate analysis of marker gene abundances showed that most of the variation in microbiome composition could be explained by the individual (permutational analysis of variance [PERMANOVA] for main individual effect, R² = 0.78, P < 0.001) (see Table S2 in the supplemental material). The diet intervention shifted the microbial composition of the entire study cohort significantly (main intervention effect, R² = 0.014, P < 0.001). Within samples from each individual, the pre- and postdiet intervention samples explain significant variation in the community composition (intervention-by-individual effect, R² = 0.083, P < 0.001). A linear mixed-effects (LME) model confirmed these results, which identified diet as a significant determinant of an individual’s microbiome composition (LME, P < 0.01). Individual gut microbiome samples grouped together in nonmetric multidimensional space (NMDS; Fig. 2C), further providing support that each individual is associated with a unique microbiome. Some individual gut microbiomes (i.e., individual 13) were more distinct from others (Fig. 2C, inset). Additionally, we used Eubacterium rectale (due to its high coverage in our data) to ask whether the diet intervention had an impact at the strain level. Strains were highly individual specific (PERMANOVA main individual effect, R² = 0.99, P < 0.001) and did not change in response to increased fiber intake (P > 0.05; Fig. 2D).

We next parsed the taxonomic data to assess which microbial taxa increased or decreased in response to the diet intervention. One species in the family Lachnospiraceae was significantly negatively associated with increasing fiber intake (Spearman r = −0.43, q = 0.01) (Fig. S2A and B). Coprococcus sp. and Anaerostipes hadrus were both positively associated with increasing fiber intake, but this association was not significant when P values were false discovery rate (FDR) corrected for multiple comparisons (r = 0.32, q = 0.33 for both species) (Fig. S2A). Furthermore, positive linear coefficients of a PERMANOVA model, which detect differences between community compositions due to the diet intervention, included genera such as Bifidobacterium, Bacteroides, and Prevotella (Fig. 3A). Conversely, Blautia and Ruminococcus contributed negative linear coefficients to the PERMANOVA model (Fig. 3A).

Of the 105 microbial genera detected in this study, Bifidobacterium was the strongest predictor genus for the postintervention microbiomes (Fig. 3A). Indeed, taxonomic analysis of the metagenomic samples identified Bifidobacterium abundances increasing, on average, 1.4-fold between the pre- and postintervention periods (Fig. S3A). Further, we identified several species of Bifidobacterium present within and across individuals, with B. adolescentis being the most abundant species on average (Fig. 3B). When we investigated the taxonomic profiles at the species level, we found that B. adolescentis, B. biavatii, B. breve, B. longum, and B. ruminantium all increased in mean abundance on a high-fiber diet, whereas the other, less abundant species exhibited no change or decreased in abundance (Fig. S3B).

Given that Bifidobacterium was the strongest predictor genus in the postfiber gut microbiomes, we employed a targeted analysis of the diversity within Bifidobacterium to examine species-level patterns. Specifically, we applied targeted amplicon approaches to amplify the groEL gene, a conserved phylogenetic marker gene to track Bifidobacterium diversity (Fig. S4). Using phylogenetic inference of the groEL gene, we compared the observed Bifidobacterium diversity observed at the community level to our targeted analysis of the groEL gene. Similar to the metagenomic analysis, we found that individuals were largely comprised of B. adolescentis and B. longum, with six other abundant species of Bifidobacterium (Fig. 3C). This analysis also revealed extensive Bifidobacterium diversity within the human gut, detecting 22 species across all individuals.

Since Bifidobacterium species are known to participate in cross-feeding with other gut microbes (reviewed in reference 29), we next assessed the cooccurrence of Bifidobacterium with other genera. Bifidobacterium was positively correlated (r = 0.43, q = 0.001) with an increasing abundance of Lactobacillus and negatively correlated with Roseburia (r = −0.49, q = 0.0002) and Ruminococcus (r = −0.38, q = 0.007) (Fig. S5), suggesting species interactions between these taxa.

Our results demonstrate that a shift in dietary fiber consumption influenced compositional changes in the gut microbial community. As such, we sought to correlate the observed taxonomic shifts to functional shifts, particularly the enrichment of genes related to carbohydrate degradation. Despite taxonomic shifts at the individual level, we observed no changes in the overall abundance (average number of normalized reads) mapping to gene families for glycoside hydrolases (GH) (Wilcoxon, P = 0.42), carbohydrate esterases (P = 0.58), glycoside transferases (P = 0.73), and polysaccharide lyases (P = 0.77) as a result of the intervention (Fig. S6A). No individual families of GH and polysaccharide lyase CAZy classes changed in abundance during the intervention when corrected for multiple comparisons (Wilcoxon, P > 0.05) (Fig. S6B). Further, the diversity (Fig. 4A) (n = 106 families; ANOVA, P > 0.05) and composition (PERMANOVA, P > 0.05) of GHs detected were indistinguishable between the pre- and postintervention samples. Compositional analysis of all genes identified by HUMAnN revealed no individual signature (PERMANOVA main individual effect, R² = 0.017, P > 0.05) (Fig. 4B) and no shifts in response to a high-fiber diet (intervention by individual effect, R² = 0.015, P > 0.05). We performed a linear discriminant analysis to determine if there were pathways that were differentially abundant due to the diet intervention and found inositol degradation (in addition to several unintegrated pathways) to be increased in abundance on a high-fiber diet (Fig. 4C, Fig. S7). For the pathways involved in SCFA metabolism, we found no significant (Wilcoxon, P > 0.05) changes as a result of the high-fiber diet (Fig. 4D).

While the presence of genes related to SCFA production provides insights into the functional changes of the microbiome, these genes only reflect the genomic potential to process these pathways. Therefore, we applied a targeted gas chromatography-mass spectrometry (GC-MS) analysis to 149 samples from 18 individuals for the presence of SCFA molecules. Across the intervention period, the average abundance of acetate, propionate, butyrate, and valerate increased (Table 1); however, these increases were not statistically significant (LME, P > 0.05) (Fig. 5). For the eight individuals with samples run in duplicate, three had biological differences that were greater than the technical variation seen in the duplicates (Fig. S8). Acetate had the least technical variation (mean coefficient of variation [CV], 45%), followed by propionate (mean CV, 51%) (Table 1, Fig. S8).

Past studies have shown that an increase in the diversity of dietary foods could lead to an increase in microbial diversity (24). Moreover, individuals living in rural societies often harbor far greater gut microbial diversity than individuals from western societies (43–45), which may in part be linked to a greater proportion of plant-based polysaccharide intake. However, we did not measure an increase in species diversity (alpha diversity) after subjects consumed >40 g of fiber from a diverse set of foods (Fig. 2B). These results could be attributed to the brevity of the intervention, as the rapid change in dietary composition may result in the loss of microbes poorly adapted to recalcitrant carbohydrates. Similarly, other studies have reported finding no increases in alpha diversity as a result of fiber intake (32, 46–49), which indicates a trade-off where fiber degraders increased while other taxa decreased. Although alpha diversity was unaffected, we did observe a significant impact of the high-fiber diet on microbial community composition (beta diversity) (Fig. 2). The composition of microbial communities within individuals shifted ∼8% during the intervention period. We found changes in communities to be at broader taxonomic levels than the strain level. We could examine strains of E. rectale due to its high coverage in our data, and we showed these strains stayed constant and individual specific during the intervention (Fig. 2D). Future work should determine if this pattern holds up for other species. While we suspect the high-fiber diet treatment played an instrumental role in shifting the microbial composition, we cannot rule out other factors, such as host genetics or nondietary behaviors. As discussed, many food-based fiber interventions have shown mixed results in changing the microbial communities (6, 50). The drastic increase in fiber from a variety of foods may lead to rapid shifts in community composition over the 2-week period. Changes in community composition pre- and postintervention were largely driven by shifts in known fiber degraders, such as Bifidobacterium, Bacteroides, and Prevotella (Fig. 3A).

We expected the taxonomic shifts in the microbiota would be associated with changes in the functional potential of the microbial communities (Fig. 4). While we initially hypothesized that a high-fiber diet would increase the abundance or diversity of carbohydrate active enzymes, we did not detect changes associated with the intervention (see Fig. S6 in the supplemental material). Our findings support a similar result showing no difference in CAZy abundance due to increased fiber intake (49). We acknowledge that sequencing depth is an important consideration in the detection of genes; increasing reads beyond our ∼1.3 million paired-end reads (average per sample; Table S1) may allow for greater detection. However, we did find a notable increase in the abundance of genes mapping to the inositol pathway (Fig. 4C). We suspect that the increased consumption of fiber-fortified cereals and legumes, which contain higher levels of inositol, during the diet intervention allowed for an expansion in organisms capable of breaking down this sugar. There is substantial interest in the role of inositol (specifically phytic acid) in its protective role against colon cancer and other metabolic disorders (51, 52). Next, we assessed whether genes involved in SCFA metabolism changed in abundance during the intervention. Although appreciable cross-feeding between lactate-producing Bifidobacterium spp. and butyrogenic bacteria has been shown (53), we did not find significant increases in genes involved in various SCFA metabolic pathways (Fig. 4D). This further supports our results showing no clear correlations between Bifidobacterium spp. and butyrate producers within our diet intervention (Fig. S5). Indeed, we would not be the first to suggest that these complex trophic interactions require more time to establish (17). Rather, our results suggest that while broad taxonomic shifts occur, these do not correspond to changes in functional potential, and fine-scale (intraspecies) shifts are less susceptible to dietary shifts on short-term timescales.

Many studies have indicated that bifidobacteria (often identified as the genus Bifidobacterium by fluorescent in situ hybridization probes, PCR, or DNA sequencing) are highly abundant in the gut following increased fiber intake (meta-analysis of 51 studies [30]). The increased abundance of Bifidobacterium is somewhat unsurprising, as they harbor numerous genetic components, such as carbohydrate active enzymes, that make them especially adapted to a fiber-rich diet (54). In one study, both resistant potato starch and inulin increased the relative abundance of Bifidobacterium spp.; however, the 16S amplicon sequencing in this study did not have the resolving power to identify which species of Bifidobacterium were increasing (17). Using a targeted amplicon approach, the groEL gene has been shown to delineate species of Bifidobacterium that otherwise share >99% sequence identity in the 16S rRNA gene, making it a robust marker gene for analyzing within-genus species diversity (55). In our study, the most abundant species of Bifidobacterium were B. adolescentis and B. longum, both of which are efficient degraders of plant-based fructooligosaccharides (FOS) and produce acetate and lactate in the process (56). Mirroring our results, other studies have found selective increases in certain species of Bifidobacterium as a result of carbohydrate intake; for example, in one study, the intake of inulin resulted in a greater increase of B. adolescentis (57). We speculate that on a high-fiber diet, bifidobacteria are the initial members of the community accessing fiber substrates, easily adapted to utilize various FOS, and pivotal to the creation of the initial metabolic cross-feeding networks. Future studies should extend the intervention period to examine the dynamics of longer-term trophic interactions in response to increased dietary fiber intake.

While SCFAs did generally increase during the diet intervention, trending toward their naturally occurring gut ratio of 3:1:1 (acetate-propionate-butyrate) (12, 58, 59), we did not observe a statistically significant increase in SCFAs postintervention. Static fecal concentrations of SCFAs may not reflect the total pool of molecules fluxing through a given individual, as the molecules are preferred substrates of the cells lining the gut epithelia (15). It is also possible that the intervention period was too short to observe increases in SCFA abundances.

It should be noted that accurate SCFA measurements are notoriously difficult. Our examination of technical variability within 44 samples from eight individuals showed that technical variation between preintervention replicates or postintervention replicates was greater than the average difference between pre- and postintervention for any given SCFA. One study reported high intrafecal variability of butyrate quantification (coefficient of variation [CV], 38%) before optimizing a freeze-drying method (60). Numerous studies have indicated the benefit of SCFAs to human health (5, 61), yet the heterogeneity in reported acetate, propionate, and butyrate abundances remains high. In one meta-analysis of fiber studies, only butyrate was generally found to increase with fiber intake, yet the heterogeneity of reported results was 70% (I²), similar to other SCFAs analyzed (30). Outside of technical limitations, shifts in microbial community structure are not predictive of changes in static measurements of fecal SCFA abundances (62). The difficulty of finding meaningful correlations between microbiome composition and SCFA abundances likely reflects a failure to measure both circulating and fecal SCFAs across time in conjunction with microbial abundances. Indeed, it has been observed that fecal levels of acetate are inversely related with the rate of its absorption (63). Future studies are needed to confirm whether correlation analysis between fecal SCFAs and microbiome composition is a useful tool to understand the interplay between microbiome, SCFAs, and health.

In sum, our results indicate that gut microbial communities are malleable to an influx in recalcitrant carbohydrates, contributing to significant community and functional shifts in certain metabolic pathways. However, these compositional changes did not correspond to broad functional changes, at least over the short-term timescales for this intervention. Further studies exploring the impact of timing and composition of dietary fiber interventions, particularly while taking into account the starting composition of the gut microbiomes of study participants, are critical for understanding the generalizability of fiber interventions for engineering microbiomes. Increasing fiber intake could have the most impact in contexts where low gut microbial diversity increases the risk of C. difficile infection, such as for nursing home residents and cancer patients or after antibiotic treatment.

Twenty-six UC Irvine students and instructors volunteered for a 3-week high-fiber diet intervention study (Fig. 1A); only 22 individuals elected to provide stool samples for microbiome analyses (from 20 of whom we recovered enough sequence data for analysis; see Table S1 in the supplemental material). The dietary intervention was approved by UC Irvine IRB number 2018-4297. For the first week of the study, all participants consumed their normal diets, tracking all nutritional information using the smartphone application MyFitnessPal (MyFitnessPal, Inc.). Prior to the end of week one, each subject provided three fecal samples from 3 days within the first week. The intervention commenced in week two, when participants were instructed to raise their dietary fiber intake to approximately 40 g per day. To assist with the dietary shifts, we provided 10 meals per week, with ∼15 g of fiber and ∼5.8 unique fruits or vegetables per meal, from the food delivery service Thistle (San Francisco, CA, USA). During week three, subjects were encouraged to further increase fiber intake to ∼50 g of fiber per day. Subjects provided three fecal samples from 3 days during week three, concluding the intervention period. As part of the CURE course, students were educated on human health, dietary information on high-fiber meals, the human gut microbiome, and the quantitative methods for microbiome analyses (from DNA extraction and library preparation to metagenome and statistical analyses), as previously described (25).

Subjects were given materials to collect fecal samples at home. Each stool sample was split into three 2-ml tubes by the individual and immediately stored in the freezer. When convenient, students transported their anonymized and coded samples using cold packs and insulated boxes to a common laboratory freezer. Upon the conclusion of the intervention period (week 1 or 3), all samples were transported to a −20°C freezer.

To characterize the bacterial community composition of the samples, DNA was extracted with the ZymoBIOMICS 96 DNA kit (product D4309) from Zymo Research using the manufacturer’s suggested protocol. Sequencing libraries were prepared using the Illumina Nextera kit and methods described in Baym et al. (64). Briefly, DNA was diluted to 0.5 ng/μl and added to 0.25 μl of Nextera enzyme and 1.25 μl of tagmentation buffer. This mixture was incubated at 55°C for 10 min and then placed on ice for the remainder of the protocol. Barcodes were added using Phusion polymerase (New England Biolabs), and excess adaptors were cleaned using AMPure XP (Beckman Coulter Life Sciences) magnetic beads. Quality and concentration were assessed using a PicoGreen assay (ThermoFisher), and the distribution of fragment sizes was determined using a Bioanalyzer. These libraries were loaded onto the Illumina Next-Seq 500 at 1.8 pM concentrations and sequenced using Illumina’s mid-output kit for 75-bp paired-end sequencing, resulting in a total of 144,023,583 reads and an average of 1,425,976 reads/sample (maximum, 5,902,966; minimum, 7) (Table S1).

To characterize the genetic diversity of Bifidobacterium at a finer genetic scale than could be assayed by metagenomics, we used genus-specific primers to target this group for sequencing (65). Sequencing libraries were prepared by setting up initial 25-μl PCRs with AccuStart II PCR ToughMix (2×), the groEL forward primer (5′-TCGTCGGCAGCGTCAGATGTGTATAAGAGACAGTCCGATTACGAYCGYGAGAAGCT-3′, 20 μM), and the groEL reverse primer (5′-GTCTCGTGGGCTCGGAGATGTGTATAAGAGACAGCSGCYTCGGTSGTCAGGAACAG-3′, 20 μM). The initial PCR ran for 28 cycles 95°C for 30 s, 60°C for 30 s, and 72°C for 50 s, followed by the addition of 0.5 μl of dual Nextera XT index (Illumina) to each sample, proceeding with an additional 8 cycles of 95°C for 30 s, 60°C for 30 s, and 72°C for 50 s. Amplicons were pooled based on visual quantification of the bands on an agarose gel and purified using magnetic Speed Beads. The pool was run on a MiSeq PE 300 at the University of California Irvine's Genetic High-Throughput Facility, resulting in a total of 20,052,935 reads and an average of 185,675 reads/sample (maximum, 6,815,601; minimum, 155).

SCFA extractions were done by following the methods of Zhao et al. (66). One hundred milligrams of fecal material was added to 1 ml of high-performance liquid chromatography (HPLC)-grade water and vortexed for 2 min. Ten microliters of 6N HCl was added to the fecal slurry and vortexed briefly. This mixture was incubated at room temperature for 10 min with occasional shaking. Afterward, the mixture was centrifuged at 14,000 × g for 1 min, and 400 μl of the supernatant was transferred to a new tube, which was then filtered through a 0.22-μm filter. An aliquot (200 μl) of this suspension was then transferred to a glass vial with a 0.2-ml vial insert and stored at −20°C. When running the sample, 10 μl of an internal standard of 10 mM ethyl butyrate was added to the extraction prior to the run. Before running each sample, the instrument was calibrated using a standard comprising 100 mg/liter acetate, propionate, isobutyrate, butyrate, isovalerate, valerate, and ethyl butyrate. Six samples were run on an Agilent 7890A gas chromatograph with dual-column flame ionization detectors. Two microliters per extracted sample was hand injected on a stainless steel column (2 m by 3.2 mm) containing 10% SP-1000 and 1% H₃PO₄ on 100/120 Chromosorb W AW (Supelco, Inc., Bellefonte, PA, USA). The flow rate of the N₂ carrier gas was 26.14 ml/min. Between sets of six samples, the instrument was washed using water and phosphoric acid. Peaks were autointegrated using ChemStation v1.0 on a PC running Windows 2000 (Microsoft). A subset of samples (n = 44 from 8 individuals) were run in duplicate to examine technical variation (see the coefficient of variation [CV] in Table 1), and the average CV was 55%.

Raw shotgun metagenome sequences were filtered using Prinseq v0.20.4 (67) to remove sequences that had a mean quality score of 30 or less. Reads from human DNA were also removed by aligning the filtered reads to the human genome (hg38), using Bowtie2 v2.2.7 (68), and keeping the reads that failed to align. A total of 130,755,383 paired-end reads (average of 1,294,607 nonhuman reads/sample) were retained and passed through MIDAS, which assigns taxonomy to short read data using a marker gene approach (69). Species counts per sample represent the average of 100 subsamples, rarefied to 900 sequences per sample using the EcolUtils (v0.1) package in R. Taxonomy was also assessed using IGGsearch (70). To analyze functional differences related to SCFA metabolism between high- and low-fiber treatment groups, HUMAnN3 (71) was used with default parameters. All pathways within the MetaCyc pathway class “Fermentation to Short-Chain Fatty Acids” were searched for within the HUMANnN pathway output, which resulted in nine pathways used for analysis (72). For genes related to carbohydrate breakdown, we translated reads using Prodigal (73) to predict open reading frames (ORFs) and searched all ORFs against the Pfam database (74) with hmmer/3.1b2 (75). Resulting PFAM annotations were then screened against CAZyDB.07202017 (76) with BLAST/2.8.1 (77) using alignments of >70% amino acid identity and 30% coverage. Alpha diversity and PERMANOVA analyses were performed using the Vegan v2.5-6 (78) package in R (79). Nonmetric multidimensional analysis was done using the metaMDS function in Vegan on Bray-Curtis distances. StrainPhlAn (80), under default parameters, was used to analyze strain-level variation within the metagenomes. To root the phylogenetic tree, Prosthecochloris aestuarii (accession no. GCA_000020625) was used, and two reference genomes of Eubacterium rectale (accession no. GCA_000209935 and GCA_001404855).

We downloaded 780 genomes from the genus Bifidobacterium on the PATRIC database (81). All genomes were screened for completeness by searching for 21 single-copy ribosomal marker genes using Prodigal (73) and HMMer v3.1b2 (75) with an E value of 1 × 10⁻¹⁰. The remaining 578 genomes were used to create a multilocus, concatenated phylogeny of the ribosomal marker genes with ClustalO v1.2.0 (82) to produce a 4,272-amino-acid alignment for phylogenetic analysis using RAxML v8.0.0 (83) under the PROTGAMMABLOSUM62 model for 100 replicates. Next, we parsed the filtered genomes for the groEL gene sequences by using 260 nonredundant gene sequences to build a groEL phylogeny under parameters identical to those of the whole-genome analysis. The groEL amino acid sequences, alignment, and phylogeny were used to construct BLASTp, HMMer, and pplacer reference databases for metagenomic analyses.

For each groEL amplicon library, sequences were quality trimmed and adapters were removed with BBDuk (84) (qtrim = rl trimq = 10 ktrim = r k = 25). Paired-end sequences were merged together with BBMerge (84), and, if paired reads did not overlap, only the forward read was retained. The reads were then searched against the groEL reference databases using BLAT (85) and hmmsearch, respectively. Passed reads were aligned with ClustalO to the pplacer reference package and placed on the groEL reference phylogeny using pplacer v1.1.alpha17 (86). Relative abundance was calculated from the single branch assignments and aggregated at the species level to be normalized by the total number of extracted groEL gene sequences. We show that the phylogenetic relationship between species of Bifidobacterium based on the groEL gene closely reflects a phylogeny based on 21 single-copy marker genes from 578 Bifidobacterium genomes (Fig. S4).

PERMANOVA was conducted on Bray-Curtis dissimilarities at the genus level with 999 permutations using the Adonis test in the vegan package in R (described below). We tested the effect of the intervention (pre- versus postfiber increase), the effect of the individual, and the interaction between these two factors. Genus contributions to significant results from the PERMANOVA model were determined by passing the resulting PERMANOVA object through the coefficients function found in the base Stats package in R. A similar procedure was used to analyze compositional differences between CAZy enzymes and HUMAnN gene predictions in the metagenomes, with permutations on Euclidean distances. Linear mixed-effects modeling, using the nLME package (79) in R, were also conducted for comparison, because they take repeated measures into account. Specifically, to support the PERMANOVA analysis of beta diversity, an LME was performed on the rank-transformed first principal coordinate of a principal coordinate analysis on the Bray-Curtis community dissimilarity matrix. Individual was used as the random effect and the model used the default autoregressive (Lag 1) structure (AR1) for regression across a time series. For the functional analyses, reads analyzed using HUMAnN3 were normalized by copies per million; CAZy reads were normalized to the total number of reads per metagenome and compared using Wilcoxon rank sum test. Gene features for HUMAnN were reduced by analyzing only unstratified data, for which 70% of samples had nonzero reads mapping to each feature. HUMAnN pathway abundances were analyzed in their entirety with stratification and without feature reduction. Lefse (87) was used to determine pathways which may differentiate pre- versus postintervention samples. Wilcoxon rank sum tests were also used to compare nutritional and gene differences between intervention periods when residuals were not normally distributed and reads or macronutrients were averaged out within individuals (by treatment) to account for repeated measures. When normality assumptions of residuals were met (tested using the Shapiro-Wilk test), ANOVAs were used. To assess which taxa were correlated with changing amounts of fiber, all species within each sample (the rarefied species abundance matrix) and fiber were correlated using the Corrr package v0.4.2 (88) in R. To analyze which genera cocorrelate with the genus Bifidobacterium, Spearman correlations were used, and, where appropriate, P values were corrected (q value) for multiple comparisons using a false discovery rate cutoff of 0.05. To assess significance of strains between individuals, cophenetic distances were calculated on the RAxML tree output from StrainPhlAn and passed into the above-described PERMANOVA model.

All scripts are stored on GitHub (https://github.com/aoliver44/Fiber-Analysis↗). All metagenomic and amplicon sequences are available from NCBI under the BioProject no. PRJNA647720↗. Metadata linking the shotgun metagenomes and groEL sequences with the appropriate sample ID and intervention can be found in Table S1.

All scripts are stored on GitHub (https://github.com/aoliver44/Fiber-Analysis↗). All metagenomic and amplicon sequences are available from NCBI under the BioProject no. PRJNA647720↗. Metadata linking the shotgun metagenomes and groEL sequences with the appropriate sample ID and intervention can be found in Table S1.