Main

Necrotizing enterocolitis (NEC) is a common and devastating disease that predominantly affects premature infants. Experimental and clinical research over four decades has revealed much about NEC risk factors (prematurity, enteral feeding, and exposure to antibiotics) and contributing mechanisms (dysbiosis, altered apoptosis, translocation of bacteria through an inadequate barrier, exuberant pro-inflammatory cascade, and in severe cases coagulopathy, gut necrosis, and shock) (1,2,3). The most promising preventative approaches to date include provision of human milk and probiotics (4,5). The observations that many Enterobacteriaceae outcompete commensal organisms in the inflamed intestine by utilizing an alternative respiratory pathway (6), together with the recent descriptions of a bloom of Enterobacteriaceae associated with NEC (7) shed new light on a possible central role of the intestinal microbiota in this disease.

The neonatal rat model of NEC is an invaluable experimental tool for examining the pathogenesis of NEC and potential mechanisms of protection (8,9,10). The strength of this model is its inclusion of stressors and enteral feeding, both of which are factors associated with human NEC. In the rat model, the stressors include separation from the dam, tube feeding, hypoxia, hypothermia, and enteral nourishment with bovine-based rat milk substitute (11).

Previous studies with this model have demonstrated a protective effect of probiotic Bifidobacterium bifidum with decreased NEC, decreased apoptosis, and decreased inflammation (8,9,10). Mouse and piglet studies have demonstrated alterations of the intestinal microbiota in NEC (12,13); however, changes in the microbiota in the rat NEC model and the impact of probiotic bifidobacteria on this microbial community have not yet been investigated.

In this study, we chose to analyze a different probiotic strain, Bifidobacterium longum subsp. infantis (B. infantis), for the following reasons. First, this bacterial strain has evolved the genetic capacity to thrive in the breast-fed human infant (14). Encoded in its genome are several unique glycosidases capable of deconstructing human milk oligosaccharides that are not found in other bifidobacteria (15). Second, B. infantis grown in culture media containing human milk oligosaccharides shows increased epithelial cell adhesion compared with B. bifidum grown under identical conditions (16). Third, this strain was shown in a phase 1 trial to be a superior colonizer of the premature infant intestinal tract (17). Finally, in a population of breast-fed infants in Bangladesh, relative percentages of B. infantis were associated with improved growth and increased responsiveness to several routine vaccines (18). We hypothesized that B. infantis in the rat NEC model would decrease the severity of NEC, decrease markers of inflammation, alter intestinal antimicrobial peptides and mediators of mucus production, and alter the cecal microbiota by decreasing Enterobacteriaceae and increasing bifidobacteria.

Results

B. infantis Decreased the Incidence and Severity of NEC

Figure 1 presents the incidence of NEC, the ileal histology scores, and the villus length and width for all the animals exposed to asphyxia and cold stress (n=50). Consistent with previous investigations (8,9,10), the formula-fed (FF) group showed an increased incidence of NEC ( Figure 1a ), an increased histologic score for NEC ( Figure 1b ), a decreased mean villous length ( Figure 1c ), and a decreased mean villous width compared with dam-fed (DF) controls ( Figure 1d ). The administration of B. infantis was protective in the FF group as manifested by a decreased incidence of NEC, a decreased histologic score, and an increased mean villous length ( Figure 1a c ).

Figure 1
figure 1

B. infantis decreased the incidence and severity of necrotizing enterocolitis (NEC). (a) Incidence of NEC P < 0.01; (b) NEC histology scores, horizontal line represents the median value for each group; (c) mean length, and (d) width of villi. ANOVA of all three groups P < 0.01 in b, c, and d. **P < 0.01, * P < 0.05 for between-group comparisons. Error bars represent SEM. White bars and DF = dam fed, gray bars and FF = formula fed, hatched bars and FF+Binf = formula with added B. infantis. Note that the incidence of NEC in the DF group was zero.

PowerPoint slide

B. infantis Decreased Inflammation

