Skip to main content
Download PDF

The fungal microbiota modulate neonatal oxygen-induced lung injury

Microbiome volume 13, Article number: 24 (2025) Cite this article

Abstract

Background

The immature lungs of very preterm infants are exposed to supraphysiologic oxygen, contributing to bronchopulmonary dysplasia (BPD), a chronic lung disease that is the most common morbidity of prematurity. While the microbiota significantly influences neonatal health, the relationship between the intestinal microbiome, particularly micro-eukaryotic members such as fungi and yeast, and lung injury severity in newborns remains unknown.

Results

Here, we show that the fungal microbiota modulates hyperoxia-induced lung injury severity in very low birth weight premature infants and preclinical pseudohumanized and altered fungal colonization mouse models. Instead of fungal communities dominated by Candida and Saccharomyces, the first stool microbiomes of infants who developed BPD had less interconnected community architectures with a greater diversity of rarer fungi. After using a pseudohumanized model to show that transfer to the neonatal microbiome from infants with BPD increased the severity of lung injury, we used gain and loss of function approaches to demonstrate that modulating the extent of initial neonatal fungal colonization affected the extent of BPD-like lung injury in mice. We also identified alterations in the murine intestinal microbiome and transcriptome associated with augmented lung injury.

Conclusions

These findings demonstrate that features of the initial intestinal fungal microbiome are associated with the later development of BPD in premature neonates and exert a microbiome-driven effect that is transferable and modifiable in murine models, which suggests both causality and a potential therapeutic strategy.

Video Abstract

Background

Bronchopulmonary dysplasia (BPD) is the most common morbidity of prematurity [1, 2]. Disproportionately affecting the most premature newborns, BPD exacts a devastating toll on the immature lung, impairing alveolarization and disrupting vascularization [3]. Despite the considerable, life-long morbidity and staggeringly high mortality of BPD [1, 2], meaningful interventions to modulate lung injury in preterm newborns remain limited [4]. While environmental triggers, such as exposure to hyperoxia, have been described [1, 2], modulators of the disease remain elusive.

The gut-lung axis is a recently described phenomenon whereby intestinal microbes may project their influence on lung pathophysiology [5, 6]. Examples of this connection are described in models of neonatal pneumonia [7, 8] and asthma [9, 10]. Intestinal microbes can influence the lung directly via the production of metabolites and indirectly by informing immune mechanisms that modulate lung injury and recovery [7, 9]. While previous work has focused on bacteria, other commensal organisms, such as fungi, may be critical immunity stimulators [9, 11]. In early life, gut fungal communities (the mycobiome) appear causally linked to pulmonary immune development [9]. This may make the mycobiome a key modulator of lung injury in BPD. However, the composition and alterations of fungal communities in BPD remain unexplored.

Therefore, we hypothesized that the composition of the intestinal mycobiome would predict BPD in very low birthweight (VLBW, < 1500 g) infants and that this effect would be transferable and modifiable in murine preclinical models of BPD. To test this hypothesis, we performed a prospective observational cohort study of VLBW preterm newborns to determine if the gut mycobiome composition predicted BPD development. To interrogate causality, we performed fecal microbiota transplant (FMT) on antibiotic-pseudohumanized mice and tested the modulation of the degree of fungal colonization in two additional murine models. Here, we show that the gut-lung axis may influence the development of BPD, that dysbiotic gut mycobiota can augment BPD-like illness via FMT, and the extent of lung injury is modifiable based on the extent of fungal colonization. These results suggest a potential therapeutic strategy for BPD.

Results

The gut mycobiota is unique in infants who will later develop BPD

To determine if the eukaryotic microbiota may influence hyperoxia-induced neonatal lung injury, we studied a cohort of very premature infants at high risk of BPD. Since we were interested in determining if differences in the intestinal microbiome preceded the later diagnosis of lung disease, we collected the first available true stool sample, which preceded the development of BPD by several weeks. In this prospective, observational cohort of 144 VLBW infants in neonatal intensive care, 102 infants produced their first non-meconium stool during the second week of life, from which we characterized the intestinal microbiota. Cohort characteristics are described in Tables 1 and S1. Of note, the data is compiled at the time of official BPD diagnosis, 36 weeks corrected gestational age, and that during the second week of life, when the stools we analyzed were produced, the infants had similar exposures. We excluded one infant for inadequate sequencing quality and one for incomplete clinical data. The remaining 42 excluded infants either did not stool or only passed a meconium stool during the second week of life. Characteristics of excluded infants or those who failed to produce a true stool during the second week of life, which did not differ from those included in the study cohort, are described in Table S2.

Table 1 Cohort characteristics and demographics
Full size table

After using source contamination modeling to identify and eliminate environmental contamination (Fig. S1–2), we first compared the intestinal microbiomes of preterm infants who later developed BPD to infants who did not develop BPD (hereafter, NoBPD). As shown in Fig. 1, the fungal intestinal microbiome (mycobiome) of BPD infants differed in community diversity, composition, and interconnectivity from NoBPD infants. We did not identify significant differences in the bacterial microbiome (bacteriome).

Fig. 1
figure 1

The fungal gut microbiota is altered in infants that later develop bronchopulmonary dysplasia. Intestinal microbiome analysis from the first true stool from 102 very low birthweight premature infants (birthweight < 1500 g). A Fungal alpha diversity of infants with moderate to severe bronchopulmonary dysplasia or death prior to 36 weeks of gestation (BPD) to those without (NoBPD) as quantified by the Shannon (p = 0.0010, Mann–Whitney) and Simpson indices (p = 0.0021, Mann–Whitney). Values represent means ± SEM. B Bacterial alpha diversity. C Fungal beta diversity. Principal coordinates analysis (PCoA) of Bray–Curtis dissimilarity displaying differences in fungal beta diversity (p = 0.040, permutational multivariate analysis of ANOVA, PERMANOVA. p = 0.44, permutational multivariate analysis of dispersion, PERMDISP). Ellipses indicate the 95% confidence interval. D Biplot of principal components analysis (PCA) with fungal community composition driven by amplicon sequence variants (ASVs) aligning to the genera Candida and Aureobasidium. PCA is shown in the inset. E Bacterial beta diversity. PCoA of Bray–Curtis dissimilarity displaying bacterial beta diversity (p = 0.06109, PERMANOVA. p = 0.37, PERMDISP). F Biplot of PCA of bacterial community composition. PCA is shown in the inset. G and H SPIEC-EASI network analysis shows that BPD infants formed sparser multikingdom interaction networks than NoBPD infants (fewer fungal than bacterial edges but a similar number of nodes, and lower edge density). See also Fig. S2–3

Full size image

For alpha diversity, the within-group variation exceeded the between-group differences for both the mycobiome and the bacteriome. The alpha diversity of the mycobiome was greater in BPD as quantified by the Shannon and Simpson diversity indices (p = 0.0010 and 0.0021, Mann–Whitney, Fig. 1A). The alpha diversity of the bacteriome was not different (Fig. 1B).

We used several complementary methods to quantify beta diversity. As visualized by principal coordinates analysis of Bray–Curtis dissimilarity matrices, disease status produced divergent clustering of the mycobiome (Fig. 1C). The mycobiome is differently composed in BPD infants (p = 0.04, R2 = 0.03 permutational multivariate ANOVA (PERMANOVA)). We verified this was not a false positive due to unequal dispersion with permutational multivariate analysis of dispersion (PERMDISP, p = 0.44). To identify microbial drivers of community composition, we performed biplot analysis of principal components analysis (PCA) to identify which amplicon sequence variants (ASVs) produced these compositional differences (Fig. 1D, inset on biplot shows PCA). We determined that multiple ASVs in the genera Candida and an ASV aligning to Aureobasidium made the most significant contributions to the structure of the mycobiome. However, the bacteriome was not different (p = 0.061, R2 = 0.017 PERMANOVA, p = 0.37, PERMDISP, Fig. 1D). Biplot analysis showed Enterococcus and multiple Klebsiella ASVs contributed to the overall bacterial community structure (Fig. 1F) but are not differently distributed by disease status.

To quantify the interconnectivity of these microbial communities, we used a novel cross-domain extension of SPIEC-EASI[11] to build multikingdom interaction networks (Fig. 1G–H, Fig. S3g–h). Compared to NoBPD infants, BPD infants have fewer fungal than bacterial edges but a similar number of nodes (Fig. S3e). Network analysis indicated that the networks of BPD infants had a lower edge density, a measure of sparsity (absolute difference, 0.017). Together, this indicates that the multikingdom interaction networks are sparser for BPD infants, with fewer fungi-to-fungi and fungi-to-bacterial interactions than NoBPD infants.

We also used negative binomial regression (DESeq2) to identify differentially abundant genera. For the mycobiome, the genus Candida was more abundant in NoBPD infants, while Mortierella was increased in BPD infants (Fig. S3a–b). For the bacteriome, we detected modest increases in the family Staphylococcaceae and in the genus Veillonella (Fig. S3c–d), similar to our earlier observations in neonatal mice [12]. We also built random forest machine learning models using both the fungal and microbiota data and found that the abundance of fungal taxa predicted the later development of BPD while the bacterial taxa did not (auROC 0.74 versus 0.48, Fig. S5).

Altogether, these analyses indicate the gut mycobiome of infants who will later develop BPD is unique several weeks to months before diagnosis. We identified global differences in alpha and beta diversity and found these communities were characterized by a lower abundance of Candida and a higher number of rarer fungal taxa, including Mortierella, arranged in more sparse and less interconnected multikingdom networks. Comparable differences in alpha and beta diversity were not observed in the gut bacteriome or the combined multikingdom dataset (Fig. S4).

Sex-specific differences in the composition of the mycobiome

Having identified that differences in the intestinal mycobiome preceded the later development of BPD, we examined if sex-specific features exist in the neonatal mycobiome (Fig. 2 and Fig. S6–7). We had predetermined to test if sex-specific differences existed in the microbiome because sex-specific differences in BPD pathogenesis are known to produce worse outcomes in male infants [13, 14]. When categorized by sex and disease state, the Shannon and Simpson diversity indices of the mycobiome and bacteriome did not differ between infants assigned female (AFAB) or male at birth (AMAB, Fig. 2A–B). Fungal community composition differed for all groups (p = 0.041, R2 = 0.071, PERMANOVA; p = 0.60, PERMDISP, Fig. 2C). AFAB infants with BPD and NoBPD clustered closer than infants AMAB (Fig. 2C). DESeq2 identified differences in the abundance of the fungal genera Mortierella and Pseudeurotium in AMAB infants who later developed BPD (Fig. S7a–b). SPIEC-EASI network analysis showed greater interconnectivity of infants with NoBPD regardless of sex assigned at birth (Fig. 2E–H, Fig. S7e–f). Infants AMAB had fewer fungal edges than those AFAB, but the number of edges was reduced only in infants AMAB who would later develop BPD (Fig. S7c). Infants AMAB had a lower edge density (absolute difference, 0.004). These differences produced sparser networks in AMAB infants, with the least interconnected networks in AMAB who later developed BPD. From these analyses, we conclude that sex-specific differences in the composition of the mycobiome are the strongest in infants AMAB who will develop BPD.

