Sphingolipid metabolism potential in fecal microbiome and bronchiolitis in infants: a case–control study

Abstract

Objective

Emerging evidence demonstrated that the structure of fecal microbiome is associated with the likelihood of bronchiolitis in infants. However, no study has examined functional profiles of fecal microbiome in infants with bronchiolitis. In this context, we conducted a case–control study. As a part of multicenter prospective study, we collected stool samples from 40 infants hospitalized with bronchiolitis (cases). We concurrently enrolled 115 age-matched healthy controls.

Results

First, by applying 16S rRNA gene sequencing to these 155 fecal samples, we identified the taxonomic profiles of fecal microbiome. Next, based on the taxonomy data, we inferred the functional capabilities of fecal microbiome and tested for differences in the functional capabilities between cases and controls. Overall, the median age was 3 months and 45% were female. Among 274 metabolic pathways surveyed, there were significant differences between bronchiolitis cases and healthy controls for 37 pathways, including lipid metabolic pathways (false discovery rate [FDR] <0.05). Particularly, the fecal microbiome of bronchiolitis cases had consistently higher abundances of gene function related to the sphingolipid metabolic pathways compared to that of controls (FDR <0.05). These pathways were more abundant in infants with Bacteroides-dominant microbiome profile compared to the others (FDR <0.001). On the basis of the predicted metagenome in this case–control study, we found significant differences in the functional potential of fecal microbiome between infants with bronchiolitis and healthy controls. Although causal inferences remain premature, our data suggest a potential link between the bacteria-derived metabolites, modulations of host immune response, and development of bronchiolitis.

Introduction

Bronchiolitis is a common acute respiratory infection and the leading cause of hospitalizations in US infants [1, 2]. Although bronchiolitis has been considered virus-induced inflammation of small airways [3], recent studies demonstrate that the pathobiology involves complex interrelations among respiratory viruses, host immune response, and human microbiome [4,5,6,7,8,9,10]. Emerging evidence also indicates the existence of “gut-lung axis” in which the gut microbiome conditions immunologic responses in the lungs to environmental challenges (e.g., viral infection) [11]. Indeed, we have previously demonstrated, in a case–control study of infants hospitalized for bronchiolitis and healthy controls [11], that the taxonomy profiles of the fecal microbiome were associated with the likelihood of bronchiolitis—e.g., infants with the Bacteroides-dominant profile were more likely to have bronchiolitis. Although previous studies suggest that the gut microbiome-derived metabolites (e.g., sphingolipids) may play an important role in the host immune development [12, 13], the functional profiles of fecal microbiome in infants were not examined in the earlier study. To address this knowledge gap, we determined the predicted function of fecal microbiome in infants with bronchiolitis and healthy infants.

Main text

Methods

This study was a secondary analysis of the data from a case–control study of infants hospitalized for bronchiolitis and healthy controls. The study design, setting, participants, and methods of data collection have been reported previously [11]. In brief, as a part of a multicenter prospective cohort study, called the 35th Multicenter Airway Research Collaboration (MARC-35) [4,5,6,7, 9], we enrolled 40 infants (aged <12 months) hospitalized for an attending physician diagnosis of bronchiolitis from November 2013 through April 2014. Bronchiolitis was diagnosed according to the American Academy of Pediatrics guidelines [14]. Exclusion criteria were a transfer to a participating hospital >48 h after the original hospitalization, delayed consent (>24 h after hospitalization), gestational age ≤32 weeks, and known comorbidities (cardiopulmonary disease, immunodeficiency, immunosuppression). In addition, during the same period, we also enrolled 115 healthy infants as the controls (age-matched within 1.5 months of cases) [11, 15,16,17]. We excluded infants with current fever, respiratory illness, or gastrointestinal illness, antibiotic treatment in the preceding 7 days, gestational age ≤32 weeks, or known comorbidities. Taken together, a total of 155 infants were eligible for the current analysis. From these infants, by using a standardized protocol [11, 15, 17], investigators conducted a structured interview and medical record review, and collected fecal specimens at the time of hospitalization (cases) or at home before the clinic visit (controls). The fecal samples were immediately stored at −80 °C. The institutional review board at each of the participating hospitals approved the study. Written informed consent was obtained from the parent or guardian.