Figure 2a summarizes the relative mRNA expression encoding pro-inflammatory (Il6, Cxcl1, Tnfa, and Il23) and anti-inflammatory (Il10) cytokines important in the pathogenesis of NEC. The significantly increased mRNA of pro-inflammatory cytokines in the FF group was attenuated by administration of B. infantis. Il10 mRNA was suppressed by formula feeding, but was not significantly rescued by B. infantis.

Figure 2
figure 2

B. infantis attenuated inflammatory responses. (a) B. infantis attenuated expression of several pro-inflammatory cytokines. ANOVA comparison of all three groups P < 0.01 for Il6, Cxcl1, and Il10 and P < 0.05 for Tnfa and Il23. (b) B. infantis attenuated nitric oxide synthase (a marker of inflammation) and Toll-like receptor 2. ANOVA comparison of all three groups P < 0.01 for Nos2, Tlr2, Tlr4 (mRNA), and TLR4 (protein). (c) B. infantis attenuated altered expression of antimicrobial peptides. ANOVA comparison of all three groups P < 0.01 for Reg3b and P < 0.05 for Reg3g. (d) Experimental NEC altered expression of trefoil factor 3 and MUC2; B. infantis altered expression of trefoil factor 3. ANOVA comparison of all three groups P < 0.01 for TFF3 and MUC2. **P < 0.01, *P < 0.05 for between-group comparisons. Error bars represent SEM. White bars = dam fed (DF), gray bars = formula fed (FF), and hatched bars = formula with added B. infantis (FF+Binf).

PowerPoint slide

Inducible nitric oxide synthase (iNOS or Nos2) is a marker of inflammation and Toll-like receptors (TLR) 2 and 4 are important sensors of microbial patterns that trigger inflammatory responses. Figure 2b summarizes iNOS, TLR4, and TLR2 expression at the protein and mRNA levels in the three treatment groups. Both Nos2 and Tlr2 were significantly increased in the FF group and attenuated by B. infantis at the mRNA level, but these changes were not significant at the protein level in a small subset of specimens. Tlr4 was unexpectedly decreased in the FF and FF+Binf group and not consistent with previous observations (11,19). We currently have no explanation for this observation.

B. infantis Attenuated Increased Expression of Antimicrobial Peptides

Reg3B and Reg3G, antimicrobial proteins of the Reg family of C-type lectins, are produced by both Paneth cells and enterocytes, and are secreted into the intestinal lumen in even very young rat pups (20,21). Reg3b and Reg3g were increased in the FF group and attenuated in the FF+Binf group ( Figure 2c ). These changes were similar to previous observations in this model with administration of B. bifidum (10).

B. infantis Altered Mucus Production

TFF3 is a peptide secreted by goblet cells and serves to stabilize the mucus layer and trigger regeneration of the injured epithelial layer, while MUC2 is the predominant glycoprotein component of the mucus layer. A decrease in TFF3-positive cells in the FF group was rescued by B. infantis ( Figure 2d ). MUC2 was increased in both the FF and FF+Binf groups. These changes occurred in the opposite direction from those reported with administration of B. bifidum (8), suggesting the possibility that bifidobacterial strains have differing impact on intestinal mucus production.

Changes in the Cecal Microbiota Were Inconsistent

The microbiota of the cecal lumen was determined in animals identically challenged in two different experiments performed 4 mo apart. At the phylum level, there was significant heterogeneity in the cecal microbiota within a given animal group and between experiments in spite of identical diet and experimental conditions. In all animals the dominant phyla were Firmicutes and Proteobacteria. In experiment A ( Figure 3a ), the cecal microbiota was dominated by Firmicutes in the DF group and by Proteobacteria in the FF group. This pattern was not seen in experiment B ( Figure 3b ). Although the number is too small for statistical analysis, there were no obvious differences in the microbiota when comparing animals exposed or not exposed to asphyxia and cold stress within each group (in experiment A, Figure 3a ; all animals in experiment B were exposed to asphyxia and cold stress).

Figure 3
figure 3