Fig. 2
figure 2

Sex-specific differences in premature infants that will develop BPD. Pre-determined analysis of 102 very low birthweight premature infants. A Fungal alpha diversity in infants in infants assigned male (AMAB) or female at birth (AFAB) with NoBPD or BPD, as quantified by the Shannon and Simpson indices. Values represent means ± SEM. Significance testing by two-way Mann–Whitney. B Bacterial alpha diversity. Significance testing by two-way Mann–Whitney. C Fungal beta diversity. Principal coordinates analysis (PCoA) of Bray–Curtis dissimilarity showing differences in the beta diversity of the mycobiome for each sex and disease state (p = 0.041, PERMANOVA, p = 0.5957, PERMDISP). Ellipses indicate the 95% confidence interval. D Bacterial beta diversity. PCoA of Bray–Curtis dissimilarity of the bacterial microbiome (p = 0.32, PERMANOVA. P = 0.64, PERMDISP). E–H SPIEC-EASI network analysis shows that BPD infants AMAB formed the sparsest networks (fewer edges, but a similar number of nodes, lower edge density). See also Fig. S6–7

Full size image

Microbiome-driven clustering identifies fungal clusters linked to clinical outcomes

Having shown that fungal microbiota features are associated with BPD development, we sought to identify other clinical characteristics that may contribute to heterogeneity in the composition of the neonatal mycobiome. To do this, we used microbiome-driven clustering to explore heterogeneity in fungal community features using Dirichlet-Multinomial Mixtures (DMM) modeling [15] (Fig. 3 and Fig. S8). In contrast to our previous approach, this orthogonal method identifies clusters based on their microbial composition, which can then be associated with clinical characteristics. This permits the identification of microbial heterogeneity within disease states that might otherwise remain undetected.

Fig. 3
figure 3

Microbiome-driven clustering analysis of heterogeneity in the mycobiome. A Principal Coordinates Analysis (PCoA) plot showing a five-cluster Dirichlet-Multinomial Mixture (DMM) model clustering of the mycobiome, n = 102. In Cluster 1, 3 infants developed BPD, and 18 did not (NoBPD). In Cluster 2: 6 BPD, 7 NoBPD, Cluster 3: 9 BPD and 3 NoBPD, Cluster 4: 5 BPD and 6 NoBPD, Cluster 5: 1 BPD and 5 NoBPD. B Fungal community diversity analyses indicate the greatest microbial richness in Cluster 3, as quantified by the Shannon diversity index. Values represent means ± SEM. Significance testing by ANOVA. C Average relative abundance of fungi in each DMM cluster. D–F A two-cluster DMM model in NoBPD infants, with the greatest richness in Cluster 2, n = 64. G–I A two-cluster DMM model in BPD infants with the highest microbial richness in Cluster 2, n = 38. Clust1–5, Cluster 1–5. Ellipses indicate the 95% confidence interval. See also Fig. S8–10

Full size image

For this study, Laplace analysis for the whole mycobiome identified 5 clusters as the optimal fit point (Fig. S8a). Membership in Cluster 1 of the fungal DMM model (Fig. 3A–C) was driven by the relative abundance of an Aureobasidium ASV and a Saccharomyces ASV. A Candida ASV drove membership in Cluster 2. Cluster 3 was driven primarily by an Aureobasidium ASV. Cluster 4 was driven by Candida and Aureobasidium ASVs and, to a lesser extent, by an Issatchenkia ASV. Finally, Cluster 5 was almost entirely driven by a Candida ASV. We then used canonical correspondence analysis (CCA) to map clinical characteristics that may be associated with the fungal composition profiles identified by the DMM components (Fig. S8c). Clusters 1 and 3 shared similar vectors and were closely associated with BPD and increased illness severity (Critical Risk Index for Babies II, CRIB-II). Cluster 2 was in opposition and, therefore, associated with opposing markers of illness severity, such as increasing gestational age, birth weight, and female sex-assigned at birth (AFAB). The limited number of samples in Cluster 5 was at right angles to most of these clinical variables. It was associated with low fungal evenness (arbitrarily defined one ASV contributing greater than half of the relative abundance).

Laplace analyses performed independently within the BPD or NoBPD cohorts identified 2 clusters as the best fit for each of these mycobiome subsets (Fig. S8d and g). In the NoBPD DMM model (Fig. 3D–F), membership in Cluster 1 was driven by Saccharomyces, Aureobasidium, and, to a lesser extent, a Candida ASV (Fig. S8h). Membership in Cluster 2 was driven by three different Candida ASVs and a Saccharomyces ASV. CCA suggests that membership in Cluster 1 is associated with the opposite clinical variables as NoBPD Cluster 1 (BW, GA, and AMAB), except for fungal abundance (Fig. S7I). In the BPD DMM model (Fig. 3G–I), membership in Cluster 1 was driven by Saccharomyces and Candida ASVs (Fig. S8e). An Aureobasidium ASV primarily drove membership in Cluster 2. CCA closely associated Cluster 1 with birthweight and, to a lesser extent, with gestational age and death (Fig. S8f). Together, these analyses identify microbial taxa that influence beta diversity and identify clinical covariates, such as illness severity, birthweight, and gestational age, that introduce heterogeneity in disease state-driven clustering. DMM modeling of the bacteriome is presented in the supplement (Fig. S9–10).

Fecal microbiota transfer transmits a severe phenotype of BPD

To test if alterations in the gut mycobiome that preceded the development of BPD in premature infants could augment lung injury, we transferred the microbiome from infants with either BPD or NoBPD into mice. Because the effects of bacterial and fungal members of the microbiome are interdependent and we could not transfer only fungal or bacterial communities without disrupting their interrelated ecology and function, we transferred intact microbial communities. Despite this consideration, this experiment tests the novel hypothesis that the intestinal microbiome can transfer an effect on lung injury from human newborns to mice. We used a pseudohumanized approach with pregestational colonization, rather than using germ-free mice, in order to avoid side effects on prenatal immune development [16]. Following antibiotic preparation, recipient female mice were orally inoculated with anoxic fecal samples from donor infants. These mice were bred, and the resulting newborn mice were exposed to hyperoxia to induce BPD-like lung injury in the standardized model of BPD[17] (“Methods” section, Fig. 4A).

Fig. 4
figure 4

Fecal microbiota transfer from infants with BPD augments lung injury in mice. A Experimental schematic showing pregestational fecal microbiota transfer to colonize neonatal mice with the intestinal microbiota from human infants with BPD or without (NoBPD). Hyperoxia (HO) or normoxia (NO) exposure was performed from day 3 to day 14 of life (P3–14). Tissue harvest and subsequent analyses were performed at P14. Schematic created using BioRender.com. Total sample size is 18–24 neonatal mice from three independent litters per experimental group. B Gut fungal microbiome analysis showing alpha diversity as quantified by the Shannon diversity index and a principal coordinates analysis (PCoA) plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 12 mice/group). C Gut bacterial microbiome analysis showing alpha diversity as quantified by the Shannon diversity index and a principal coordinates analysis (PCoA) plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 12 mice/group). D Lung fungal microbiome analysis showing alpha diversity as quantified by the Shannon diversity index and a principal coordinates analysis (PCoA) plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 12 mice/group). E Lung bacterial microbiome analysis showing alpha diversity as quantified by the Shannon diversity index and a principal coordinates analysis (PCoA) plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 12 mice/group). F Representative hematoxylin and eosin confocal micrographs of the distal lung of 14-day-old neonatal mice born to pseudohumanized mice with fecal microbiota transport from BPD-associated microbiomes (BPD-FMT) or NoBPD-associated microbiomes (NoBPD-FMT). Scale bars, 100 µm. G Histomorphological analysis and forced oscillometry. Morphometry is quantified by mean linear intercept and radial alveolar count. Forced oscillometry is quantified by resistance and compliance. Significance testing by two-way ANOVA. H Representative confocal microscopy images of immunofluorescence of alpha-smooth muscle actin and von Willebrand’s factor showing vascular changes consistent with early pulmonary hypertension and reduced vascular density of 20–50 µm vessels. Scale bars, 20 and 50 µm. Significance testing by two-way ANOVA. I–K Hyperoxia exposure induces different transcriptomic profiles in the lungs in BPD-FMT or NoBPD-FMT mice by RNA-seq (n = 4 mice/group). I Volcano plots of the hyperoxia-induced transcriptomic alterations in BPD-FMT or NoBPD-FMT mice. J Venn diagram showing differential gene expression in BPD-FMT and NoBPD-FMT as compared to SPF controls. K Heatmaps of the top 10 most differentially expressed genes, either unique BPD-FMT, shared or unique to NoBPD. FiO2, fraction of inspired oxygen. Values represent means ± SEM. See also Fig. S11–12

Full size image

In the neonatal mouse gut, where we sampled the microbiome from the intact terminal ileum due to our prior demonstration of its biological relevance in BPD [12], beta but not alpha diversity of both the fungal (p = 0.046, R2 = 0.094 PERMANOVA, Fig. 4B) and bacterial (p = 0.0461, R2 = 0.090, PERMANOVA, Fig. 4C) was significantly different in mice that received stool from infants with BPD (BMD-FMT) than those from infants without BPD (NoBPD-FMT). Their respective lung microbiomes were not altered (Fig. 4D–E). The FMT-transferred microbiome also differed from the control conventionally raised specific pathogen-free microbiome (SPF, Fig. S11).

Neonatal BPD-FMT mice had more marked histological lung injury when exposed to hyperoxia than NoBPD-FMT mice (Larger mean linear intercept (MLI) and lower radial alveolar count (RAC), p < 0.0001 and p = 0.0005, two-way ANOVA, Fig. 4F–G). Functionally, BPD-FMT mice had increased pulmonary resistance (p = 0.0015) and decreased compliance (p = 0.020). We also used immunofluorescence and morphometry to measure neointimal changes consistent with pulmonary hypertension and identified BPD-FMT mice had reduced vascular density of 20–50 µm diameter vessels (p = 0.0065, Fig. 4H).