16S rRNA gene sequencing was performed based on the methods adapted from the NIH Human Microbiome Project. Briefly, bacterial genomic DNA was extracted using MO BIO PowerMag DNA Isolation Kit (Mo Bio Lab; Carlsbad, CA). The 16S rDNA V4 regions were amplified by PCR and sequenced in the MiSeq platform (Illumina; SanDiego, CA) using 2 × 250 bp paired-end protocol. Sequencing read pairs were demultiplexed based on the unique molecular barcodes, and reads were merged using USEARCH v7.0.1090, allowing no mismatches and a minimum overlap of 50 bases. We trimmed the merged reads at the first base with a Q5 quality score. We calculated the expected error after taking into account all Q scores across all the bases of a read and the probability of an error occurring. We also applied a quality filter to the resulting merged reads, discarded the reads containing >0.05 expected errors. We constructed rarefaction curves of bacterial operational taxonomic units (OTUs) using sequence data for each sample to ensure coverage of the bacterial diversity present (Fig. 1). 16S rRNA gene sequences were clustered into OTUs at a similarity cutoff value of 97% using the UPARSE algorithm; OTUs were mapped to the SILVA Database to determine taxonomies. Abundances were recovered by mapping the demultiplexed reads to the UPARSE OTUs.

Fig. 1

Rarefaction curves for bacterial operational taxonomic units of the fecal microbiome. The horizontal axis indicates sequence depth while the vertical axis indicates the number of bacterial operational taxonomic units (OTUs). All 155 fecal specimens had sufficient depth to obtain high degree of sequence coverage (rarefaction cutoff, 1470 reads/specimen)

To infer the functional capabilities of the fecal microbiome based on the OTU (taxonomy) data, we used a bioinformatic approach, Tax4Fun [18]. This approach links the 16S rRNA gene sequences with the functional annotation of sequenced bacterial genomes by identifying a nearest neighbor based on a minimal 16S rRNA gene sequence similarity. Next, the predicted metagenomes were categorized by function at the Kyoto Encyclopedia of Genes and Genomes (KEGG) ortholog and pathway levels [19]. We tested for significant differences in the functional category abundances between cases and controls using Welch’s unequal variances t test. Resulting P values were adjusted for multiple hypothesis testing by converting to false discovery rate q values using the Benjamini–Hochberg procedure, with q values of <0.05 considered statistically significant. To validate the findings, we performed random permutation testing with 1000 permutations for each of the pathways of interest, which corresponds to the situation when the abundance of pathways is randomly assigned to cases and controls contained in dataset. Once the dataset was permuted, we tested for the differences in abundances between cases and controls. We repeated the randomization 1000 times and recorded the squared error of the models averaged for every repetition. Additionally, to further examine the differences in the pathways of interest, we constructed multivariable linear regression models adjusting for potential confounders (age, sex, race/ethnicity, maternal antibiotic use during pregnancy, history of prematurity, mode of delivery, feeding status, and lifetime history of antibiotic and corticosteroid use), based on a priori knowledge [5, 6, 9]. Furthermore, to determine the relationship between the abundance of bacteria genus and metabolic pathways of interest, we examined their correlations with the use of scatterplots fitting locally weighted scatterplot smoothed (LOWESS) curve as well as Spearman’s correlation. The analyses used R version 3.3 with phyloseq package [20] and STAMP version 2.1 [21].

Results

At the four participating hospitals, a total of 40 infants hospitalized for bronchiolitis (cases) and 115 age-matched healthy infants (controls) were enrolled (Table 1). Overall, the median age was 3 months (IQR, 2–5 months) and 55% were male. All 155 fecal specimens had sufficient depth to obtain high degree of sequence coverage (rarefaction cutoff, 1470 reads/specimen; Fig. 1). The fecal microbiome were dominated by four genera: Escherichia (22%), Bifidobacterium (19%), Enterobacter (15%), and Bacteroides (13%). The characteristics of the fecal microbiome differed between cases and controls (Table 2). For example, infants with bronchiolitis had a higher proportion of Bacteroides-dominant profile and lower proportion of Enterobacter/Veillonella-dominant profile, compared to healthy controls (P = 0.01).