Alterations in the cecal microbiota. (a) Cecal microbiota at the phylum level in experiment A. (b) Cecal microbiota at the phylum level in experiment B. Each bar represents an individual animal. Stress = hypoxia and cold stress. DF = dam fed, FF = formula fed, FF+Binf = formula with added B. infantis. In experiment B, all of the animals were exposed to hypoxia and cold stress. Green = Firmicutes, purple = Proteobacteria, the following represented < 3% of total bacteria identified and are visible in the lower aspect of the bars for only a few animals: red = Actinobacteria, yellow = Bacteroidetes, light blue = Deferribacteres. (c) Linear discriminate analysis (LDA) scores for taxa differing between treatment groups in experiment A (only the animals exposed to cold stress and hypoxia are included). LDA score correlates roughly with P values compared with the other groups (e.g., LDA score of 1 ~ P value of 0.1, LDA score of 2 ~ P value of 0.01, etc.). Blue bars represent taxa significantly increased in the FF+Binf group (from top to bottom four unidentified taxa, Turicibacteraceae, Turicibacterales, Bacillales, Staphylococcaceae, Clostridiaceae). Green bars represent taxa significantly increased in the FF group (from top to bottom γ-Proteobacteria, Enterobacteriales, Enterobacteriaceae, Proteobacteria, Exiguobacteraceae, Exiguobacterales, Enterococcaceae). Red bars represent taxa significantly increased in the DF group (from top to bottom Lactobacillaceae, Lactobacillales, Bacilli, Firmicutes, Veillonellaceae, Actinobacteria, Actinobacteria, Actinomycetales, Micrococcaceae, Streptococcaceae). (d) LDA scores from experiment B. The color scheme is the same as in panel c with significantly increased taxa (top to bottom) in the FF+Binf group (Lactobacillaceae, Bacilli, Lactobacillales, Firmicutes, Actinomycetales, Actinobacteria, Actinobacteria, Micrococcaceae), the FF group (Enterococcaceae), and the DF group (γ-Proteobacteria, Enterobacteriales, Enterobacteriaceae, Proteobacteria). (e) Principal coordinate analysis from experiment A. Each symbol represents one animal (red circles = DF, blue triangles = FF, green squares = FF+Binf). (f) Principal coordinate analysis from experiment B (the color/symbol scheme is the same as panel e).

PowerPoint slide

Linear discriminant analysis (LDA) is a method for demonstrating differences in complex datasets. Bars represent taxa that are dominant in a given group; the length of each bar represents the level of significance with LDA scores > 2 considered statistically significant. In experiment A ( Figure 3c ), there were significantly higher numbers of seven taxa in the FF group (green bars), ten taxa in the DF group (red bars), and nine taxa in the FF+Binf group (blue bars). In experiment B ( Figure 3d ), the dominant taxa in each group were dramatically different from those observed in experiment A.

Cladograms and principal coordinate analyses are complementary methods for illustrating differences between microbial communities. In the cladogram (Supplementary Figure S1 online), the innermost circle represents phylum, the next circle outward representing class, etc. Taxa that were more prominent in pups from a given group are color coded; the remaining taxa were present in at least one group and not significantly different between groups (tallow nodes in the cladograms). In experiment A, the FF group was dominated by γ-Proteobacteria, the DF group was dominated by Firmicutes (mostly Lactobacillales), and the FF+Binf group had features of both. The microbial composition was markedly different in experiment B with all groups showing a mix of γ-Proteobacteria and Firmicutes plus some Bacteroidetes (DF group) and Actinobacteria (FF+Binf group). Principal coordinate analyses based on Unifrac distances demonstrate the pronounced differences between groups in experiment A ( Figure 3e ) that is less apparent in experiment B ( Figure 3f ).

