Mother's milk microbiota is associated with the developing gut microbial consortia in very-low-birth-weight infants

Abstract

Mother's milk contains diverse bacterial communities, although their impact on microbial colonization in very-low-birth-weight (VLBW, <1,500 g) infants remains unknown. Here, we examine relationships between the microbiota in preterm mother's milk and the VLBW infant gut across initial hospitalization (n = 94 mother-infant dyads, 422 milk-stool pairs). Shared zero-radius operational taxonomic units (zOTUs) between milk-stool pairs account for ∼30%-40% of zOTUs in the VLBW infant's gut. We show dose-response relationships between intakes of several genera from milk and their concentrations in the infant's gut. These relationships and those related to microbial sharing change temporally and are modified by in-hospital feeding practices (especially direct breastfeeding) and maternal-infant antibiotic use. Correlations also exist between milk and stool microbial consortia, suggesting that multiple milk microbes may influence overall gut communities together. These results highlight that the mother's milk microbiota may shape the gut colonization of VLBW infants by delivering specific bacteria and through intricate microbial interactions.

Keywords: antibiotics; direct breastfeeding; donor milk; human milk; microbiome; microbiota; mother’s milk; nutrient fortification; preterm infant; very-low-birth-weight infant.

Graphical abstract

Figure 1

Microbial alpha diversity of paired milk-stool samples (A) Number of zOTUs and (B) Shannon index over time, stratified by sample type (n = 422 stools, n = 334 mother’s milk). Solid lines represent the mean alpha diversity over time, and shaded areas represent the 95% confidence interval. p values are from linear mixed-effects models adjusted for sample type, postnatal week, DNA extraction batch, and an interaction term between sample type and postnatal week. (C) Relationships between the number of zOTUs in paired milk-stool samples were examined for the entire cohort (in yellow) and by postnatal period and feeding variables of interest. Given the non-linear relationships observed (milk cut-point: zOTUs ≤ 200 or > 200), the number of zOTUs in paired milk-stool samples was investigated separately for each segment. p values (colored according to models stratified by each feeding variable of interest) are from unadjusted linear mixed-effects models due to sample size constraints within each segment. (D) Relationships with Shannon index in paired milk-stool samples were assessed using linear mixed-effects models adjusted for postnatal week, DNA extraction batch, infant sex assigned at birth, birth weight stratum, and feeding variables of interest. Interaction terms were removed from final models if they were not statistically significant (p > 0.05; denoted with NS). Abbreviations: zOTUs, zero-radius operational taxonomic units; MOM, mother’s milk; HMBF, human milk-based fortifier; BMBF, bovine milk-based fortifier.

Figure 2

Microbial beta diversity of paired milk and stool samples (A) Principal coordinate analysis plots of beta diversity metrics. R² and p values are from adonis models adjusted for postnatal week, DNA extraction batch, and participant identification. (B) UniFrac distances between paired milk-stool samples were assessed over time for the entire cohort (in yellow) and by feeding variables of interest. Solid lines represent the mean distance between paired samples, while shaded areas represent 95% confidence intervals. p values are from linear mixed-effects models adjusted for postnatal week, DNA extraction batch, infant sex assigned at birth, birth weight stratum, and feeding variables of interest. Interaction terms were tested between the feeding variables and postnatal week but were removed from final models if they were not statistically significant (p > 0.05; denoted with NS). Abbreviations: PC, principal component; MOM, mother’s milk; HMBF, human milk-based fortifier; BMBF, bovine milk-based fortifier.

Figure 3

Shared microbial taxa between paired milk-stool samples (A) Venn diagram showing the total number of zOTUs unique to milk and stool, and the total number of zOTUs that were shared in paired milk-stool samples. (B) Number of zOTUs shared between paired milk-stool samples over time. Solid lines represent the mean number of shared zOTUs over time, while the shaded areas represent 95% confidence intervals. p values are from linear mixed-effects models adjusted for postnatal week, DNA extraction batch, infant sex assigned at birth, birth weight stratum, and feeding variables of interest. (C) Shared zOTU ratio was calculated as the number of shared zOTUs between paired milk-stool samples divided by the total number of zOTUs in each corresponding stool sample. This model was adjusted as described in (B). (D) Relationships between postnatal week, feeding variables of interest, and the likelihood of sharing a greater number of zOTUs between milk-stool pairings were assessed using an adjusted repeated measures Poisson regression model as described in (B). (E) Alluvial diagram depicting the taxonomy of shared zOTUs between paired milk-stool samples. Node sizes represent the number of milk-stool pairings that shared a zOTU mapping back to the specified taxa. For visual clarity, only taxa that were shared in approximately 10% of milk-stool pairings are listed. Abbreviations: zOTUs, zero-radius operational taxonomic units; IRR, incidence rate ratio; MOM, mother’s milk; HMBF, human milk-based fortifier; BMBF, bovine milk-based fortifier.