Table 1 Patient characteristics of 40 cases (infants with bronchiolitis) and 115 controls (healthy infants) at enrollment

Table 2 Richness, alpha-diversity, and relative abundance of fecal microbiome in infants with bronchiolitis and healthy controls

Between the infants with bronchiolitis and healthy controls, we compared the functional potential of fecal microbiome inferred from the 16S rRNA gene sequencing data. Of 6402 KEGG orthologs (orthologous genes) surveyed, the abundances of 319 genes were significantly different (q < 0.05; Table 3). The functional differences involved genes with diverse metabolic functions—e.g., carbohydrate, amino acid, and lipid metabolism. To make the data presentation and interpretation more meaningful, the genes were further consolidated into 274 KEGG pathways. Among these, there were significant differences between bronchiolitis cases and healthy controls for 37 pathways, including lipid metabolic pathways (q < 0.05; Table 4; Fig. 2). Particularly, the fecal microbiome of bronchiolitis cases had consistently higher abundances of gene function related to the sphingolipid metabolic pathways compared to that of controls (all q < 0.05)—i.e., sphingolipid (ko00600) and glycosphingolipid (ko00603, ko00604) metabolic pathways (Fig. 3). For each of these 3 pathways, the permutation test was significant (all random permutation P < 0.05), supporting the validity of the observed between-group differences. In the multivariable models adjusting for 9 patient-level factors (age, sex, race/ethnicity, maternal antibiotic use during pregnancy, history of prematurity, mode of delivery, feeding status, and lifetime history of antibiotic and corticosteroid use), the difference in the abundances of 3 sphingolipid metabolic pathways remained significant (all P < 0.05). Additionally, these pathways were more abundant in infants with Bacteroides-dominant microbiome profile compared to the other microbiome profiles (all q < 0.001; Fig. 4). Likewise, there was a positive correlation between the abundance of Bacteroides genus and each of the 3 sphingolipid metabolic pathways (all P < 0.001; Fig. 5; Table 4).

Table 3 Predicted KEGG orthologs with significant differences in relative abundance between infants with bronchiolitis and healthy controls

Table 4 Predicted KEGG pathways with significant differences in relative abundance between infants with bronchiolitis and healthy controls

Fig. 2

Predicted KEGG pathways with significant differences in relative abundance between infants with bronchiolitis and healthy controls. Of 274 KEGG pathways surveyed, the relative abundance of 37 genes was significantly different (false discovery rate, q < 0.05) between infants with bronchiolitis and healthy controls

Fig. 3

Box-whisker plots of the three sphingolipid metabolic pathways that distinguish between the fecal microbiome of infants with bronchiolitis and that of healthy controls. The predicted metagenome of fecal microbiome in infants with bronchiolitis had a higher abundance of the a ko00600 (q = 0.03), b ko00603 (q = 0.03), and c ko00604 (q = 0.048) pathways compared to that in healthy controls. The horizontal line represents the median; the bottom and the top of the box represent the 25th and the 75th percentiles; whiskers represent 5 and 95% percentiles

Fig. 4

Box-whisker plots of the three sphingolipid metabolic pathways that distinguish four fecal microbiome profiles. The relative abundance of a ko00600, b ko00603, and c ko00604 pathways were consistently higher in infants with Bacteroides-dominant microbiome profile compared to the others (all q < 0.001). The four fecal microbiota profiles were derived using partitioning around medoids clustering method with Bray–Curtis distance. The optimal number of clusters was identified by the use of gap statistic. The horizontal line represents the median; the bottom and the top of the box represent the 25th and the 75th percentiles; whiskers represent 5 and 95% percentiles. BCP Bacteroides-dominant profile, BFP Bifidobacterium-dominant profile, ESP Escherichia-dominated profile, EVP Enterobacter/Veillonella-dominant profile

Fig. 5

Correlations between the abundance of Bacteroides and the three sphingolipid metabolic pathways. There was a positive correlation between the abundance of Bacteroides and each of the three sphingolipid metabolic pathways. a ko00600 (Spearman’s r = 0.77; P < 0.001), b ko00603 (Spearman’s r = 0.73; P < 0.001), and c ko00604 (Spearman’s r = 0.74; P < 0.001). The fitted line represents locally weighted scatterplot smoothed (lowess) curve