Figure 4 presents relative abundances of three taxa of interest for each group (panel a from experiment A and panel b from experiment B). The differences in the DF groups between the two experiments are dramatic. While the differences between the FF and the FF+Binf groups did not reach statistical significance, in both experiments there were consistent trends: a relative increase in Lactobacillaceae and decreases in Enterobactericeae and Enterococcaceae in the FF+Binf group. Bifidobacteria were not present in the ceca in significant numbers in any group. Rarefaction curves demonstrate a trend towards differing species richness and diversity between groups in experiment A, which is not seen in experiment B (Supplementary Figure S2 online).

Figure 4
figure 4

Relative abundance of Enterobacteriaceae, Lactobacillaceae, and Enterococcaceae for each animal exposed to hypoxia and cold stress (n=50). Solid lines represent mean values. Differences between all three groups were significant at P < 0.05 for Enterobacteriaceae in both experiments and for Lactobacillaceae and Enterococcaceae in experiment A. Differences between FF and FF+Binf were not significant for any of these three families in either experiment A or B. Circles = dam fed (DF), squares = formula fed (FF), triangles = formula with added B. infantis (FF+Binf).

PowerPoint slide

Linear regression of percentage of two dominant taxa (Enterobacteriaceae and Lactobacillaceae) with IL6, IL8, IL23, and iNOS expression for all 50 pups exposed to stress in both experiments demonstrated a weak correlation between iNOS and Enterobacteriaceae ( Figure 5a ). This association was stronger when analyzing just the stressed pups in experiment A ( Figure 5b ), but nonsignificant when analyzing the pups in experiment B (data not shown). A negative association between Lactobacillaceae and iNOS was observed in experiment A ( Figure 5c ) but not experiment B. To assess the impact of one outlier with very high iNOS expression, the analyses were repeated without including that pup ( Figure 5 , right panels).

Figure 5
figure 5

Ileal iNOS and cecal microbes. (a) Correlation of cecal Enterobacteriaceae with ileal iNOS for all animals exposed to hypoxia and cold stress (left panel R2 = 0.14, P < 0.01) and for all animals excluding one outlier (right panel R2 = 0.11, P = 0.02). (b) Correlation of cecal Enterobacteriaceae and ileal iNOS in all animals exposed to hypoxia and cold stress in experiment A (left panel R2 = 0.28, P=0.01) and excluding one outlier (right panel R2 = 0.43, P < 0.01). (c) Correlation of cecal Lactobacillaceae and ileal iNOS in all animals exposed to hypoxia and cold stress in experiment A (left panel R2 = 0.23, P=0.03) and excluding one outlier (right panel R2 = 0.41, P < 0.01).

PowerPoint slide

Discussion

Probiotics decrease the risk of NEC in premature infants, but the mechanisms of protection are unclear. Previous studies of probiotics in animal models of NEC show efficacy and support several possible mechanisms, including decreased apoptosis (9,22), improved intestinal integrity (8), and altered expression of antimicrobial peptides (10). In vitro studies have demonstrated that factors secreted by B. infantis suppress expression of innate immune response genes (including Tlr2, Tlr4, Il6, and Il8), resulting in decreased inflammation in immature, but not mature intestinal cells and tissues (23). A recent comparison of three Bifidobacterium and two Lactobacillus strains in the rat model of NEC found differences in protective ability among strains and no clear benefit of combinations over single strains (24). In vitro, B. infantis was superior to L. acidophilus in reducing innate immune gene expression and inflammation (23).

Our findings here support the hypothesis that B. infantis decreases the incidence and severity of NEC, and attenuates the hyper-inflammatory immune response of the immature intestinal mucosa as evidenced by decreases in expression of pro-inflammatory cytokines and Nos2 (iNOS). The results of experiment A are consistent with the hypothesis that a pro-inflammatory milieu stimulates an overgrowth of several species of Enterobacteriaceae (6,25); however, this association was not seen in experiment B.