Next, we used RNA-seq to quantify transcriptomic changes in the lungs of these mice. Analysis confirmed unique global transcriptomic changes in the lungs of mice by both FMT and hyperoxia exposure (Fig. 4I–K). Volcano plots show a clear difference in the response of NoBPD-FMT and BPD-FMT (Fig. 4I) to hyperoxia exposure. To identify genes uniquely regulated in mice exposed to FMT, we compared their lung transcriptome to hyperoxia or normoxia-exposed conventional mice (SPF controls without any microbiome alterations, Fig. 4J). We identified 210 differentially expressed genes (DEGs) in BPD-FMT mice and 363 genes in NoBPD-FMT mice as compared to 793 DEGs in the conventional SPF mice. Pathway analysis identified upregulated cell death, interleukin(IL)-4 and IL-13 signaling, and extracellular matrix organization in BPD-FMT mice (Fig. S12). In contrast, NoBPD-FMT mice upregulated oxidation, Cyp2e1 reactions, and glutathione conjugation. Both groups upregulated the Keap1-Nfe2l2 cellular stress pathway but to a greater extent in BPD-FMT mice. Ingenuity pathway analysis uncovered increased acute phase response and Nrf2-mediated oxidative stress response in BPD-FMT but not in NoBPD-FMT mice (Fig. S13), with corresponding decreases in the dilated cardiomyopathy signaling pathway, LXR/RXR activation, and CREB signaling. In contrast, SPF controls had upregulated pathogen-activated cytokine storm, angiogenic, IL-6, and IL-17 signaling that was not observed in either FMT group. From this experiment, we conclude that a BPD-associated gut microbiome is transferable from human newborns to neonatal mice and augments lung injury during hyperoxia exposure.

Inhibition of perinatal fungal colonization reduces lung injury

We have previously shown that mice without a microbiome (germ-free) are protected from hyperoxia-induced lung injury [18], so next, we used a loss-of-function approach to test if limiting fungal colonization reduced lung injury. To do this, we exposed pregnant mice to the antifungal fluconazole prior to delivery (Fig. 5A). The gut mycobiomes of their offspring were indeed suppressed (Fig. S14), and we did not observe antifungal-driven augmentation of selected fungal species (Fig. S15). Fluconazole-induced suppression of the mycobiome also produced clear differences in gut beta diversity in samples where fungi were detectable (Fig. 5B–C). We did not observe a difference in the lung microbiome (Fig. 5D–E).

Fig. 5
figure 5

Inhibiting fungal colonization reduces hyperoxia-induced lung injury in mice. A Experimental schematic showing prenatal fluconazole exposure to inhibit neonatal fungal colonization. Hyperoxia (HO) or normoxia (NO) exposure was performed from day 3 to day 14 of life (P3–14). Tissue harvest and subsequent analyses were performed at P14. Schematic created using BioRender.com. Total sample size is 18–24 neonatal mice from three independent litters per experimental group. B Gut fungal microbiome analysis, showing alpha diversity as quantified by the Shannon diversity index and a PCoA plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney and PERMANOVA (n = 8 mice/group). C Gut bacterial microbiome analysis, showing alpha diversity as quantified by the Shannon diversity index and a principal coordinates analysis (PCoA) plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney and PERMANOVA (n = 8 mice/group). D Lung fungal microbiome analysis, showing alpha diversity as quantified by the Shannon diversity index and a PCoA plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 8 mice/group). E Lung bacterial microbiome analysis, showing alpha diversity as quantified by the Shannon diversity index and a PCoA plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 8 mice/group). F Representative hematoxylin and eosin confocal micrographs of the distal lung of 14-day-old neonatal mice born to mouse dams exposed to fluconazole or conventional specific pathogen-free (SPF), vehicle-exposed, control mice. Scale bars, 100 µm. Histomorphological analysis and forced oscillometry. Morphometry is quantified by mean linear intercept and radial alveolar count. Forced oscillometry is quantified by resistance and compliance. Significance testing by two-way ANOVA. G Representative immunofluorescence of alpha-smooth muscle actin. Quantification of vascular density of 20–50 µm vessels, and echocardiography quantification of right ventricular systolic pressure. Scale bars, 100 µm. Significance testing by two-way ANOVA. H Flow cytometric quantification of group 2 innate lymphoid cells (ILC2), expressed as a percentage of CD45+ cells. ILC2 were defined as live, lineage negative, CD45+, CD90.2+, CD127+, IL-33R (ST2)+, and GATA3+ cells. Significance testing by two-way ANOVA (n = 5–7 mice/group). I–K Hyperoxia exposure alters the lung transcriptome in the lungs in fluconazole-exposed or SPF mice by RNA-seq (n = 4 mice/group). I Volcano plot of the hyperoxia-induced transcriptomic alterations in fluconazole-exposed mice and SPF controls. J Venn diagram showing differential gene expression in fluconazole-exposed neonatal mice compared to SPF controls. K Heatmaps of the top 10 most differentially expressed genes either unique to SPF controls, shared, or unique to fluconazole-exposed mice. FiO2, fraction of inspired oxygen. Values represent means ± SEM. See also Fig. S14–17, Table S4

Full size image

Reduced fungal colonization was associated with improvement in both lung structural maturity after hyperoxia exposure (MLI p < 0.001, two-way ANOVA; RAC p < 0.001) and pulmonary mechanics (resistance p = 0.0026, compliance p = 0.0005, Fig. 5F). Using echocardiography, we also detected a functionally meaningful decrease in right ventricular systolic pressure, which is the most reliable measure of pulmonary hypertension in neonatal mice (p = 0.011, Fig. 5G). Flow cytometry showed that group 2 innate lymphoid cell (ILC2) populations were increased in the lung (p < 0.001, Fig. 5H, S16).

We performed RNA-seq on the lungs of these mice to gain insight into how suppressing fungal colonization might influence their response to lung injury. Volcano plots showed a clear difference in the response to hyperoxia exposure (Fig. 5I) that resulted from the differential expression of 759 genes (DEGs) in SPF controls and 944 DEGs in fluconazole mice (Fig. 5J). Heat maps of the top 10 most differentially regulated genes unique to either SPF controls, fluconazole mice, or shared between these groups (Fig. 5K). The most increased genes in SPF mice notably included matrix metallopeptidases, Inhba (which encodes for a TGF-β family member) Il6 (interleukin-6), and Fosl1 (which is an upstream regulator of several developmental pathways including vasculogenesis). In fluconazole mice, the most increased genes notably included several metabolism-related genes (Adig, Agrp, Trib3, Cox6a3), transcription regulation (Nupr1), and TGF-β signaling (Trib3 and Wfikkn1). Pathway analysis identified that the most upregulated disease ontology terms in hyperoxia-exposure SPF controls were obstructive and interstitial lung disease, arteriosclerosis and pulmonary fibrosis (Fig. S17). We also found several upregulated cytokine signaling-related pathways, in particular IL-4, IL-13, and IL-17 signaling, extracellular matrix integrity, and HIF1 (hypoxia-inducible factor-1) signaling pathways. In contrast, fluconazole upregulated pathways were related to mitochondrial metabolism, reactive oxygen species, xenobiotic metabolism, and transcription regulation. From this experiment, we conclude that inhibiting initial fungal colonization reduces the extent of hyperoxia-induced lung injury and improves lung development and mechanics.

Augmentation of perinatal fungal colonization amplifies BPD severity

Antibiotic-facilitated fungal colonization is an established technique to augment fungal colonization in mice [19]. Therefore, to understand if increased fungal colonization worsened disease severity, we used antibiotic-assisted fungal colonization to amplify the colonization of the mouse fungal commensal Candida tropicalis (Ctrop, Fig. 6A). We did not use Candia albicans, because it is not a natural mouse commensal like it is in humans. Indeed, the average number of fungal reads was markedly increased (230-fold). As intended, this was driven by an increase in the relative abundance of Candida ASVs in the small intestines of these mice (Fig. S18). This was accompanied by alterations in both the lung and gut microbiome (Fig. 6B–E).

Fig. 6
figure 6

Amplifying fungal colonization augments hyperoxia-induced lung injury in mice. A Experimental schematic showing cefoperazone-facilitated intestinal colonization with the mouse commensal yeast, Candida tropicalis (C. trop). Hyperoxia (HO) or normoxia (NO) exposure was performed from day 3 to day 14 of life (P3-14). Tissue harvest and subsequent analyses were performed at P14. Total sample size is 18–24 neonatal mice from three independent litters per experimental group. B Gut fungal microbiome analysis, showing alpha diversity as quantified by the Shannon diversity index and a principal coordinates analysis (PCoA) plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 8 mice/group). C Gut bacterial microbiome analysis, showing alpha diversity as quantified by the Shannon diversity index and a principal coordinates analysis (PCoA) plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 8 mice/group). D Lung fungal microbiome analysis, showing alpha diversity as quantified by the Shannon diversity index and a PCoA plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 8 mice/group). E Lung bacterial microbiome analysis, showing alpha diversity as quantified by the Shannon diversity index and a PCoA plot showing Bray–Curtis dissimilarity. Significance testing by Mann–Whitney or PERMANOVA (n = 8 mice/group). F Representative hematoxylin and eosin confocal micrographs of the distal lung of 14-day-old neonatal mice born to mouse dams colonized by C. trop or conventional specific pathogen-free (SPF), vehicle-exposed, control mice. Scale bars, 100 µm. Histomorphological analysis and forced oscillometry. Morphometry is quantified by mean linear intercept, radial alveolar count. Forced oscillometry is quantified by resistance and compliance. Significance testing by two-way ANOVA. G Representative immunofluorescence of alpha-smooth muscle actin. Quantification shows vascular density of 20–50 µm vessels and echocardiography quantification of right ventricular systolic pressure. Scale bars, 100 µm. Significance testing by two-way ANOVA. H Flow cytometric quantification of group 2 innate lymphoid cells (ILC2), expressed as a percentage of CD45+ cells. ILC2 were defined as live, lineage negative, CD45+, CD90.2+, CD127+, IL-33R (ST2)+, and GATA3+ cells. Significance testing by two-way ANOVA (n = 5–7 mice/group). I–K Hyperoxia exposure alters the lung transcriptome in the lungs in fluconazole-exposed or SPF mice by RNA-seq (n = 4 mice/group). I Volcano plots of the hyperoxia-induced transcriptomic alterations in Ctrop mice and specific pathogen-free (SPF) controls. J Venn diagram showing unique gene regulation in Ctrop mice compared to SPF controls. K Heatmaps of the top 10 most differentially expressed genes either unique to SPF controls, shared, or unique to Ctrop mice. FiO2, fraction of inspired oxygen. Schematic created using BioRender.com. Sample size is 18–24 neonatal mice from three independent litters per experimental group. Values represent means ± SEM. See also Fig. S18–19 and Table S4