Discussion

By predicting the functional potential of the fecal microbiome from 40 infants with bronchiolitis and 115 healthy age-matched controls enrolled in a case–control study, we found significant differences in the abundance of genes related to multiple metabolic pathways. Of these, the gene function related to sphingolipid metabolic pathways was consistently more abundant in the fecal microbiome of bronchiolitis cases compared to that of healthy controls. The current study extends the previously identified association of Bacteroides-dominated fecal microbiome profile with higher likelihoods of bronchiolitis by demonstrating the functional potential of the gut microbiome in infants.

Sphingolipids are a class of complex lipids containing a backbone of sphingoid bases. These lipids have long been known as structural components of human cell membranes and as a component of surfactant, but have more recently emerged as signaling molecules that modulate the host immune response and contribute to the pathogenesis of respiratory diseases, such as bronchiolitis, pneumonia, and asthma [7, 22]. While sphingolipids production is ubiquitous in eukaryotes, it is also produced by several bacteria genera such as Bacteroides, Prevotella, and Porphyromonas [12]. Recently, experimental models reported that Bacteroides-derived sphingolipids (e.g., α-galactosylceramide) play an important role in host immunomodulation similar to lipopolysaccharide (LPS), another family of bacteria-derived glycolipid. For example, Wieland Brown et al. demonstrated that Bacteroides-derived α-galactosylceramide binds to CD1d and activates mouse and human invariant natural killer T (iNKT) cells both in vitro and in vivo [12]. In contrast, An et al., using neonatal mouse models, found that treatment with a different Bacteroides-derived glycosphingolipids (GSL-Bf717) reduces the number of colonic iNKT cells and subsequent colonic inflammation [13]. Although reverse causation—e.g., bronchiolitis per se or treatments for bronchiolitis result in perturbation of the fecal microbiome—is also possible, these prior studies, coupled with our findings, collectively suggest that Bacteroides-dominant microbiome in the gut, through their sphingolipid production, may contribute to inappropriate immune responses and bronchiolitis pathogenesis in infants. Our data should encourage future investigations into the mechanisms linking the individual gut microbiome-derived metabolites to the host immune response in the gut and respiratory tract (the gut-lung axis).

In sum, on the basis of the predicted metagenome in this case–control study, we found significant differences in the functional potential of fecal microbiome between infants with bronchiolitis and healthy controls. Particularly, the fecal microbiome in infants with bronchiolitis had consistently higher abundances of gene function related to the sphingolipid metabolic pathways. Although causal inferences remain premature, our data may suggest a potential link between the bacteria-derived metabolites, modulations of host immune response, and development of bronchiolitis. Our findings should facilitate further metagenomic, metatranscriptomic, and metabolomic (including Bacteroides-derived galactosylceramide [13]) investigations into the role of gut microbiome in the bronchiolitis pathogenesis. Our data also encourage researchers to integrate these “omics” approaches with mechanistic evaluations in experimental models in order to develop new preventive and therapeutic strategies (e.g., microbiome modification) for infants with bronchiolitis.

Limitations

Our study has several potential limitations. First, the location of fecal sample collection differed between cases and controls. However, in both populations, the fecal samples were refrigerated immediately after collection and the literature reported that refrigeration is associated with no significant alteration in fecal microbiota composition [23]. Second, the functional potential of fecal microbiome was inferred from the 16S rRNA gene sequencing data rather than measured by metabolomics or metatranscriptomics, or from metagenomic sequencing. However, a study has shown a strong correlation between the predicted metagenome and metagenome sequencing data (r > 0.85) in the NIH Human Microbiome Project samples (including fecal samples) [18]. Third, the concentration of metabolites was not measured in the fecal samples. This is an important area for examination in our future work. Lastly, the study design precluded us from examining the succession of fecal microbiome and its relation to the development of respiratory disease in early childhood. To address this question, the study populations are currently being followed longitudinally to age 6 years, with fecal sample collections at multiple time-points.