Patterned progression of gut microbiota associated with necrotizing enterocolitis and late onset sepsis in preterm infants: a prospective study in a Chinese neonatal intensive care unit

Ethics

This study was approved by the joint committee of ethics of Shanghai Children’s Medical Center, Shanghai Jiao Tong University School of Medicine (SCMCIRB-K2013022). Detailed written informed consent was obtained from parents before enrolment.

Patients

Newly born preterm infants with a gestational age less than 33 weeks and birth weight over 950g were enrolled from Neonatal Intensive Care Unit (NICU) at Shanghai Children’s Medical Center from July 2013 to December 2014. The exclusion criteria were (1) diagnosed with early-onset sepsis, (2) hepatic diseases, (3) renal impairment (Cr > 88 µM), 4) diagnosed with intestinal obstruction, (5) in foreseeable need of major cardiovascular or abdominal surgeries (except for male circumcision or PDA ligation), (6) estimated parenteral support to supply over 50% of daily caloric intake for more than four days, (7) given intravenous antibiotics administration (except prophylactic regimen of cefotaxime, piperacillin-tazobactam and/or metronidazole), (8) history of oral antibiotics administration, (9) grossly bloody stools at admission, and (10) over five days old.

NEC cases were defined as infants who met the criteria for Stage II and Stage III NEC diagnosis (Bell et al., 1978), including radiographic intestinal dilation, ileus, pneumatosis intestinalis, and/or absent bowel sounds with or without abdominal tenderness, and/or mild metabolic acidosis and thrombocytopenia. An LOS case was defined if an infant (1) had a positive hemoculture or other suspicious loci of infection after 72 h of life, or (2) presented with septic signs/symptoms reviewed and diagnosed independently by at least two neonatologists, and had been responding well with advanced antibiotics (e.g., Meropenem) after diagnosis. Infants with no infectious complications were regarded as controls.

Sample collection and handling

Fecal sample collection started from neonatal meconium until death or discharge, whichever came first. Although we intended to collect fecal samples every day, due to working shifts and flexible clinical scheduling, we set seven days as the maximum interval between two collections from every infant. Every sample was collected from infants’ diaper with a sterile spatula into cryogenic vials within 30 min of defecation. Then the sample was immediately placed on dry ice and stored at −80 °C within 30 min without additives. All samples were collected and stored before knowing the diagnosis of respective patients.

DNA extraction and quality control amplification and 16s rRNA gene sequencing

Microbial genomic DNA was isolated from each fecal specimen using the E.Z.N.A. Soil DNA Kit (Omega Bio-Tek, Norcross, GA, U.S.) according to the manufacturer’s protocols. The concentration and purity of the DNA were determined by NanoDrop 2000 UV-vis spectrophotometer (Thermo Scientific, Wilmington, DE, USA), and the DNA quality was checked by 1% agarose gel electrophoresis.

Broad-range PCR and High-throughput Sequencing of 16s rRNA gene amplicons

The V3-V4 hypervariable regions of the bacterial 16S rRNA gene were amplified by PCR from each sample using bacterialarchaeal primers 338F (5′-ACTCCTACGGGAGGCAGCAG-3′) and 806R (5′-GGACTACHVGG GTWTCTAAT-3′) using a thermocycler PCR system (GeneAmp 9700, ABI, USA). The PCR reactions were as follows: 3 min of denaturation at 95 °C, 27 cycles of 30 s at 95 °C, 30 s annealing at 55 °C and 45 s elongation at 72 °C, and a final extension at 72 °C for 10 min. The PCR reactions were performed in triplicate, with each 20 µL mixture containing four µL 5X FastPfu Buffer, two µL 2.5 mM dNTPs, 0.8 µL of each primer (five µM), 0.4 µL FastPfu Polymerase (TransGen Biotech, Beijing, China) and 10 ng template DNA. PCR products were separated from impurities and genomic DNA by running in 2% agarose gels. The PCR bands were further purified using the AxyPrep DNA Gel Extraction Kit (Axygen Biosciences, Union City, CA, USA), and quantified using QuantiFluor-ST (Promega, USA) according to the manufacturer’s protocols. Equimolar amounts of purified amplicons were pooled and paired ended sequenced (2 × 300) on an Illumina MiSeq platform (Illumina, San Diego, USA) according to the standard protocols of Majorbio Bio-Pharm Technology Co. Ltd. (Shanghai, China). The reads were de-multiplexed using the Illumina software and separate FASTQ files were generated for each specimen and deposited to the Sequence Read Archive NCBI under the BioProject accession PRJNA470548. Another public archive repository is available at Figshare doi: https://doi.org/10.6084/m9.figshare.7205102.

Raw data processing

