Gram-negative Microbiota Blooms in Premature Twins Discordant for Parenteral Nutrition-associated Cholestasis

Abstract

Parenteral nutrition-associated cholestasis (PNAC) causes serious morbidity in the neonatal intensive care unit. Infection with gut-associated bacteria is associated with cholestasis, but the role of intestinal microbiota in PNAC is poorly understood. We examined the composition of stool microbiota from premature twins discordant for PNAC as a strategy to reduce confounding from variables associated with both microbiota and cholestasis. Eighty-four serial stool samples were included from 4 twin sets discordant for PNAC. Random Forests was utilized to determine genera most discriminatory in classifying samples from infants with and without PNAC. In infants with PNAC, we detected a significant increase in the relative abundance of Klebsiella, Veillonella, Enterobacter, and Enterococcus (P < 0.05). Bray-Curtis dissimilarities in infants with PNAC were significantly different (P < 0.05) from infants without PNAC. Our findings warrant further exploration in larger cohorts and experimental models of PNAC to determine if a microbiota signature predicts PNAC, as a basis for future interventions to mitigate liver injury.

What Is Known

• Parenteral nutrition-associated cholestasis causes significant morbidity.
• Infection with certain gut-associated bacteria has been associated with cholestasis.
• There is a paucity of data investigating the role of intestinal microbiota in parenteral nutrition-associated cholestasis.

What Is New

• Significant alterations in gut microbiota diversity and overrepresentation of certain Gram-negative bacteria were associated with parenteral nutrition-associated cholestasis in twin analyses.
• Findings may form the basis for further exploration to determine if a microbiota signature predicts parenteral nutrition-associated cholestasis.

Parenteral nutrition (PN) is an essential form of nutritional support for premature infants, yet is frequently complicated by its association with liver disease. Preterm infants have increased risk of developing PN-associated cholestasis (PNAC) because of immaturity of the hepatobiliary system compounded by lack of enteral feeding ⁽¹⁾, with cholestasis often persisting after PN cessation ⁽²⁾. Thus, PNAC causes significant morbidity in the neonatal intensive care unit (NICU) and beyond ⁽³⁾. Clinical factors associated with PNAC include low birth weight, longer duration of PN, delayed enteral nutrition, antimicrobial administration, Gram-negative infections, necrotizing enterocolitis, and gastrointestinal malformations ⁽⁴⁾.

Altering the composition of PN and lipid-sparing strategies have shown partial success to prevent or alleviate the development of cholestasis ⁽⁵⁾. More recently, studies have focused on the role of the intestinal microbiota in the development of PNAC. PN has been associated with alteration of the intestinal microbiota ⁽⁶⁾, although this can be difficult to separate from the effect of lack of enteral nutrition on the microbiota. Evidence also suggests that cholestasis may be related to the composition of the microbiota ^(7,8). Bile acids, diminished in cholestasis, normally have an antimicrobial effect on the microbiota; conversely, the intestinal microbiota metabolizes and biotransforms bile acids ⁽⁸⁾. In adults, perturbations to the intestinal microbiota have been associated with worsening liver function, cirrhosis, and primary biliary cholangitis, with dysbiosis worsening as liver disease progresses ^(9,10).

In a prior study of infants with short bowel syndrome receiving PN, infants who developed PN-associated liver disease displayed less diverse microbiota and a higher abundance of pathogenic bacteria ⁽¹¹⁾. Moreover, infection with Gram-negative bacteria, including Escherichia coli urinary tract infections or Klebsiella septicemia, has been associated with the development of cholestasis ^(9–13).

Early enteral feeding has been shown to be the most effective maneuver to prevent the development of cholestasis ⁽¹⁴⁾; however, initiation of enteral feeding in the NICU is often restricted by clinical factors. Other targets to decrease the burden of PNAC in the NICU are needed and the intestinal microbiota provides a compelling target. Due to methodological challenges in comparing infants with different genetics, diets, and environments, little is known about changes in the intestinal microbiota associated with PNAC in NICU patients. Therefore, the goal of this study was to investigate intestinal microbiota composition from serial stool samples of premature twins discordant for PNAC.

METHODS