Circumstantial evidence strongly supports a role for dysbiosis in the pathogenesis of NEC. In the mouse, experimental NEC is associated with increased Enterobacteriaceae (particularly Citrobacter and Klebsiella species) and decreased microbial diversity (12), whereas in the piglet experimental NEC is associated with increased surface-associated Clostridia species (13). Surveys of the fecal microbiota in premature infants demonstrate that dysbiosis due to increased Enterobacteriaceae is associated with NEC (7,26). Early dysbiosis may predispose premature infants to subsequent development of NEC (27), and increased administration of antibiotics (28,29) and acid-blocking agents (30,31) are both associated with dysbiosis and with increased susceptibility to NEC. In experiment A, we found marked dysbiosis in the FF group with high levels of Enterobacteriaceae that was partially attenuated by B. infantis; however, this pattern was not seen in experiment B. While the numbers are too small for conclusions, the data suggest that hypoxia and cold stress might not significantly alter the cecal microbiota, and that in the FF group the observed changes reflect either an effect of formula feeding or the stress of being housed without a surrogate dam (or both).

Despite the protective effects of probiotic clearly evident in the current study, the numbers of bifidobacteria in the cecum were very small in both experiments, suggesting that either the protective effects of this strain occurred “upstream” from the cecum, or perhaps more likely the beneficial effects require only minimal colonization. It would have been ideal to examine the microbiota of the distal ileum; however, the limited ileal tissue was used to assess NEC pathology, mRNA, and protein. Future studies in this model comparing the microbiota of the distal ileum and the cecum at days 1–4 of life would be valuable in further investigating the development of dysbiosis.

The marked differences between the two identical experiments underscore a variability in the cecal microbial community that remains unexplained, but worthy of comment. The two experiments were performed 4 mo apart with animals from the same facility receiving identical diets and interventions. We have previously noted that transportation of mice can cause a transient but significant alteration in the intestinal microbiota (loss of Lactobacilli) (32). We speculate that some unknown stressor increased the Enterobacteriaceae and decreased the Lactobacillaceae in the DF group in experiment B; it is striking that in spite of high numbers of Enterobacteriaceae, none of the pups in the DF group in experiment B developed NEC. Although we have previously demonstrated that infection with Helicobacter hepaticus (often asymptomatic in laboratory animals) worsens the severity of NEC in this model (11), we did not find significant numbers of H. hepaticus (phylum Proteobacteria, class ɛ-Proteobacteria) in either experiment A or B.

In spite of the marked differences in the microbiota between experiments A and B, the incidence and severity of NEC and the changes in markers of inflammation, antimicrobial peptides, and mucus-related molecules were very similar between experiments. This observation may suggest that changes in the microbiota in this model of experimental NEC are secondary to inflammation rather than directly causative of disease. An alternative explanation is that in the presence of significant stress, even subtle changes in the microbiota increase the risk of NEC.

Two key questions will require further investigation. First, is dysbiosis a cause or consequence of NEC disease? One scenario is that formula feeding and environmental stress trigger tissue inflammation, which in turn leads to NEC and to alterations in the gut microbiota, but the dysbiosis is of minimal consequence. Another scenario is that perturbation in host–microbe homeostasis induced by formula feeding and environmental stress causes a dysbiosis, and that the altered composition of the microbiota causes inflammation and mucosal disease. The data from this study tend to favor the first scenario, but are not definitive. Furthermore, we do not view the two scenarios as mutually exclusive, i.e., dysbiosis and inflammation may escalate each other, ultimately resulting in NEC.

The second key question is how does B. infantis prevent NEC? B. infantis might exert its effects primarily on host mucosa by suppressing detrimental inflammation (23) or by altering mucus or antimicrobial peptide expression. On the other hand, B. infantis might mediate its protective properties by attenuating dysbiosis via microbe–microbe interactions. Answers to these questions will likely help in developing more effective preventative and therapeutic interventions for NEC.

Methods