Raw data were processed according to the standard protocols provided by Majorbio Bio-Pharm Technology Co. Ltd. (Shanghai, China) as previously described (Liu et al., 2018; Wang et al., 2018). In short, raw sequencing data was first de-multiplexed. Sequence reads were then subjected to quality filtering utilizing Trimmomatic software (Bolger, Lohse & Usadel, 2014) and were truncated at any site with a Phred score <20 over a 50 bp-sized window. Barcode matching with the primer mismatch from 0 to 2 nucleotides was adopted and reads containing ambiguous characters were removed. After trimming, FLASh (Fast Length Adjustment of Short Read) (Magoč & Salzberg, 2011), a read pre-processing software, assembled and merged the paired-end reads from fragments and generated >10 bp overlapped, with the dead match ratio of 0.2. Unassembled reads were discarded. From the 192 fecal samples sequenced, a total of 7,472,400 optimized V3-V4 tags of 16s rRNA gene sequences were generated (Table S1).

To unbiasedly compare all the samples at the same sequencing depth, the ”sub.sample” command of mothur program (version1.30.1) (Schloss et al., 2009) was used for normalization to the smallest sample size. Chimera was detected and removed by UCHIME Algorithm. The effective reads were then sorted by cluster size and processed using Operational Taxonomic Units (OTUs) with 97% similarity cutoff UPARSE-OTU algorithm (implementing “cluster_otus” command) (Edgar, 2013) in USEARCH(v10)(UPARSE version 7.1). The taxonomy of each 16S rRNA gene sequence was analyzed by RDP Classifier algorithm (Wang et al., 2007) against the Silva (SSU128) 16S rRNA database (Quast et al., 2012) using confidence threshold of 70%. Each sequence was assigned the taxonomy by QIIME (Caporaso et al., 2010). The representative sequences were allocated phylogenetically down to the domain, phylum, class, order, family, and genus levels (Table S2). The relative abundance of a given taxonomic group was calculated as the percentage of assigned sequences over total sequences.

Within-sample diversity (alpha diversity), including Shannon index and observed species richness (Sobs), was obtained using the “summary.single” command of mothur program (version1.30.1) (Schloss et al., 2009). Between-sample diversity (beta diversity) was obtained by calculating weighted UniFrac distances between samples.

Statistical and bioinformatics analyses

Patients characteristics

From July 2013 to December 2014, a total of 130 preterm infants admitted to the neonatal intensive care unit (NICU) of Shanghai Children’s Medical Center met the criteria of our study and a total of 1698 samples were collected. 192 fecal samples from 24 well-sampled preterm infants were sequenced. Four subsequently developed NEC (2 in stage IIA and 2 in stage IIB) and three developed LOS (2 with positive hemoculture of Klebsiella pneumoniae; the other one was diagnosed upon sepsis-related signs and symptoms, lab test of white blood cells > 20 cells/microL and her effective reaction to vancomycin). The remaining 17 served as matched controls (Fig. 1, Table S3). Fecal samples were collected between days 1 and 69 of life. Numbers of samples collected and interval of sampling varied among patients but met our preset criteria of less than 7 days between sampling. The average number of sample collected for NEC, LOS and control patients was 11, 14 and six respectively. The number of samples per patient was higher in the NEC and LOS groups because the severity of the disease required longer hospitalization (p = 0.046).

All 24 infants profiled were delivered by Cesarean section, fed on infant formula and prescribed with prophylactic antibiotics regimen (cefotaxime, piperacillin-tazobactam and/or metronidazole) right after they were admitted to our NICU. No infant was prescribed probiotics during the study. There was no significant difference in gestational age (p = 0.074), birth weight (p = 0.111) or gender proportions (p = 0.822) among the three groups. The average age at diagnosis for both disease groups was 16 days and there was no statistical difference between the groups (p = 0.629) (Table 1). Therefore, we assigned day 16 to discharge as early disease interval, day 4–8 as early pre-onset interval and day 9–15 as late pre-onset interval for the control group (Table S4).

Longitudinal Microbiome Diversity of NEC and LOS patients

To get an overview of gut microbiota in patients, we analyzed the microbial richness of the NEC and LOS patients over time. Similar to the control group, the case groups showed a decreasing trend in observed species (Sobs) from early post-partum stage (EPP) to early disease (ED) stage (Fig. 2 (A) control group, p < 0.01; (B) NEC group, p = 0.044; (C) LOS group, p = 0.013; Dataset S1, Sheet “Sobs” two way RM ANOVA, p < 0.0001). The greatest decline in sobs was from early pre-onset (EPO) to late pre-onset (LPO). However, the decrease in the disease groups was less significant than the control group (control group p = 0.0004, NEC group p = 0.18, LOS group p = 0.066). The Sobs then stabilized from LPO onward with no significant difference between adjacent time intervals.

Next, we analyzed gut microbiome evenness over time. Similar to Sobs, the Shannon indices decreased significantly from the early post-partum (EPP) to early disease (ED) stage (Fig. 3A control group 2.768 to 1.004, p = 0.04; (B) NEC group, 3.141 to 0.578, p = 0.01; (C) LOS group, 2.641 to 0.470, p = 0.01).