Full size image

Augmented fungal exposure amplified structural and functional lung injury (MLI p = 0.0037, two-way ANOVA; RAC p = 0.0002, resistance p = 0.0026, and compliance p = 0.0005, Fig. 6F). We also observed a decrease in the number of vessels under 50 microns in diameter (p < 0.0001, Fig. 6G) in Ctrop mice but without a clear increase in right ventricular systolic pressure by echocardiography. Finally, flow cytometry also showed an increase in ILC2s (p < 0.0001, Fig. 6H), which is the opposite of the fluconazole mice.

Overall, RNA-seq revealed more modest differences between SPF controls and Ctrop mice (Fig. 6I) than between fluconazole mice and controls, with a greater number of overlapping DEGs (483, Fig. 6J). However, we still identified 340 unique DEGs in Ctrop mice. The most differentially expressed genes included several C–C motif chemokine ligand-related genes (Ccl3, Ccl4, Ccl18, Xcl1, and Xcl2) and granzyme B (Gzmb, Fig. 6K). Pathway analysis revealed increases in several cytokine signaling pathways in Ctrop mice, including Toll-like receptor signaling, IL-10 signaling, and neutrophil degranulation (Fig. S16). From this experiment, we conclude that augmenting fungal colonization can increase the severity of lung injury. Altogether, this suggests that modulating fungal colonization modulates the extent of neonatal hyperoxia-induced lung injury.

Discussion

Here, we show the fungal members of the gut microbiome contribute to the severity of BPD. In a cohort of premature infants, using an unbiased approach that considered both bacterial and microeukaryotic members of the intestinal microbiome, we sought to identify prognostic features in the multi-kingdom microbiome that could determine which infants would later develop BPD. Because we aimed to identify early prognostic features, we selected samples as close to birth as possible, starting with the first true non-meconium stool samples from these infants. We found that the gut mycobiome, but not the bacterial microbiome, differed several weeks before diagnosis from neonates who later developed BPD from those who did not develop BPD (NoBPD).

In adults, fungi are vital members of the human microbiota [20]. Unlike bacterial commensals, fungal non-pathological and non-parasitic functions remain poorly understood [21], particularly in newborns, primarily because of technical limitations. Motivated by these differences in the mycobiome that precede the development of BPD in premature infants, we turned to preclinical models of BPD to test if the composition of the mycobiome was linked to the severity of lung injury. However, fungi are integrated members of the microbiome and are unlikely influence the host completely independently of bacteria. So first, we used fecal microbiota transplant to show that the transfer of the intact multikingdom gut microbiome from human infants with BPD into antibiotic-prepared mice worsened the extent of neonatal hyperoxia-induced lung injury as compared to mice that received their microbiota from infants without BPD. Guided by our previous observations that animals without a microbiome are protected from hyperoxia-induced lung injury [18, 22], but mice with an antibiotic-suppressed microbiome have augmented lung injury [23], as well as our exploration of the impact of neonatal hyperoxia on the neonatal small intestine[12], next we used a loss-of-function approach and depleted the perinatal mycobiome by exposing the maternal mice to fluconazole. These mice were also partially protected from lung injury. Finally, we used a gain-of-function approach to show that lung injury was increased by augmenting fungal colonization with the mouse commensal C. tropicalis. Our findings suggest the existence of a gut-lung axis, wherein altering the extent of fungal colonization modulates the extent of lung injury.

While the microbiome has recently emerged as a contributor to clinical heterogeneity in adult patients in intensive care [22, 24,25,26,27,28,29,30,31], our knowledge of the mycobiome in newborn respiratory disease is limited [32]. Here, in one of the most extensive studies of the multikingdom intestinal microbiome in preterm newborns, we show that the composition of gut mycobiome during the second week of life predicts the development of BPD. We also identified robust sex-specific differences between infants who would develop BPD that were not observed in infants who did not develop BPD. In light of recent sex-specific differences in the development of BPD detected by single-cell transcriptomics [13, 33,34,35], this is an interesting finding worth further investigation.

Collectively, our analyses demonstrate that the composition of the intestinal mycobiome of infants who did not develop BPD was more uniform. In contrast, those who eventually developed BPD had more disparate mycobiomes. This finding suggests that a particular pattern of mycobiome development may be necessary to impart resistance to the development of BPD, and failure to do so in various ways is associated with disease development. A more robust health-associated microbiome with a less specific disease-associated microbiome is a well-described ecological pattern termed the Anna Karenina phenomenon [36]. This paradigm holds for prior genomic, proteomic, and transcriptomic analysis of the pathogenesis of BPD [37] and, therefore, fits within our broader understanding of the disease.

A previous study of the bacterial microbiota identified a divergence over time between the bacterial community composition in the gut and lungs of preterm infants that would eventually develop BPD [38], suggesting that a gut-lung axis may influence the development of BPD. However, while we could not verify serial divergence given our study design, we also did not detect any significant differences in the bacteriome in the infants we studied. While the mechanisms underlying the gut-lung axis require further exploration in BPD, work in other childhood diseases has suggested some potential mechanisms. Mouse models of childhood asthma suggest intestinal fungi are more immunostimulatory than bacteria, while intestinal bacteria are more responsible for producing metabolites detectable in the host serum [9]. Gut intestinal fungi, mainly the yeast C. albicans, are known to provoke strong immune responses in the host [39]. Indeed, monocolonization with Candida can stimulate postnatal immune development [19]. The greater abundance of Candida-dominated communities we observed in NoBPD infants is intriguing in this regard, as it may suggest a certain level of fungal immune stimulation opposes the development of BPD.

The primary limitation of this study was that we could not simultaneously sample both the lung and gut microbiomes of these infants. The study design may have also been strengthened by longitudinal sampling. However, in this study, we aimed to sample the gut microbiota as early as possible after birth to identify early prognostic markers of BPD. Future experiments may benefit from an integrated multiomic approach to identify and validate potential candidate microbes or metabolites. Similarly, it is currently not feasible to confine the transfer of microbes to only one organ system. Quantifying how spillover from FMT might affect other host organs, such as the lungs, remains challenging. A strength of this study was the rigorous exploration of the multikingdom microbiota using the most up-to-date bioinformatic techniques [32] and the demonstration that the mycobiome effect was transferable by FMT and modifiable in other mouse models.

In conclusion, in VLBW preterm newborns, we found that baseline compositional differences in the gut mycobiome predicted the eventual development of BPD. Furthermore, we showed that this dysbiotic mycobiota could augment lung injury when transferred into newborn mice and that modifying fungal colonization modulated the extent of lung injury in preclinical mouse models of BPD. These findings suggest that the gut mycobiota may represent a therapeutic target in newborn lung disease.

Methods

Materials availability

This study did not generate unique reagents or transgenic animals.

Data availability

  • Raw RNA-seq, 16S, and ITS2 sequencing data are deposited at the NCBI Sequence Read Archive (Access Number: PRJNA982880).

  • Processed data files and original code are available at github.com/WillisLungLab/MiBPD.

Experimental model and subject details

Study subjects

University of Tennessee Health Science Center cohort

We performed a prospective observational cohort study at the regional NICU at Regional One Health (Memphis, TN, USA) from July 2017 to January 2020, pre-registered with clinicaltrials.gov (NCT03229967). The Microbial Influence on the Development of Bronchopulmonary Dysplasia (MiBPD) study approached the mothers of all inborn very low birthweight infants (VLBW, < 1500 g) without known major congenital anomalies or immunodeficiency for potential consent for enrollment during the first week of life unless they were enrolled in a concurrent trial of probiotics that enrolled < 10 infants during this period. Mothers provided written consent on their own behalf and assent for the inclusion of their newborn before enrollment. Samples were collected weekly for the first month of life. Meconium samples collected within the first 72 h have been described [40]. Analysis in this study was constrained to fecal samples collected during the second week of life. Demographic and clinical data were collected to describe the clinical course and to identify the severity of bronchopulmonary dysplasia (BPD) according to the NICHD workshop definition [41]. The risk of developing BPD was predicted using the NICHD Neonatal BPD Outcomes Estimator (neonatal.rti.org, 2022 updated formula) [42]. Human subjects data is reported in accordance with the extended STROBE guidelines [43].

University of Alabama at Birmingham cohort

Samples from this cohort were collected in oxygen-reduced PBS immediately prior to FMT to the mice. The samples were used for the generation of the pseudohumanized mouse model and were collected from August 2021 to June 2022.

Animals

An 8-week-old female C57Bl/6 J mice (Stock Number 000664) were purchased from the Jackson Laboratory and allowed to acclimate in SPF housing with 4 animals per cage, on a 12-h light/dark cycle, at 20–24 °C, with ad libitum access to irradiated rodent chow (Teklad 7904, Envigo) in the vivarium at the UAB for 2 weeks. Animal studies were performed in accordance with the recommendations in the Guide for the Care and Use of Laboratory Animals of the National Institutes of Health under protocol 22,042, approved by the Institutional Animal Care and Use Committee at UAB. Animal experimentation is reported in accordance with the ARRIVE guidelines [44].

Generation of pseudohumanized mice

