The neonatal gut microbiota: A role in the encephalopathy of prematurity

Abstract

Preterm birth correlates with brain dysmaturation and neurocognitive impairment. The gut microbiome associates with behavioral outcomes in typical development, but its relationship with neurodevelopment in preterm infants is unknown. We characterize fecal microbiome in a cohort of 147 neonates enriched for very preterm birth using 16S-based and shotgun metagenomic sequencing. Delivery mode strongly correlates with the preterm microbiome shortly after birth. Low birth gestational age, infant sex assigned at birth, and antibiotics associate with microbiome composition at neonatal intensive care unit discharge. We integrate these data with term-equivalent structural and diffusion brain MRI. Bacterial community composition associates with MRI features of encephalopathy of prematurity. Particularly, abundances of Escherichia coli and Klebsiella spp. correlate with microstructural parameters in deep and cortical gray matter. Metagenome functional capacity analyses suggest that these bacteria may interact with brain microstructure via tryptophan and propionate metabolism. This study indicates that the gut microbiome associates with brain development following preterm birth.

Keywords: brain MRI; encephalopathy of prematurity; gut microbiome; gut-brain modules; microbiome-gut-brain axis; neonate; preterm.

Graphical abstract

Figure 1

Overview of microbiota profiles in neonates based on 16S rRNA gene sequencing (A) Relative abundances of the 20 most abundant amplicon sequence variants (ASVs) identified across the dataset are visualized per sample, with all other ASVs grouped together as residuals. Samples are ordered based on hierarchical clustering of the Bray-Curtis dissimilarity matrix using average linkage (see dendrogram). (B) Non-metric multidimensional scaling (NMDS) plot based on Bray-Curtis dissimilarity between samples; data points, and ellipses are colored by sample type. The ellipses denote the standard deviation of data points belonging to each sample type, with the center points of the ellipses calculated using the mean of the coordinates per group. (C) Microbiota alpha diversity measured by Shannon index (left) and observed ASVs (right) presented as boxplots and individual data points. ∗∗∗ indicates q < 0.001 in pairwise comparisons using emmeans following linear mixed effects model comparing alpha diversity indices between the groups and time points. (D) Differentially abundant ASVs in association with preterm status at time point 1. Bar plots depict MaAsLin2 analysis results. ASVs present with at least 1% of abundance in at least 5% of samples were analyzed (10 ASVs) and significant results are shown (Benjamini-Hochberg [BH] corrected p < 0.25 as default). Lengths of the bars correspond to the MaAsLin2 model coefficient, which relates to the strength of the association. Error bars indicate the standard error (SE) of the model coefficient. MaAsLin2 models were adjusted for postnatal age at sampling. Sample sizes: term time point 1 = 12, preterm time point 1 = 58, preterm time point 2 = 103. See Figure S2 for overview of microbiome profiles in preterm neonates arising from shotgun sequencing.

Figure 2

Covariates associated with preterm infant gut microbiota (A) Univariable PERMANOVA results showing the association between perinatal variables and the gut bacterial community composition at each time point and for each data type. Left: ASV from 16S rRNA sequencing, middle: species from shotgun sequencing, and right: gut metabolic modules (GMMs) calculated from KEGG orthologs arising from shotgun sequencing. The variance explained is estimated for each variable independently and is indicated by a percentage/blue shades. Significance of PERMANOVA was based on 1,000 permutations and was adjusted for multiple comparisons using the Benjamini-Hochberg (BH) method; asterisks denote statistical significance (ˆq ≤ 0.1, ∗q ≤ 0.05, ∗∗q ≤ 0.01). (B and C) Differentially abundant ASVs in association with perinatal factors at time point 1 (B) and 2 (C). Bar plots depict MaAsLin2 analysis results. ASVs present with at least 1% of abundance in at least 5% of samples were analyzed (14 ASVs for time point 1 and 21 for time point 2) and significant results are shown (BH-corrected p < 0.25 as default). Bars are colored according to the covariate they are associated with. Lengths of the bars correspond to the MaAsLin2-model coefficient, which relates to the strength of the association. Error bars indicate the standard error (SE) of the model coefficient. In baseline models, we adjusted for postnatal age (time point 1), or GA at birth and sample collection (time point 2); in full adjusted models, all covariates with q value <0.1 from univariable PERMANOVAs were tested simultaneously. Here, GA at birth was dichotomized to group the infants into extremely (GA at birth <28 completed weeks) and very (GA at birth <32 completed weeks) preterm. Sample sizes for 16S rRNA sequencing: preterm time point 1 = 58, preterm time point 2 = 103; sample sizes for shotgun sequencing: preterm time point 1 = 23, preterm time point 2 = 97. See Tables S2, S3, S4, and S5 for detailed MaAsLin2 results and Table S6 for alpha diversity associations.