Two way RM ANOVA showed significant Shannon index divergent among three groups before disease onset (Dataset S1, sheet “Shannon”, EPP to ED, p = 0.0017). Moreover, during early disease stage, the Shannon indices were different among three groups (Fig. 4, facet early disease, p = 0.0037), suggesting that microbiota distortion may precede NEC and LOS onset. As diseases progressed, the NEC group differed significantly with the LOS group during middle disease interval but insignificantly during late disease interval (Fig. 4, facet middle disease, p = 0.034; facet late disease, p = 0.750). Upon alleviation of both diseases, the Shannon indices rose back to the early pre-onset levels (Fig. 3B) NEC group. early pre-onset at 1.925 vs. post disease at 1.320, p = 0.79; (C) LOS group, early pre-onset at 2.473 vs. post disease at 1.463, p = 0.16).

Kinetics of microbiome composition

To compare the beta-diversity of the three groups over time, we applied Principal Component Analysis (PCoA) to weighted UniFrac distance matrix (Clarke, 1993). Bacterial composition of three groups during early post-partum interval were the most similar compared with other time intervals, with the first principal coordinates accounted for 33.01% (Fig. 5A). Then beta diversity continued to separate from one another. The first principal coordinate one (PC1) increased from 33.01% at the early post-partum to 35.23% at the early pre-onset stage, 38.36% at the late pre-onset stage and eventually reaching 42.32% at the early disease stage (Figs. 5B to 5D). This continuous increase in beta-diversity suggested that the phylogenetic composition of the patients’ microbiome started to deviate from the control group before the onset of diseases. As diseases progressed, the phylogenetic similarity between the NEC and the LOS disease groups diverged further and peaked at 59.53% in middle disease stage then came down gradually to 42.8% at post disease stage (Fig. 5E to Fig. 5G). This trend in phylogenetic dissimilarity suggested that the microbiome composition of the NEC and LOS patients might have deviated from normal even before the onset of diseases. Also, the further separation between the NEC and the LOS groups could be a result of different treatment strategies.

Colonization trend at the genus level

In the analyses of intestinal microbiome alpha (Figs. 2–4) and beta diversity (Fig. 5), detectable differences were observed among the three groups, especially during the transition from the LPO to ED stage. This indicated that the microbiota assembly differences between the case groups and control group. To further investigate which microbiota composition was correlated with the onset and/or progression of NEC and LOS, we tracked the longitudinal compositional changes in genera abundance. We filtered the genus of over 10% relative abundance among all samples and plotted relative abundance over time (Fig. 6).

At the early post-partum stage, all three groups showed high proportion of Lactococcus, Bacillus and Pseudomonas. However, ZIBR model the disease groups showed significantly higher OTUs that matched to Bacillus (NEC 15.05% and LOS 15.97% compared to 6.02% of control, p = 0.032) and Solibacillus (8.88% in NEC and 9.61% in LOS compared to 3.65% of control, p = 0.047) from the case groups (Dataset S2). Moreover, Enterococcus proportion (Fig. 6B, purple area) was much higher in LOS patients (20.72%) than the normal controls (6.66%, Fig. 6A, purple area) but almost absent in NEC patients (0.51%) (Fig. 6B). While all three groups showed increases in Klebsiella and Escherichia-Shigella and decrease in Lactococcus from EPP to ED, the rates of change were different among the three groups. The LOS group exhibited the most drastic changes, with a rapidly increase of Klebsiella (from 4.71% to 58,90%), Escherichia-Shigella (from 2.02% to 18.16%) and Streptococcus (from 1.22% to 12.68%) (Fig. 6C). Together, these three genera accounted for almost 100% of all bacteria (Fig. 6C). In addition, Lactococcus decreased more rapidly than the other groups, from 24.54% at EPP to 0.94% before LPO (Fig. 6C magenta area).

Besides, the increase of Klebsiella was the most minimal in NEC patients (Fig. 6B grey area, from 7.17% at EPP to 35.63% at ED). Moreover, a rapid surge of Enterococcus, Staphylococcus and Streptococcus from EPO to ED was only observed in NEC patients (Fig. 6B, purple, dark and light blue area).

As NEC and LOS progressed with medical intervention, the genus in case groups underwent another round of drastic changes. Most notably, the fluctuation of Enterococcus, Klebsiella, Staphylococcus and Peptoclostridium during the disease stages (Figs. 6B and 6C stage ED to LD), which might be resultant from different healthcare strategies applied in two groups. Interestingly, as patients approached remission, the composition became more balanced and resembled more to that of the normal control, except for a higher level of Clostridium. In summary, relative to patients in the control group, we observed different patterns of temporal alterations in bacterial composition among NEC and LOS patients. Rapid changes in relative abundance of certain genera were revealed as early as early pre-onset of stages and were the most notable in LOS patients.