Pseudohumanized mice were generated by fecal microbiota transplant (FMT) from human infant donors with or without severe BPD, following antibiotic preconditioning, as previously described [45,46,47], generating mice we named NoBPD-FMT and BPD-FMT, respectively. Briefly, 9-week-old dams were subjected to antibiotic administration for 2 weeks in their drinking water, allowing them to access ad libitum. The first week consisted of systemically absorbable antibiotics: penicillin VK, cefoperazone sodium, and clindamycin sodium (all at a concentration of 1 mg/mL). This was followed by a bowel cleansing regimen by oral-gastric gavage with 200 μL of polyethylene glycol 3350 at 425 g/L. This was repeated twice. Mice were placed on the second week of antibiotic therapy in the drinking water, consisting of the non-absorbable antibiotics ertapenem sodium, vancomycin sodium, and neomycin sodium (all at 1 mg/mL). Cage bottoms were changed daily and after each gavage of polyethylene glycol, as mice are coprophagous and could reinoculate themselves with stool. After 2 weeks of antibiotic treatment, female mice were given stool from neonatal human donors with or without BPD. The stool was preserved in oxygen-reduced sterile PBS. We then performed FMT. Mice in each group were gavaged twice, using a sample of 200 μL from the same human donor. Stool samples were collected after antibiotic therapy, FMT engraftment, or at harvest time. After engraftment, dams inoculated with the same donors FMT were bred simultaneously to generate neonatal mice for further experimentation. On day 2 after birth, the pups were inoculated with the same donor FMT to create BPD-FMT and NoBPD-FMT mice, respectively. Each experiment was conducted with at least two litters arising from the same human donor and a minimum of three human donors in each experimental group (6 litters = 3 donors/group).

Perinatal antifungal exposure model

Cages of 4 dams were randomized to receive either fluconazole (250 mg/L) or standard water (control) ad libitum, (approximately 0.04–0.5 mg/g/day), on estimated gestational day 15 (E15). The exposure was continued until the resulting offspring were euthanized on postnatal day 14 (P14).

Perinatal antibiotic-assisted fungal colonization model

Similarly, cages of 4 dams were randomized to receive either cefoperazone (500 mg/mL) or standard water (control) ad libitum on E15. Pregnant mice were then split into single cages before giving birth. At birth, cefoperazone was discontinued, and the antibiotic-preconditioned dams were gavaged with 10,000 CFU of Candida tropicalis in PBS (5). Repeat exposure was performed on P2 and P3 to both pups and dams, with an additional inoculum administered to the maternal abdomen (6).

Neonatal hyperoxia exposure model

We utilized our established, standardized hyperoxia-exposure model of BPD [17, 23]. We time-mated the mice to produce multiple simultaneous birth cohorts with the same perinatal exposure. On the day of birth, pups from two consecutively born litters with the same perinatal exposure were pooled, and the pups were redistributed evenly between the two dams. Litter sizes for all experiments were adjusted to 5–8 pups per treatment group to minimize nutritional effects on lung development. One pooled litter was randomly assigned to hyperoxia [fraction of inspired oxygen (FiO2) 0.85]) and the other to air (FiO2 0.21). Pups were exposed continuously for 11 days, starting on the third day after birth (P3-14). Dams were rotated daily between the two pooled litters to limit the complications of hyperoxia. All animals were housed under specific-pathogen-free conditions with a 12-h night-day cycle and fed Envigo 7904 diet ad libitum.

Fungal culture

C. tropicalis strain (Castellani) Berkhout, a mouse fungal commensal (5, 9), was used for in vivo antibiotic-facilitated colonization. Before oral inoculation, yeast were stab-inoculated into Sabouraud Dextrose Broth 4% (pH 5.6, Thermofisher Scientific) and grown for 48 h at 28 °C. Yeasts were enumerated and resuspended in PBS prior to oral gavage of 10,000 cfu in 0.05 mL to each mouse.

Method details

Data collection

Clinical data were collected from the medical record that described clinical illness severity, infant characteristics, maternal and infant demographics, the extent of prematurity, respiratory parameters, and known microbiome modulators when the infant reached 36 weeks’ corrected gestational age or was discharged. BPD was adjudicated at 36 weeks corrected gestational age using the NICHD workshop definition [41]. We quantified clinical illness severity with the Critical Risk Index for Babies II (CRIB II) [48]. The risk of developing BPD was predicted using the NICHD Neonatal BPD Outcomes Estimator (neonatal.rti.org, 2022 updated formula) [42].

Bacterial and fungal amplicon sequencing

For the human cohorts, the first non-transitional stool produced during the second week of life was collected in oxygen-reduced PBS with a sterile spatula immediately following evacuation into a clean diaper. For mouse samples, we collected the whole left lung and terminal ileum [49]. Microbial DNA was extracted, amplified, and sequenced as described [40, 50]. Negative and positive controls were collected and passed through the DNA extraction and sequencing steps. Sequence data were processed with QIIME 2 [51] and DADA2 [52], decontaminated using PERFect [53], and analyzed in R. Alpha diversity using Shannon and Simpson indices, and beta diversity using Bray–Curtis dissimilarity was performed using the vegan package in R. Network analysis was performed using SPIEC-EASI. Details regarding data processing and ecological analyses are available in the online supplement.

RNA sequencing

RNA sequencing was performed on isolated tissue from the left middle lobe using established methods [54]. Briefly, 25 mg of tissue was collected in RLT buffer (Qiagen) with 1% (v/v) 2-Mercaptoethanol (BME, Sigma Aldrich), homogenized using a QIAshredder (Qiagen), and isolated with a RNeasy Kit (Qiagen) with on-column removal of DNase. DNA purity and concentration were validated with the Bioanalyzer High Sensitivity DNA Analysis chip (Agilent Technologies). Samples were sent to Novogene for RNA library preparation and transcriptome sequencing on the HiSeq 2500 platform (Illumina) using a 150-nucleotide paired-end read length [54].

FlexiVent pulmonary function testing

The P14 mice used for this analysis were not used for other analyses. Forced oscillometry was performed as described [55, 56]. Briefly, we performed a tracheostomy on the sedated mice with the appropriate cannula and secured it with two 3–0 sutures. We performed the neonatal pulmonary function assessment program on the flexiVent system (SCIREQ). We calibrated the system with the tracheal cannula before each mouse.

Histology and immunofluorescence

We gravity-inflated the right lung of each mouse to 20 cm H2O with 4% formalin (Thermo Fisher Scientific) and placed them into 100% ethanol within 16 h. After paraffin embedding, we cut 4–5 m slides at 3 random levels through the block stained with hematoxylin and eosin for morphometry. At each level, three slides were made from adjacent sections for immunofluorescence. These slides were deparaffinized and rehydrated with graduated concentrations of ethanol. We performed antigen retrieval with heated 10 µM Sodium Citrate Buffer with 20% Tween. Slides were fenestrated, blocked in 1% BSA and anti-CD16/32 (BD Biosciences) for 30 min, then incubated overnight with either anti-vWF and or anti-αSMA at 4 °C in incubation buffer (1% BSA, 1% normal donkey serum, 0.3% Triton X-100, and 0.01% sodium azide in PBS). The next day, they were washed in PBS, blocked again, counterstained with the appropriate secondary antibodies, washed, and finally counterstained with ProLong Diamond Anti-Fade Mounting Media with DAPI (Thermofisher Scientific). Images were captured using an A1R HD inverted confocal microscope (Nikon).

Bacterial and fungal microbiome sequencing

The first non-transitional stool sample collected during the second week of life was used to assess the human gut microbiota. In mice, liquid nitrogen snap-frozen homogenized whole left lung or a 1-cm section of terminal ileum were used because of the advantages in quantifying the adherent mucosal microbiota [49]. We used ZymoBIOMICS whole organism and DNA microbial community standards (Zymo Research) and sterile-filtered PBS to create pairs of positive and negative controls appropriate for each step of the sample collection and isolation procedure. Microbial DNA was extracted using the ZymoBIOMICS DNA Miniprep Kit (Zymo Research) with the addition of a lyticase (200 Units per 1 mL sample volume) predigestion with 10-min bead bashing on a Fisher Vortex Genie 2 (Thermo Fisher Scientific) at RT, then 37 °C with gentle agitation for 20 min, as previously validated to facilitate extraction of fungal DNA [50]. Extracted DNA was sequenced using the MiSeq platform (Illumina) at Novogene (human samples) or Argonne National Laboratory (mouse samples) using 16S V3-4 (bacterial) rRNA with primers 515F-806R and the internal transcribed spacer (ITS) 2 (fungal) rRNA gene primers. Primer design was chosen based on prior published optimization for mycobiome amplification [57, 58]. Quantitative RT-PCR of the ITS1, ITS2, and 16S genes was performed as described [59], using TaqMan probes (Thermo Fisher Scientific) validated with ZymoBIOMICS controls (Zymo Research) and additional human and mouse samples. QC was performed at the UAB Genomics Core.

Quantification and statistical analysis

Bioinformatics

Microbiomics

We used QIIME 2 [51] to import and trim ITS2 and 16S rRNA gene amplicon sequences. DADA2 [52] was used to denoise and cluster amplicon sequence variants (ASVs). ASV taxonomy was assigned with the RDP naive Bayesian Classifier using the SILVA 138.1 and UNITE 9.0 databases for bacteria and fungi, respectively. ASVs unassigned at the phylum level or assigned to the phylum Chytridiomycota were removed, as they had bootstrap confidence scores well below the minimum threshold of 50%, and confirmatory analysis with the eukaryotic UNITE 9.0 database revealed these to be mostly host mitochondrial or plant contaminants. Contamination source modeling was performed using the R package PERFect [53], permutation filtering with p-value ordering at an alpha value of 0.05. Decontaminated ASV tables were then exported into R where they were normalized, transformed, and analyzed with the following packages: vegan, phyloseq, DESeq2, SpiecEasi (Sparse Inverse Covariance Estimation for Ecological Association Inference), NetCoMI (Network Construction and comparison for Microbiome data), tidyverse, igraph, mia, and caret. Alpha diversity was quantified with the Shannon and Simpson Indices. We displayed beta diversity with both principal components analysis of Euclidean distances and principal coordinates analysis of Bray–Curtis dissimilarity matrices. Significance testing was performed using PERMANOVA, and confounding was evaluated by differential dispersion using PERMDISP. Feature selection was performed using DESeq2 and Random Forest using ranger and caret. Unsupervised clustering was performed using Dirichlet-Multinomial Mixture modeling using the packages mia and DirichletMultinomial. Feature selection was performed using the fitted function. Clinical data was mapped to DMM clusters using canonical correspondence analysis (CCA) using the vegan package.

RNAseq

Raw sequencing data was aligned, processed, and normalized using the CLC Genomic Workbench (Qiagen) using the mouse reference genome GRC38-mm10, and downstream analysis was performed in R. Differential expression analysis was performed using DESeq2. Differentially expressed genes were input into Ingenuity IPA (Qiagen) to identify significantly regulated networks and pathways. Additional pathway analysis was performed in R using Disease Ontology, Reactome, and Kyoto Encylopedia of Genes and Genomes (KEGG).

Statistical analysis