This protocol was approved by the Animal Care and Use Committee of the University of Arizona. Neonatal Sprague–Dawley rats (Harlan Laboratories, Madison, WI) were collected by Caesarean section 1 d before scheduled birth and their first feeding started 2 h after delivery. Rat pups were hand-fed five times daily with a total volume of 850 µl of rat milk substitute formula (11) (FF, n = 20), the identical formula plus 5 × 106 CFU per day of Bifidobacterium longum subsp. infantis ATCC 15697 (FF+Binf, n = 19), or dam-fed by surrogate mothers (DF, n = 20). The dose is the same as in previous experiments with a different bifidobacterium (8) and is comparable in CFU/g body weight to doses given to premature infants. Fifty animals (FF 17, FF+Binf 17, and DF 16) were exposed to asphyxia (breathing 100% nitrogen gas for 60 s) and cold stress (4 °C for 10 min) twice daily; the remaining nine animals were not exposed to asphyxia or cold stress. After 96 h, all of the surviving animals were killed by decapitation. Animals that developed signs of distress or imminent death before 96 h were killed and included in the study (n=1, DF).

NEC Evaluation

After euthanasia, a 2-cm piece of distal ileum was removed and fixed in 70% ethanol, paraffin embedded, sectioned at 4–6 µm, and stained with hematoxylin and eosin (H&E) for histological evaluation of NEC. Pathological changes in intestinal architecture were evaluated using our previously published NEC scoring system (8,33). Histological changes in the ileum were scored by a blinded evaluator and graded as follows: 0 (normal)—no damage; 1 (mild)—slight submucosal and/or lamina propria separation; 2 (moderate)—moderate separation of submucosa and/or lamina propria, and/or edema in submucosal and muscular layers; 3 (severe)—severe separation of submucosa and/or lamina propria, and/or severe edema in submucosa and muscular layers, regional villous sloughing; 4 (necrosis)—loss of villi and necrosis. Intermediate scores of 0.5, 1.5, 2.5, and 3.5 were also utilized to more accurately assess levels of ileal damage when necessary (34). Experimental NEC was defined as a histologic score of 2 or greater (8,35).

Morphometric Measurements in the Ileum

A 2-cm section of distal ileum stained with H&E was used for morphometric measurements as previously described (34). Briefly, 20 villi were measured in each histological sample and 8–10 animals were evaluated per experimental group. Sections from animals with a NEC score of 3 and higher were not included in analyses because of the lack of intact tissue to evaluate. Villi were measured from the tip to the crypt base using an image analysis system (Image-Pro Plus; Media Cybernetics, Silver Spring, MD) in a blind manner to prevent observer bias.

RNA Preparation and Real-Time Polymerase Chain Reaction

Total RNA was isolated from ileal tissue using the RNeasy Mini Kit (Qiagen, Santa Clarita, CA). RNA concentration was quantified by ultraviolet spectrophotometry at 260 nm using a NanoDrop (Thermo Fisher Scientific, Wilmington, DE). cDNA synthesis and real-time PCR were performed as previously described (10). The PCR primers for Reg3b (NM053289, also referred to as Pap1) and Reg3g (NM 173097, also referred to as Pap3) have been previously reported (10). TaqMan primers and probes were used for the detection of Nos2, Cxcl1, Il6, Il23, Muc2, Tff3, Tlr2, and Tlr4. Reporter dye emission was detected by an automated sequence detector combined with ABI Prism 7700 Sequence Detection System software (Applied Biosystems, Foster City, CA). Real-time PCR quantification was then performed using TaqMan 18S controls.

Western Blot Analysis

Ileal tissue was examined for protein expression of TLR2 and TLR4 as well as iNOS as previously described (9). The following antibodies (all from Santa Cruz Biotechnology, Dallas, TX) were used: goat anti-TLR-4 polyclonal antibody (1:500, sc-3002), goat anti-TLR-2 polyclonal antibody (1:500, sc-16237), or mouse anti-NOS2 monoclonal antibody (1:500, sc-7271).

Microbiota Analysis

To ensure accurate detection of bifidobacteria, preliminary analyses were performed on unrelated specimens with and without the additional steps of bead beating and addition of lysis buffer (Supplementary Figure S3 online). Consistent with a previous report, these results underline the necessity of both steps to ensure detection of bifidobacteria (36). For all microbiota studies, intact ceca were placed in 5 ml of RNAlater solution (Qiagen) (37), kept at 4 °C overnight, and then stored at −80 °C until processing. Ceca samples were thawed and DNA was extracted using the QIAamp DNA Stool Mini Kit (Qiagen) with the additional enzymatic lysis and homogenization by bead beating (36).