Figure 3

Dimensionality reduction of the microbiota community composition data Bacterial ASV correlations with the first four orthogonal principal coordinates (PCo-s), showing the top 20 strongest correlations for each PCo. The percentage refers to the variance explained by each of the PCo-s. Red indicates positive and blue negative correlations between the PCo-s and relative abundances of ASVs. Sample size n = 79 (linked MRI and microbiome data). See Figure S3 for bacterial species (shotgun sequencing) correlations with the PCo-s.

Figure 4

Microbiota associations with MRI features of encephalopathy of prematurity (A) Regression results for brain volumetric measures. (B) Regression results for brain microstructural measures. Models are adjusted for gestational age at birth and at scan; microbiota PCo-s and alpha diversity metrics were adjusted for gestational age at sampling via linear regression, retaining the residuals. Points correspond to the standardized model coefficient. Error bars indicate the SE of the model coefficient. Full color points indicate nominal p value <0.05; asterisks (∗) indicate Benjamini-Hochberg (BH) method-adjusted p value < 0.25. Red indicates positive and blue negative associations. Relative volumes are calculated by normalizing to total tissue volume (the sum of the volumes of cortical gray matter, white matter, deep gray matter, cerebellum, brainstem, hippocampi, and amygdalae). FA, fractional anisotropy; RD, radial diffusivity; NDI, neurite density index; ODI, orientation dispersion index; ISO, isotropic volume fraction; cGM, cortical gray matter; dGM, deep gray matter, CB, cerebellum; sulc, sulcal depth; GI, gyrification index; g, general factor; SE, standard error. Sample sizes (total n = 79): volumetric and cortical structural complexity analysis = 76, white matter microstructure analysis = 74, and cortical and deep gray matter and cerebellar microstructural diffusion analysis = 74. See Table S8 for contextualization of the image features in respect to GA at birth and at scan, Tables S9 and S10 for ASV- and species-level MaAsLin2 results, and Figure S5 for representative brain maps.

Figure 5

Taxa-level analyses correlating brain microstructural features with the relative abundances of ASVs Analyses were conducted using MaAsLin2, testing for differences in ASVs present with at least 1% of abundance in at least 10% of samples (n = 13 ASVs). ASVs are ordered by the strength of association with each brain imaging feature. Lengths of the bars correspond to the MaAsLin2 model coefficient, which relates to the strength of the association. Error bars indicate the standard error (SE) of the model coefficient. Full color bars and asterisks (∗) indicate Benjamini-Hochberg (BH) method-adjusted p value <0.25. Red indicates positive and blue negative associations. Sample size: white matter microstructure analysis = 74, and cortical and deep gray matter microstructure analysis = 74. ASV, amplicon sequence variant; MaAsLin, Microbiome Multivariate Association with Linear Models; FA, fractional anisotropy; RD, radial diffusivity; NDI, neurite density index; ODI, orientation dispersion index; cGM, cortical gray matter; dGM, deep gray matter; g, general factor; SE, standard error. See Tables S9 and S10 for ASV- and species-level MaAsLin2 results.

Figure 6

Gut-brain modules in association with brain microstructure in preterm infants (A) Mean relative abundance of the GBMs reflecting functional potential of the metagenome; bars are colored by the prevalence of the modules. 42 out of 56 GBMs were detected in the microbiome-MRI matching dataset (n = 77); all are present in at least two samples. (B) GBM correlations with the first four orthogonal principal coordinates (PCo-s) calculated from 16S rRNA beta diversity data, showing the top 20 strongest correlations for each PCo. The percentage refers to the variance explained by each of the PCo. Red indicates positive and blue negative correlations between the PCo-s and GBMs. (C) GBMs in correlation with brain microstructural features. Analyses were conducted using MaAsLin2, testing for differences in GBMs present in at least 10% of samples (n = 34 modules). Modules are ordered (left to right) by the prevalence in the dataset. Color corresponds to MaAsLin2 model coefficient, which relates to the strength of the association, with blue indicating negative and red positive correlations. Asterisks (∗) indicate Benjamini-Hochberg (BH) method-adjusted p value <0.25. Sample size n = 77 (linked MRI and metagenomic shotgun data). FA, fractional anisotropy; RD, radial diffusivity; NDI, neurite density index; ODI, orientation dispersion index; cGM, cortical gray matter; dGM, deep gray matter; g, general factor; GBM, gut-brain module. See Table S11 for MaAsLin2 results and Figure S4 for species contribution to GBMs.