The human data were analyzed using SPSS Statistics (IBM). The murine data were analyzed using R and Prism 10 (GraphPad, Inc.). Animal experiments are representative of 2–3 pooled mouse litters per human donor with at least 3 human donors per experimental group and analyzed by two-way ANOVA with a Šidák correction.

Data availability

Raw RNA-seq, 16S, and ITS2 sequencing data are deposited at the NCBI Sequence Read Archive (Access Number: PRJNA982880).

Processed data files and original code are available at github.com/WillisLungLab/MiBPD.

References

  1. Lapcharoensap W, Lee HC, Nyberg A, Dukhovny D. Health care and societal costs of bronchopulmonary dysplasia. NeoReviews. 2018;19:e211–23.

    Article  PubMed  PubMed Central  Google Scholar 

  2. Lui K, Lee SK, Kusuda S, Adams M, Vento M, Reichman B, et al. Trends in outcomes for neonates born very preterm and very low birth weight in 11 high-income countries. J Pediatr. 2019;215:32-40.e14.

    Article  PubMed  Google Scholar 

  3. Kinsella JP, Greenough A, Abman SH. Bronchopulmonary dysplasia. Lancet. 2006;367:1421–31.

    Article  PubMed  Google Scholar 

  4. Thébaud B, Goss KN, Laughon M, Whitsett JA, Abman SH, Steinhorn RH, et al. Bronchopulmonary dysplasia. Nat Rev Dis Primers. 2019;5:78.

    Article  PubMed  PubMed Central  Google Scholar 

  5. Dang AT, Marsland BJ. Microbes, metabolites, and the gut–lung axis. Mucosal Immunol. 2019;12:843–50.

    Article  CAS  PubMed  Google Scholar 

  6. Willis KA, Stewart JD, Ambalavanan N. Recent advances in understanding the ecology of the lung microbiota and deciphering the gut-lung axis. Am J Physiol-lung C. 2020;319:L710–6.

    Article  CAS  Google Scholar 

  7. Stevens J, Steinmeyer S, Bonfield M, Peterson L, Wang T, Gray J, et al. The balance between protective and pathogenic immune responses to pneumonia in the neonatal lung is enforced by gut microbiota. Sci Transl Med. 2022;14:eabl3981. https://doiorg.publicaciones.saludcastillayleon.es/10.1126/scitranslmed.abl3981.

  8. Deshmukh HS, Liu Y, Menkiti OR, Mei J, Dai N, O’Leary CE, et al. The microbiota regulates neutrophil homeostasis and host resistance to Escherichia coli K1 sepsis in neonatal mice. Nat Med. 2014;20:524–30.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Bernardes EVT, Pettersen VK, Gutierrez MW, Laforest-Lapointe I, Jendzjowsky NG, Cavin J-B, et al. Intestinal fungi are causally implicated in microbiome assembly and immune development in mice. Nat Commun. 2020;11:2577–616.

    Article  Google Scholar 

  10. Arrieta M-C, Stiemsma LT, Dimitriu PA, Thorson L, Russell S, Yurist-Doutsch S, et al. Early infancy microbial and metabolic alterations affect risk of childhood asthma. Sci Transl Med. 2015;7:307ra152.

    Article  PubMed  Google Scholar 

  11. Tipton L, Müller CL, Kurtz ZD, Huang L, Kleerup E, Morris A, et al. Fungi stabilize connectivity in the lung and skin microbial ecosystems. Microbiome. 2018;6:12.

    Article  PubMed  PubMed Central  Google Scholar 

  12. Abdelgawad A, Nicola T, Martin I, Halloran BA, Tanaka K, Adegboye CY, et al. Antimicrobial peptides modulate lung injury by altering the intestinal microbiota. Microbiome. 2023;11:226.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  13. Cantu A, Gutierrez MC, Dong X, Leek C, Sajti E, Lingappan K. Remarkable sex-specific differences at single-cell resolution in neonatal hyperoxic lung injury. Am J Physiol-lung C. 2022;324(1):L5–31.

  14. Trembath A, Laughon MM. Predictors of bronchopulmonary dysplasia. Clin Perinatol. 2012;39:585–601.

    Article  PubMed  PubMed Central  Google Scholar 

  15. Holmes I, Harris K, Quince C. Dirichlet multinomial mixtures: generative models for microbial metagenomics. PLoS ONE. 2012;7:e30126.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  16. Koren O, Konnikova L, Brodin P, Mysorekar IU, Collado MC. The maternal gut microbiome in pregnancy: implications for the developing immune system. Nat Rev Gastroenterol Hepatol. 2024;21:35–45.

    Article  PubMed  Google Scholar 

  17. Nardiello C, Mižíková I, Silva DM, Ruiz-Camp J, Mayer K, Vadász I, et al. Standardisation of oxygen exposure in the development of mouse models for bronchopulmonary dysplasia. Dis Model Mech. 2017;10:185–96.

    CAS  PubMed  PubMed Central  Google Scholar 

  18. Dolma K, Freeman AE, Rezonzew G, Payne GA, Xu X, Jilling T, et al. Effects of hyperoxia on alveolar and pulmonary vascular development in germ free mice. Am J Physiol- Lung Cell Mol Physiol. 2019;295:L86.

    Google Scholar 

  19. Jiang TT, Shao T-Y, Ang WXG, Kinder JM, Turner LH, Pham G, et al. Commensal fungi recapitulate the protective benefits of intestinal bacteria. Cell Host Microbe. 2017;22:809-816.e4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Iliev ID, Funari VA, Taylor KD, Nguyen Q, Reyes CN, Strom SP, et al. Interactions between commensal fungi and the C-type lectin receptor Dectin-1 influence colitis. Science (New York, NY). 2012;336:1221789–31317.

    Article  Google Scholar 

  21. Laforest-Lapointe I, Arrieta M-C. Microbial eukaryotes: a missing link in gut microbiome studies. Msystems. 2018;3:e00201-e217.

    Article  PubMed  PubMed Central  Google Scholar 

  22. Ashley SL, Sjoding MW, Popova AP, Cui TX, Hoostal MJ, Schmidt TM, et al. Lung and gut microbiota are altered by hyperoxia and contribute to oxygen-induced lung injury in mice. Sci Transl Med. 2020;12:eaau9959.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Willis KA, Siefker DT, Aziz MM, White CT, Mussarat N, Gomes CK, et al. Perinatal maternal antibiotic exposure augments lung injury in offspring in experimental bronchopulmonary dysplasia. Am J Physiol- Lung Cell Mol Physiol. 2019;3:21.

    Google Scholar 

  24. Bongers KS, Chanderraj R, Woods RJ, McDonald RA, Adame MD, Falkowski NR, et al. The gut microbiome modulates body temperature both in sepsis and health. Am J Resp Crit Care. 2023;207(8):1030–41.

  25. Moutsoglou DM, Tatah J, Prisco SZ, Prins KW, Staley C, Lopez S, et al. Pulmonary arterial hypertension patients have a proinflammatory gut microbiome and altered circulating microbial metabolites. Am J Resp Crit Care. 2023;207(6):740–56.

  26. Tashiro H, Shore SA. The gut microbiome and ozone-induced airway hyperresponsiveness mechanisms and therapeutic prospects. Am J Resp Cell Mol. 2021;64:283–91.

    Article  CAS  Google Scholar 

  27. Wang Z, Locantore N, Haldar K, Ramsheh MY, Beech AS, Ma W, et al. Inflammatory endotype–associated airway microbiome in chronic obstructive pulmonary disease clinical stability and exacerbations: a multicohort longitudinal analysis. Am J Resp Crit Care. 2021;203:1488–502.

    Article  CAS  Google Scholar 

  28. Dickson RP, Schultz MJ, van der Poll T, Schouten LR, Falkowski NR, Luth JE, et al. Lung microbiota predict clinical outcomes in critically ill patients. Am J Resp Crit Care. 2020;201:555–63.

    Article  CAS  Google Scholar 

  29. O’Dwyer DN, Ashley SL, Gurczynski SJ, Xia M, Wilke C, Falkowski NR, et al. Lung microbiota contribute to pulmonary inflammation and disease progression in pulmonary fibrosis. Am J Resp Crit Care. 2019;199:1127–38.

    Article  Google Scholar 

  30. Dickson RP, Singer BH, Newstead MW, Falkowski NR, Erb-Downward JR, Standiford TJ, et al. Enrichment of the lung microbiome with gut bacteria in sepsis and the acute respiratory distress syndrome. Nat Microbiol. 2016;1:16113.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Combs MP, Wheeler DS, Luth JE, Falkowski NR, Walker NM, Erb-Downward JR, et al. Lung microbiota predict chronic rejection in healthy lung transplant recipients: a prospective cohort study. Lancet Respir Medicine. 2021;9:601–12.

    Article  CAS  Google Scholar 

  32. Alcazar CG-M, Paes VM, Shao Y, Oesser C, Miltz A, Lawley TD, et al. The association between early-life gut microbiota and childhood respiratory diseases: a systematic review. Lancet Microbe. 2022;3:e867–80.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Xia S, Ellis LV, Winkley K, Menden H, Mabry SM, Venkatraman A, et al. Neonatal hyperoxia induces activated pulmonary cellular states and sex-dependent transcriptomic changes in a model of experimental bronchopulmonary dysplasia. Am J Physiol-lung C. 2023;324(2):L123–40.

  34. Lao JC, Bui CB, Pang MA, Cho SX, Rudloff I, Elgass K, et al. Type 2 immune polarization is associated with cardiopulmonary disease in preterm infants. Sci Transl Med. 2022;14:eaaz8454.

    Article  CAS  PubMed  Google Scholar 

  35. Hurskainen M, Mižíková I, Cook DP, Andersson N, Cyr-Depauw C, Lesage F, et al. Single cell transcriptomic analysis of murine lung development on hyperoxia-induced damage. Nat Commun. 2021;12:1565.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Zaneveld JR, McMinds R, Thurber RV. Stress and stability: applying the Anna Karenina principle to animal microbiomes. Nat Microbiol. 2017;2:17121.

    Article  CAS  PubMed  Google Scholar 

  37. Lal CV, Ambalavanan N. Biomarkers, early diagnosis, and clinical predictors of bronchopulmonary dysplasia. Clin Perinatol. 2015;42:739–54.

    Article  PubMed  PubMed Central  Google Scholar 

  38. Gallacher D, Mitchell E, Alber D, Wach R, Klein N, Marchesi JR, et al. Dissimilarity of the gut-lung axis and dysbiosis of the lower airways in ventilated preterm infants. Eur Respir J. 2020;1901909. https://doiorg.publicaciones.saludcastillayleon.es/10.1183/13993003.01909-2019.

  39. Underhill DM, Iliev ID. The mycobiota: interactions between commensal fungi and the host immune system. Nat Rev Immunol. 2014;14:405–16.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  40. Willis KA, Purvis JH, Myers ED, Aziz MM, Karabayir I, Gomes CK, et al. Fungi form interkingdom microbial communities in the primordial human gut that develop with gestational age. FASEB J. 2019;33:12825–37.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Higgins RD, Jobe AH, Koso-Thomas M, Bancalari E, Viscardi RM, Hartert TV, et al. Bronchopulmonary dysplasia: executive summary of a workshop. J Pediatrics. 2018;197:300–8.

    Article  Google Scholar 

  42. Laughon MM, Langer JC, Bose CL, Smith PB, Ambalavanan N, Kennedy KA, et al. Prediction of bronchopulmonary dysplasia by postnatal age in extremely premature infants. Am J Resp Crit Care. 2011;183:1715–22.

    Article  Google Scholar 

  43. von Elm E, Altman DG, Egger M, Pocock SJ, Gøtzsche PC, Vandenbroucke JP, et al. The Strengthening the Reporting of Observational Studies in Epidemiology (STROBE) statement: guidelines for reporting observational studies. Plos Med. 2007;4:e296.

    Article  Google Scholar 

  44. Kilkenny C, Browne WJ, Cuthill IC, Emerson M, Altman DG. Improving bioscience research reporting: the ARRIVE guidelines for reporting animal research. 2010.

  45. Staley C, Kaiser T, Beura LK, Hamilton MJ, Weingarden AR, Bobr A, et al. Stable engraftment of human microbiota into mice with a single oral gavage following antibiotic conditioning. Microbiome. 2017;5:87.

    Article  PubMed  PubMed Central  Google Scholar 

  46. Chen X, Li P, Liu M, Zheng H, He Y, Chen M-X, et al. Gut dysbiosis induces the development of pre-eclampsia through bacterial translocation. Gut. 2020;513–22. https://doiorg.publicaciones.saludcastillayleon.es/10.1136/gutjnl-2019-319101

  47. Wrzosek L, Ciocan D, Borentain P, Spatz M, Puchois V, Hugot C, et al. Transplantation of human microbiota into conventional mice durably reshapes the gut microbiota. Sci Rep. 2018;8:6854.

    Article  PubMed  PubMed Central  Google Scholar 

  48. Ezz-Eldin ZM, Hamid TAA, Youssef MRL, Nabil HE-D. Clinical Risk Index for Babies (CRIB II) scoring system in prediction of mortality in premature babies. J Clin Diagnostic Res. 2015;9:SC08-11.

    Google Scholar 

  49. Baker JM, Hinkle KJ, McDonald RA, Brown CA, Falkowski NR, Huffnagle GB, et al. Whole lung tissue is the preferred sampling method for amplicon-based characterization of murine lung microbiota. Microbiome. 2021;9:99.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Pierre JF, Phillips GJ, Chandra LC, Rendina DN, Thomas-Gosain NF, Lubach GR, et al. Lyticase facilitates mycobiome resolution without disrupting microbiome fidelity in primates. J Surg Res. 2021;267:336–41.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  51. Bolyen E, Rideout JR, Dillon MR, Bokulich NA, Abnet CC, Al-Ghalith GA, et al. Reproducible, interactive, scalable and extensible microbiome data science using QIIME 2. Nat Biotechnol. 2019;37:852–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Callahan BJ, McMurdie PJ, Rosen MJ, Han AW, Johnson AJA, Holmes SP. DADA2: High-resolution sample inference from Illumina amplicon data. Nat Methods. 2016;13:581–3.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Smirnova E, Huzurbazar S, Jafari F. PERFect: PERmutation Filtering test for microbiome data. Biostatistics. 2018;20:615–31.

    Article  PubMed Central  Google Scholar 

  54. Wang W, Yan T, Guo W, Niu J, Zhao Z, Sun K, et al. Constitutive GLI1 expression in chondrosarcoma is regulated by major vault protein via mTOR/S6K1 signaling cascade. Cell Death Differ. 2021;28:2221–37.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  55. Nicola T, Ambalavanan N, Zhang W, James ML, Rehan V, Halloran B, et al. Hypoxia-induced inhibition of lung development is attenuated by the peroxisome proliferator-activated receptor-γ agonist rosiglitazone. Am J Physiol- Lung Cell Mol Physiol. 2011;301:L125–34.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  56. Nicola T, Hagood JS, James ML, MacEwen MW, Williams TA, Hewitt MM, et al. Loss of Thy-1 inhibits alveolar development in the newborn mouse lung. Am J Physiol-lung C. 2009;296:L738–50.

    Article  CAS  Google Scholar 

  57. Tiew PY, Dicker AJ, Keir HR, Poh ME, Pang SL, Aogáin MM, et al. A high-risk airway mycobiome is associated with frequent exacerbation and mortality in COPD. Eur Respir J. 2021;57:2002050.

    Article  CAS  PubMed  Google Scholar 

  58. Ali NABM, Aogáin MM, Morales RF, Tiew PY, Chotirmall SH. Optimisation and benchmarking of targeted amplicon sequencing for mycobiome analysis of respiratory specimens. Int J Mol Sci. 2019;20:4991.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  59. Jilling T, Ren C, Yee A, Aggarwal S, Halloran B, Ambalavanan N, et al. Exposure of neonatal mice to bromine impairs their alveolar development and lung function. Am J Physiol- Lung Cell Mol Physiol. 2018;314:L137–43.

    Article  PubMed  Google Scholar 