Figure 4

Likelihood of zOTU sharing in paired milk-stool samples depending on postnatal period and feeding practices Repeated measures logistic regressions were used to assess how postnatal period and feeding practices influence the likelihood of sharing a zOTU within the 12 most commonly shared genera. Separate models were run for each genus, with the outcome being whether a paired milk-stool sample had a shared zOTU mapping to the specified genera (yes/no). Models were adjusted for postnatal week, DNA extraction batch, infant sex assigned at birth, birth weight stratum, and feeding variables of interest. p values were adjusted using a Benjamini-Hochberg false discovery rate to account for multiple comparisons. Results for unclassified Enterobacteriaceae are not reported since almost all milk-stool pairs had a shared zOTU mapping back to this family, leading to model convergence issues. Abbreviations: OR, odds ratio; CI, confidence interval; MOM, mother’s milk; HMBF, human milk-based fortifier; BMBF, bovine milk-based fortifier.

Figure 5

Dose-response relationships between bacterial intakes from mother’s milk and their concentrations in infant stools (A) Relationships between cumulative 3-day milk bacterial intakes of commonly shared taxa and their concentrations in corresponding infant stools were assessed using unadjusted linear mixed-effects models. Sample pairs were included in these models if both milk and stool contained the respective taxa, and p values were adjusted using a Benjamini-Hochberg false discovery rate. Solid lines represent the mean and shaded areas represent the 95% confidence interval. (B) Relationship between the cumulative 3-day total bacterial intake from mother’s milk and the total bacterial concentration in stool, as described in (A). (C) Summary of all p values from unadjusted linear mixed-effects models in (A) and (B) in addition to models stratified by postnatal period and feeding variables of interest. p values from the main effects were FDR-adjusted to account for multiple comparisons. (D) Statistically significant relationships between cumulative 3-day intakes of total bacteria and commonly shared taxa with corresponding concentrations in infant stools, stratified by postnatal period and feeding variables of interest. Abbreviations: NS, non-significant; FDR, false discovery rate; MOM, mother’s milk; HMBF, human milk-based fortifier; BMBF, bovine milk-based fortifier.

Figure 6

Correlations between milk bacterial intakes and concentrations in infant stools (A) Chord diagram displaying statistically significant relationships between cumulative 3-day milk bacterial intakes of commonly shared genera and their concentrations infant stools. The width of the linkage is proportional to the strength of the correlation. (B) Heatmap displaying spearman rank correlations between 3-day milk bacterial intakes and concentrations in infant stools for the entire cohort. ∗p < 0.05; ∗∗p < 0.01; ∗∗∗p < 0.001. (C) Heatmaps showing the spearman rank correlations as described in (A) stratified by postnatal period. Black borders are used to highlight key findings across heatmaps.

Figure 7

Relationships between the microbiotas in mother’s milk and infant stools according to maternal-infant antibiotic exposure Relationships with (A) the number of zOTUs and (B) Shannon index in paired milk-stool samples stratified by antibiotic group: maternal low-infant low (n = 22 mother-infant dyads, n = 80 paired milk-stool samples), maternal low-infant high (n = 34 mother-infant dyads, n = 166 paired milk-stool samples), maternal high-infant low (n = 11 mother-infant dyads, n = 52 paired milk-stool samples), and maternal high-infant high (n = 27 mother-infant dyads, n = 124 paired milk-stool samples). High antibiotic use was defined as >1 day for mothers and >3 days for infants. Solid lines represent the mean alpha diversity, and shaded areas represent the 95% confidence interval. p values are from unadjusted linear mixed-effects models. (C) Heatmaps showing spearman rank correlations between cumulative 3-day milk bacterial intakes of commonly shared genera and their concentrations in infant stools stratified by antibiotic group. ∗p < 0.05; ∗∗p < 0.01; ∗∗∗p < 0.001.