Illumina Sequencing

Amplification and sequencing were performed as described previously (38). Briefly, the V4 domain of bacterial 16S rRNA genes was amplified using primers F515 (5′-NNNNNNNN GTGTGCCAGCMGCCGCGGTAA-3′) and R806 (5′-GGACTACHVGGGTWTCTAAT-3′) (39), with the forward primer modified to contain a unique 8 nt barcode (italicized poly-N section of the primer above) and a 2 nt linker sequence (bold, underlined portion) at the 5′ terminus. PCR reactions contained 5–100 ng DNA template, 1X GoTaq Green Master Mix (Promega, Madison, WI), 1 mmol/l MgCl2, and 2 pmol of each primer. Reaction conditions consisted of an initial 94 °C for 3 min followed by 35 cycles of 94 °C for 45 s, 50 °C for 60 s, and 72 °C for 90 s, and a final extension of 72 °C for 10 min. PCR amplicons were pooled at roughly equal amplification intensity ratios, purified using the Qiaquick PCR purification kit (Qiagen), and submitted to the UC Davis Genome Center DNA Technologies Core for Illumina paired-end library preparation, cluster generation, and 250 bp paired-end sequencing on an Illumina MiSeq instrument.

Sequencing Analysis

The QIIME software package (University of Colorado, Boulder, CO, version 1.7.0) was used to analyze data from the sequencing run (39,40). Sequences were quality filtered and demultiplexed, and then operational taxonomic units (OTUs) were assigned using UCLUST (drive5.com, Tiburon, CA) (41), based on 97% pairwise identity as previously described (38). OTUs went through a secondary filtration by 0.005% to remove low-abundance OTUs (38). Taxonomic classification of the filtered OTUs was based on the Ribosomal Database Project classifier (Michigan State University, East Lansing, MI) (42) against a representative subset of the Greengenes 16S rRNA database (Second Genome, South San Francisco, CA, gg_13_5 release) (43). OTU sequence alignment was carried out with PyNAST (University of Colorado) (38,44) and used to build a phylogenetic tree for β diversity analyses. β diversity was estimated by calculating abundance-weighted and unweighted UniFrac (45) distances. Samples were clustered based on between-sample distances.

The linear discriminate analysis (LDA) effect size (LEfSe) (46) module (Galaxy, Boston, MA) for biomarker discovery was used to determine taxonomically significant differences between the DF, FF, and FF+Binf groups. LEfSe uses a factorial Kruskal–Wallis sum-rank test (α = 0.05) to identify taxa with significant differential abundances between categories (using one-against-all comparisons), followed by a Wilcoxon Mann–Whitney test to assess subclass variation, and LDA to determine the effect size of each differentially abundant taxon.

Statistical Analysis

DF, FF, and FF+Binf groups were compared using ANOVA followed by Fisher PLSD and by the Student’s t-test at the 95% CI. Analysis of NEC score between groups was accomplished using the Kruskal–Wallis test for nonparametric values followed by pairwise comparison using the Mann–Whitney test. The Pearson’s chi-squared (χ2) test was used to analyze differences in incidence of NEC between groups. All statistical analyses were conducted using the statistical program StatPlus:mac LE for Macintosh computers (AnalystSoft, Alexandria, VA). All numerical data are expressed as mean ± SE.

Statement of Financial Support

UC Davis, Davis, CA: National Institutes of Health (NIH) HD059127 (M.A.U., C.L.B., D.A.M.), NIH AI32738 (C.L.B.), NIH AT007079 (D.A.M.), The Development and Promotion of Science and Technology Talents Project, Thailand (A. Kananurak). U Arizona, Tucson, AZ: NIH HD039657, Mead Johnson, and Meiji Dairies Co (B.D.).

Disclosure

There are no conflicts of interest.