Download references

Acknowledgements

We thank the infants, their parents, and the study coordinators, without whom this study would have been impossible. Confocal microscopy was performed at the UAB High-Resolution Imaging Facility with the assistance of Robert Grabski. We support the inclusive, diverse, and equitable conduct of research. Multiple authors of this paper self-identify as underrepresented ethnic or gender minorities or as members of the LGBTQA+ community.

Funding

Research reported in this article was supported by the National Institutes of Health under award numbers K08 HL151907 (KW), L40 HL159581 (KW), KL2 TR3097 (VJ), K08 HL141652 (CVL), R21 HD100917 (SG), R01 HL156275 (NA), R01 HL129907 (NA), R01 AI134796 (BP), R01 CA253329 (JP), and R21 AI163503 (JP), the Dixon Fellowship (MS), the Microbiome Center of the University of Alabama at Birmingham (KW), and the Division of Neonatology at the University of Tennessee Health Science Center (AT). The content is solely the authors' responsibility and does not necessarily represent the official views of the respective funders.

Author information

Author notes

  1. Ajay J. Talati and Kent A. Willis are both co-senior authors.

Authors and Affiliations

  1. Division of Neonatology, Department of Pediatrics, Heersink School of Medicine, The University of Alabama at Birmingham, Birmingham, AL, USA

    Isaac Martin, Mary Silverberg, Ahmed Abdelgawad, Kosuke Tanaka, Brian A. Halloran, Teodora Nicola, Samuel J. Gentle, Namasivayam Ambalavanan, Viral G. Jain, Charitharth V. Lal, Tamás Jilling & Kent A. Willis

  2. College of Medicine, The University of Tennessee Health Science Center, Memphis, TN, USA

    Erin D. Myers

  3. Division of Neonatology, Department of Pediatrics, College of Medicine, The University of Tennessee Health Science Center, Memphis, TN, USA

    Jay P. Desai, Catrina T. White & Ajay J. Talati

  4. Division of Cardiology, Department of Medicine, Wake Forest University School of Medicine, Wake Forest University, Winston-Salem, NC, USA

    Ibrahim Karabayir & Oguz Akbilgic

  5. Epidemiological Cardiology Research Center, Wake Forest University, Winston-Salem, NC, USA

    Oguz Akbilgic

  6. Departments of Biology and Mathematics & Statistics, James Madison University, Harrisonburg, VA, USA

    Laura Tipton

  7. Department of Clinical Pharmacy and Translational Science, College of Pharmacy, The University of Tennessee Health Science Center, Memphis, TN, USA

    Brian M. Peters

  8. Department of Molecular Microbiology and Immunology, The University of Texas at San Antonio, San Antonio, TX, USA

    Luan D. Vu

  9. Department of Biological Sciences, Louisiana State University, Baton Rouge, LA, USA

    Stephania A. Cormier

  10. Pennington Biomedical Research Center, Baton Rouge, LA, USA

    Stephania A. Cormier

  11. Department of Nutritional Sciences, College of Agricultural and Life Sciences, The University of WI–Madison, Madison, WI, USA

    Joseph F. Pierre

  12. Department of Obstetrics and Gynecology, College of Medicine, The University of Tennessee Health Science Center, Memphis, TN, USA

    Ajay J. Talati

Authors
  1. Isaac Martin
    View author publications

    You can also search for this author inPubMed Google Scholar

  2. Mary Silverberg
    View author publications

    You can also search for this author inPubMed Google Scholar

  3. Ahmed Abdelgawad
    View author publications

    You can also search for this author inPubMed Google Scholar

  4. Kosuke Tanaka
    View author publications

    You can also search for this author inPubMed Google Scholar

  5. Brian A. Halloran
    View author publications

    You can also search for this author inPubMed Google Scholar

  6. Teodora Nicola
    View author publications

    You can also search for this author inPubMed Google Scholar

  7. Erin D. Myers
    View author publications

    You can also search for this author inPubMed Google Scholar

  8. Jay P. Desai
    View author publications

    You can also search for this author inPubMed Google Scholar

  9. Catrina T. White
    View author publications

    You can also search for this author inPubMed Google Scholar

  10. Ibrahim Karabayir
    View author publications

    You can also search for this author inPubMed Google Scholar

  11. Oguz Akbilgic
    View author publications

    You can also search for this author inPubMed Google Scholar

  12. Laura Tipton
    View author publications

    You can also search for this author inPubMed Google Scholar

  13. Samuel J. Gentle
    View author publications

    You can also search for this author inPubMed Google Scholar

  14. Namasivayam Ambalavanan
    View author publications

    You can also search for this author inPubMed Google Scholar

  15. Brian M. Peters
    View author publications

    You can also search for this author inPubMed Google Scholar

  16. Luan D. Vu
    View author publications

    You can also search for this author inPubMed Google Scholar

  17. Viral G. Jain
    View author publications

    You can also search for this author inPubMed Google Scholar

  18. Charitharth V. Lal
    View author publications

    You can also search for this author inPubMed Google Scholar

  19. Stephania A. Cormier
    View author publications

    You can also search for this author inPubMed Google Scholar

  20. Joseph F. Pierre
    View author publications

    You can also search for this author inPubMed Google Scholar

  21. Tamás Jilling
    View author publications

    You can also search for this author inPubMed Google Scholar

  22. Ajay J. Talati
    View author publications

    You can also search for this author inPubMed Google Scholar

  23. Kent A. Willis
    View author publications

    You can also search for this author inPubMed Google Scholar