Subjects and Sample Collection

Subjects were chosen from a larger cohort in 2 ongoing neonatal microbiome studies from Level IV NICUs at Inova Fairfax Hospital, VA (IRB approval #15-1945) and the University of Virginia (IRB approval #18244). Neonates with an anticipated length of stay over 5 days were recruited. Detailed maternal, pregnancy, and delivery data were collected. While in the NICU infants had stool collected at least twice weekly and frozen at −80 °C within 12 hours. Detailed data regarding feeding, medications, and health status were collected.

Twin sets were selected for discordance for PNAC, that is, twin pairs simultaneously receiving PN but only 1 twin developed PNAC. Selected subjects were admitted to the NICU for prematurity only; no subject developed necrotizing enterocolitis. Subjects were not included from the larger cohort if they were singletons, had cholestasis from causes other than PNAC, were admitted to the NICU for causes other than prematurity, had congenital abnormalities or developed necrotizing enterocolitis. PNAC was defined as a direct bilirubin level of ≥1 mg/dL, deemed to be caused by PN by the treating physician, with no other apparent causes for cholestasis ⁽¹⁵⁾. Tests conducted to screen for other causes of cholestasis included a right upper quadrant ultrasound, Alpha-1 Antitrypsin deficiency testing, urine CMV PCR, urine culture, and thyroid-stimulating hormone; specific tests conducted varied for each cholestatic subject. Twin sets were chosen as a strategy to reduce potential confounding variables associated with both the microbiota and cholestasis including age, intrauterine environment, genetics [although the environment has been shown to dominate the shaping of the microbiome more than genetics ⁽¹⁶⁾] and maternal breast milk.

DNA Extraction and 16S Ribosomal RNA Gene Sequencing

Whole stool samples were stored at -80 °C until DNA extraction. Prepped samples were loaded on the EZ1 Advanced (Qiagen, Valencia, CA) using the EZ1 Tissue Kit and the Bacterial DNA Extraction protocol card. Samples were cleaned and concentrated using the DNeasy PowerClean Cleanup kit (Qiagen).

Sequencing libraries were prepared using a Nextera XT kit (Illumina, San Diego, CA) using a modified Illumina 16S Metagenomics Sequencing Library Preparation protocol for analysis of hypervariable region V4. Each sample was sequenced on the Illumina MiSeq with paired-end reads of 301 bp. Sequencing of negative controls of lysis buffer and positive controls of Staphylococcus aureus (Strain NCTC 8532, ATCC, VA) and E coli (Strain NCTC 9001, ATCC, VA) were included. The sequence data are being deposited in NCBI Sequence Read Archive (SRA).

16S Ribosomal RNA Gene Analysis

Sequence analysis was performed using DADA2 (version 1.6.0) ^(17,18). The forward read was truncated to 200 base pairs and reads with ambiguous ‘N’ bases, and >2 expected errors were filtered out. Chimeras were removed. Taxonomy was assigned using the Ribosomal Database Project's naïve Bayesian classifier ⁽¹⁹⁾ with RDP training set 16. Resulting sequence variant counts (SVs) and taxonomic assignments were analyzed using several R packages: Phyloseq (v. 1.22.3), Vegan (v. 2.4.6), GGplot2 (v. 2.2.1), and mice (v. 2.46.0) ^(20,21). Ultralow abundance SVs (<2 raw reads and present in <5% of samples) were removed. Missing values in the metadata were imputed using mice. Benjamini-Hochberg-corrected Wilcoxon rank sum tests were used to compare the alpha diversity and genus-relative abundances. A PERNANOVA (Adonis test in the Vegan R package) was used to compare the sample clusters displayed in the nonmetric multidimensional scaling (NMDS) ordination plot.

Random Forest Model Generation

The R package AUCRF (v. 1.1) was used to generate a Random Forest (RF) model to classify the samples according to the disease state from which the individual each sample was sourced. Relative SV counts, agglomerated to the genus taxonomic level, were used to train the RF. During the creation of an optimal RF, 5-fold cross validation with 20% of the data left out for validation was used to assess the fit of the model. The forest used has 1000 trees with a node split (mtry) of 2, and the classes were manually weighted to account for the imbalance between the PNAC and non-PNAC classes. Features were selected by maximizing model accuracy and minimizing the features required for classification. All code is available on Github (https://github.com/Tjmoutinho/PNAC_microbiome).

RESULTS

Demographics and Clinical Characteristics

A total of 84 serial stool samples were included from 4 twin sets discordant for PNAC. Twins ranged from 25 to 31 weeks’ gestational age (mean = 27 weeks). All twin sets were delivered by emergency Cesarean Section. Each twin set has 1 subject that was diagnosed with PNAC for part of the sampling period. Birth weights did not differ between twins who did and did not develop cholestasis (P > > 0.05, Wilcoxon rank sum test); see birth weights in Supplementary Table 1 (Supplemental Digital Content, https://links.lww.com/MPG/B769). There was no significant difference in antibiotic exposure between twins who did and did not develop cholestasis. Both members from 1 twin set (twin set 3) were treated for sepsis, with negative blood and urine cultures, and hence had more days of antibiotic use than the other twin sets. No other subject was diagnosed with sepsis or an infection; all blood cultures drawn on all subjects during the study were negative. Twins who developed cholestasis trended towards being on PN longer than their sibling without cholestasis (Table 1) because of slower tolerance of advancement of enteral feeds or less satisfactory weight gain, although this difference was not significant; however, cholestasis in the affected twin in all twin sets developed while both twins were still on PN. When enteral feeding was started, all subjects received maternal breast milk. Demographic and clinical characteristics of twin sets are listed in Table 1. The range of peak cholestasis in the twins with PNAC was 1.9 to 10.1 mg/dL, mean = 5 mg/dL; see hepatic function panel values for subjects in Supplementary Table 2 (Supplemental Digital Content, https://links.lww.com/MPG/B769).

TABLE 1: Demographic and clinical characteristics of twin sets

Parenteral Nutrition-associated Cholestasis Stool Samples Were Different From Nonparenteral Nutrition-associated Cholestasis Samples; 5 Significantly Differentially Abundant Genera Were Driving This Difference

All samples were classified by the clinical PNAC diagnosis at sampling. PNAC and non-PNAC samples cluster separately (P < 0.001, PERMANOVA; Fig. 1A). There were 5 statistically significant differentially abundant genera between PNAC and non-PNAC samples (P < 0.05, corrected Wilcoxon rank-sum; Fig. 1B). PNAC was associated with increased levels of Gram-negative bacteria Klebsiella, Veillonella, and Enterobacter, increased levels of Enterococcus and decreased levels of Escherichia/Shigella. Klebsiella, Enterococcus, and Veillonella have the highest median relative abundances and tended to be more abundant in PNAC samples.

FIGURE 1: (A) Microbiota 16S rRNA gene sequence data is visualized using nonmetric multidimensional scaling (Bray-Curtis) to compare samples. The PNAC sample cluster is different from the non-PNAC samples cluster (P < 0.001, PERMANOVA). (B) Differentially abundant taxa in samples with and without PNAC are starred (^∗) (Benjamini-Hochberg corrected P value <0.05). (C) Alpha diversity (Shannon Diversity Index) was not significantly different between PNAC versus non-PNAC. The sample with the greatest direct bilirubin level per subject is starred (^∗); the peak direct bilirubin was 1.9, 3.4, 4.6, and 10.1 mg/dL for the cholestatic twin of sets 1 to 4, respectively. (D) Relative abundance of microbial genera is displayed for each sample. (E) AUCRF was used to select the optimal random forest model with the lease number of predictive features. The optimal RF model generated for classifying the microbiota samples had an accuracy of 90%. Using 5-fold cross validation, with 1/5 of the samples for a validation set, the average accuracy was 88%, indicating that the optimal model was not over fit to the dataset. The 5 microbial features (genera) used in the optimal RF model are displayed. Each feature is listed in order of descending importance to the model. AUCRF, Area Under the ROC Curve Random Forest; PNAC = parenteral nutrition-associated cholestasis; RF = Random Forest.

Alpha diversity (Shannon index) was not significantly associated with PNAC. In several subjects, the alpha diversity increased with age, but this was not seen in all individuals (Fig. 1C). Five genera dominated most samples (Staphylococcus, Escherichia/Shigella, Klebsiella, Veillonella, Enterococcus) consisting of 2 phyla (Proteobacteria and Firmicutes) (Fig. 1D).

Random Forest Model Predicts Parenteral Nutrition-associated cholestasis Based on Relative Abundance of Specific Genera

Feature reductions were performed using AUCRF (Area Under the ROC Curve Random Forest) to identify the microbial features that were predictive of PNAC versus non-PNAC samples. Feature reduction involved starting with all genera with at least 1 nonzero median value (Fig. 1B). The optimal RF model consisted of 5 genera and predicted sample diagnosis with 90% accuracy (Fig. 1E). Five-fold cross validation indicated that the model was not overfit to the dataset. The high predictive capability of the optimal RF model suggests a strong microbial signature was associated with PNAC.

DISCUSSION

We detected significant differences in the stool microbiota of premature twins discordant for PNAC. Specifically, the relative abundance of genera differed, including increased levels of Gram-negative bacteria Veillonella, Klebsiella, and Enterobacter in infants with PNAC. This finding is intriguing, given well known associations of Gram-negative infections in neonates and the development of cholestasis ^(12,13) and raises the question of whether increased intestinal colonization alone, rather than frank infection, with these microbes is sufficient to provoke cholestasis. Indeed, increased intestinal Enterobacteriaceae⁽⁹⁾ has been described in adults with cirrhosis, whereas increased Klebsiella and Veillonella has been described in adults with primary biliary cholangitis ⁽¹⁰⁾. A predominance of Enterobacteriaceae has also been identified in infants with short bowel syndrome and PN-associated liver disease ⁽¹¹⁾. The actual species and strains associated with PNAC will need to be identified in future studies, using techniques, such as shotgun metagenomic sequencing.

The extent to which these changes in the relative abundance of genera occur before cholestasis versus a result of cholestasis cannot be answered in this pilot study as cholestatic subjects from 2 twin sets (1 and 2) developed cholestasis before the first stool sample was collected. In 1 twin set ⁽⁴⁾, however, the subject who developed PNAC displayed an increase in Veillonella abundance relative to their non-PNAC twin at the same time point. This question is of the utmost importance to answer, because of the significant potential utility of being able to identify, which infants might be susceptible to PNAC before clinical evidence of PNAC develops. Early interventions to mitigate the development of PNAC could then be initiated in at risk infants, for example, early enteral feeding and the use of fish oil-based lipid emulsions ⁽⁵⁾. Although this approach would be difficult to test in human neonates, gnotobiotic models of PNAC might offer 1 path forward for testing whether specific Gram-negative bacteria are predictive, or even causal, in the development of PNAC ⁽²²⁾.

Certainly, the ability to predict the development of cholestasis based on a microbiota signature would be powerful. In our study, a predictive model of PNAC was created based on 5 genera, providing a method for identifying informative microbial signatures using RF. This model needs to be further trained and tested on larger datasets to ensure that it remains valid in a broader population of at-risk infants. In addition, testing it in at-risk infants who have not yet developed PNAC will be imperative. If increased levels of certain Gram-negative bacteria are found to be causal or compounding in PNAC, targeted interventions, such as antibiotics, probiotics, or fecal microbiota transplantation could be investigated.

Limitations of this study include the small number of infant twin pairs. This limited sample size is, however, offset somewhat by the use of discordant twin pairs, potentially reducing confounding from the intrauterine environment, age, genetics, and maternal breast milk, and rigorous analysis of serial stool samples. An additional limitation is that tests used to exclude other causes of cholestasis in each cholestatic subject were not homogeneous and were dependent on the clinical situation.

In conclusion, we found marked differences in the diversity and relative abundance of stool microbiota in infants who developed PNAC versus twin controls, most notably an increase in Gram-negative bacteria Veillonella, Klebsiella, and Enterobacter. This finding warrants further exploration in larger cohorts and experimental models to determine if signals in the microbiota are predictive of PNAC, thereby providing a basis for interventions to mitigate the development of PNAC in vulnerable infants.