Contributions

Conceptualization: NA, BP, CL, SC, JP, AT, KW. Clinical Sample and Data Collection: IM, MS, EM, JD, SG, VJ, KW. Experimentation: MS, AA, BH, KT, TN, CW, SJ, KW. Analysis: IM, IK, OA, LT, LV, VJ, TJ, KW. Funding Acquisition: AT, KW Supervision: AT, KW Manuscript: IM, KW. All authors have approved the final manuscript as submitted.

Corresponding author

Correspondence to Kent A. Willis.

Ethics declarations

Ethics approval and consent to participate

The human study protocol was approved by the Institutional Review Board of The University of Tennessee under protocol 17–05311-XP and conducted in accordance with the Declaration of Helsinki. The human cohort for FMT samples was approved under protocol IRB-300003994 by the IRB at the University of Alabama at Birmingham (UAB, Birmingham, AL, USA). Mouse experiments were performed at UAB under IACUC protocol 22042.

Consent for publication

Assent for publication was obtained from the mother during study enrollment.

Competing interests

CL is the founder of ResBiotic Nutrition, Inc. and Alveolus Bio, Inc., from which he reports salary and stock options, and TN is now an employee, thereof. KW and NA are scientific advisors to these companies and report consulting fees and stock options from the same. The remaining authors report no disclosures.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

40168_2025_2032_MOESM1_ESM.pdf

Supplementary Material 1: Figure S1. Source contamination modeling. a) Rarefication plot of unfiltered ITS2 sequencing. b) Rarefication plot of ITS2 sequencing after decontamination modeling with PERFect. c) Rarefication plot of unfiltered 16S sequencing. d) Rarefication plot of 16S sequencing after decontamination modeling with PERFect. e) Fungal taxa filtering loss plots. f) Bacterial taxa filtering loss plots. Please see tables of filtered taxa at github.com/WillisLungLab/MiBPD/SourceContamination. Figure S2. Source contamination modeling. a) Fungal contaminate taxa detected by PERFect. b) Bacterial contaminate taxa detected by PERFect. c) Average abundance of filtered fungal taxa. d) Average abundance of filtered bacterial taxa. Please see tables of filtered taxa at github.com/WillisLungLab/MiBPD/SourceContamination. Figure S3. The BPD-associated mycobiome is unique. (a, b) Relative abundance of key fungi by DESeq2. (c,d) Relative abundance of key bacteria by DESeq2. (d) Number of edges in each SPIEC-EASI network. (e) Number of nodes in each SPIEC-EASI network. (e) Complete SPIEC-EASI BPD network with all genera labeled. (F) Complete NoBPD network with all genera labeled. Figure S4. The overall composition of the multikingdom microbiome. (a) Alpha diversity. (c) Beta diversity (PCoA of Bray-Curtis dissimilarity, p = 0.0743, PERMANOVA, p = 0.3749, PERMDISP). Ellipses indicate the 95% confidence interval. Figure S5. BPD status is predictable using machine-learning models of the fungal microbiota several weeks before the adjudication of the clinical diagnosis at 36 weeks corrected gestational age. (a-b) Exploratory random forest machine learning model using genera level mycobiome and clinical data (Should be considered exploratory due to high risk of overfitting). (c-d) Exploratory random forest machine learning model using genera-level bacterial microbiota and clinical data. Figure S6. Mycobiome features differences by sex assigned at birth. (a) Fungal alpha diversity. (b) Bacterial alpha diversity. Shannon diversity decreased in infants assigned female at birth (AFAB) as compared to infants assigned male at birth (AMAB, p = 0.0427, Mann-Whitney). Simpson diversity was increased in infants AFAB as compared to those AMAB (p = 0.0239, Mann-Whitney). (c) Principal coordinates analysis of BrayCurtis dissimilarity of fungal beta diversity (p = 0.2279, PERMANOVA; p = 0.9464, PERMDISP). (d) Principal coordinates analysis of Bray-Curtis dissimilarity of bacterial beta diversity (p = 0.7846, PERMANOVA; p = 0.3957, PERMDISP). Figure S7. Microbiota alterations by disease status and sex assigned at birth. (a-b) DESeq2 selected features of the mycobiota. (c) SPIEC-EASI network edges. (d) SPIEC-EASI network nodes. (e) SPIEC-EASI networks for NoBPD infants assigned male at birth (AMAB) or assigned female at birth (AFAB). (f) SPIEC-EASI networks for BPD infants. Figure S8. Dirichlet-Multinomial Mixture (DMM) modeling of the mycobiome. (a) Laplace analysis of the mycobiome indicates optimal model fit points with 5 clusters. (b) Variable importance analysis of the fungal DMM model. (c) Canonical Correspondence Analysis (CCA) for the DMM model. (d) Laplace analysis of mycobiome samples from NoBPD infants shows an optimal fit point at 2 clusters. (e) Variable importance analysis of the NoBPD DMM model. (f) CCA for the NoBPD DMM model. (g) Laplace analysis of mycobiome samples from BPD infants shows an optimal fit point at 2 clusters. (h) Variable importance analysis of the BPD DMM model. (i) CCA for the BPD DMM model. Figure S9. Dirichlet-Multinomial Mixture (DMM) modeling of the bacterial microbiome. (a) Principal Coordinates Analysis (PCoA) of the bacterial microbiome DMM model. (b) Alpha diversity of the bacterial microbiome DMM model. (c) Relative abundance of the bacterial microbiome in the DMM components. (d) PCoA of the DMM model of the bacterial microbiome of infants with NoBPD. (e) Alpha diversity of the multikingdom microbiome. (f) Relative abundance of the bacterial microbiome in the NoBPD DMM model. (g) PCoA of the DMM model of the bacterial microbiome of infants with BPD. (h) Alpha diversity of the bacterial microbiome. (i) Relative abundance of the multikingdom microbiome in the BPD DMM model. Ellipses indicated the 95% confidence interval. Figure S10. DMM modeling of the bacterial microbiome. (a) Laplace analysis of the mycobiome indicates an optimal model fit point with four components. (b) Variable importance analysis of the multikingdom microbiome DMM model. (c) Canonical Correspondence Analysis (CCA) for the multikingdom microbiome DMM model. (d) Laplace analysis of multikingdom microbiome samples from NoBPD infants shows an optimal fit point at three components. (e) Variable importance analysis of the NoBPD DMM model. (f) CCA for the NoBPD DMM model. (g) Laplace analysis of multikingdom microbiome samples from BPD infants shows an optimal fit point at two components. (h) Variable importance analysis of the BPD DMM model. (i) CCA for the BPD DMM model. Figure S11. Comparison between the fungal microbiomes of FMT-colonized and SPF control mice. Figure S12. Comparative RNA-seq pathway analysis for fecal microbiota transfer experiment. Comparative gene enrichment analysis was performed using the top UP and DOWN regulated genes in NoBPD hyperoxia vs NoBPD normoxia compared to BPD hyperoxia vs NoBPD normoxia top regulated genes. Top pathways (based on adjusted p-value) are shown. (a) Top (5) Disease Ontology terms between BPD-FMT HO and NoBPD-FMT NO mice. (b) Network of top (2) DO terms and genes. (c) Top (5) KEGG pathways between BPD-FMT HO and NoBPD-FMT NO mice. (d) Network of top (2) KEGG pathways and genes. (c) Top (5) Reactome pathways between BPD-FMT HO and NoBPD-FMT NO mice. (d) Network of top (2) Reactome pathways and genes. Figure S13. Ingenuity pathway analysis. (A) Pathways in BPD-FMT mice. (B) Pathways in NoBPD-FMT mice. (C) Pathways in specific pathogen-free (SPF) mice. Figure S14. Read counts in fluconazole-exposed neonatal mice. (a) ITS2 read counts in neonatal mouse gut samples. (b) ITS2 read counts in neonatal mouse lung samples. Figure S15. Relative abundance of fungal genera between fluconazole and control mice. Potential anti-fungal amplification of selected fungal colonization occurred in only one sample (G5). Of note, absolute fungal colonization in the fluconazole-exposed group is so low, no significant differences could be detected by DESeq2. Figure S16. Flow cytometry gating strategy. ILC2s were defined as live, lineage negative, ST2+CD127+Gata3+ cells. FMO, fluorescence minus one. Figure S17. RNA-seq pathway analysis of fluconazole-exposed mice. (a) Significant Disease Ontology terms between Fluco HO and SPF NO mice. (b) Network of significant DO terms and genes. (c) Significant KEGG pathways between Fluco HO and SPF NO mice. (d) Network of significant KEGG pathways and genes. (c) Significant Reactome pathways between Fluco HO and SPF NO mice. (d) Network of significant Reactome pathways and genes. Figure S18. (e) Gut Candida and total fungal relative abundance. (f) Lung Candida and total fungal relative abundance. Figure S19. RNA-seq pathway analysis of C.tropicalis-colonized mice. (a) Significant Disease Ontology terms between C.trop HO and SPF NO mice. (b) Network of significant DO terms and genes. (c) Significant KEGG pathways between C.trop HO and SPF NO mice. (d) Network of significant KEGG pathways and genes. (c) Significant Reactome pathways between C.trop HO and SPF NO mice. (d) Network of significant Reactome pathways and genes. Table S1. Cohort characteristics and demographics. Table S2. Characteristics and demographics of infants without a stool for analysis. Table S3. PERMANOVA analysis of potential fungal covariates. Table S4. PERMANOVA analysis of potential bacterial covariates. Table S5. Type 2 Innate Lymphoid Cell Flow Cytometry Panel.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Martin, I., Silverberg, M., Abdelgawad, A. et al. The fungal microbiota modulate neonatal oxygen-induced lung injury. Microbiome 13, 24 (2025). https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s40168-025-02032-x

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s40168-025-02032-x

Keywords

Microbiome

ISSN: 2049-2618

Contact us