Skip to main content

Advertisement

Early-life gut microbiome modulation reduces the abundance of antibiotic-resistant bacteria

Article metrics

Abstract

Background

Antibiotic-resistant (AR) bacteria are a global threat. AR bacteria can be acquired in early life and have long-term sequelae. Limiting the spread of antibiotic resistance without triggering the development of additional resistance mechanisms is of immense clinical value. Here, we show how the infant gut microbiome can be modified, resulting in a significant reduction of AR genes (ARGs) and the potentially pathogenic bacteria that harbor them.

Methods

The gut microbiome was characterized using shotgun metagenomics of fecal samples from two groups of healthy, term breastfed infants. One group was fed B. infantis EVC001 in addition to receiving lactation support (n = 29, EVC001-fed), while the other received lactation support alone (n = 31, controls). Coliforms were isolated from fecal samples and genome sequenced, as well as tested for minimal inhibitory concentrations against clinically relevant antibiotics.

Results

Infants fed B. infantis EVC001 exhibited a change to the gut microbiome, resulting in a 90% lower level of ARGs compared to controls. ARGs that differed significantly between groups were predicted to confer resistance to beta lactams, fluoroquinolones, or multiple drug classes, the majority of which belonged to Escherichia, Clostridium, and Staphylococcus. Minimal inhibitory concentration assays confirmed the resistance phenotypes among isolates with these genes. Notably, we found extended-spectrum beta lactamases among healthy, vaginally delivered breastfed infants who had never been exposed to antibiotics.

Conclusions

Colonization of the gut of breastfed infants by a single strain of B. longum subsp. infantis had a profound impact on the fecal metagenome, including a reduction in ARGs. This highlights the importance of developing novel approaches to limit the spread of these genes among clinically relevant bacteria. Future studies are needed to determine whether colonization with B. infantis EVC001 decreases the incidence of AR infections in breastfed infants.

Trial registration

This clinical trial was registered at ClinicalTrials.gov, NCT02457338.

Background

Each year, more than two million people in the United States develop antibiotic-resistant infections, and at least 23,000 die as a result [1]. Furthermore, it is estimated that over 70% of the bacteria responsible for healthcare associated infections are resistant to at least one of the antibiotics used worldwide as a first-line therapy [2, 3]. Consequently, drug-resistant bacterial infections require longer treatment periods, which cost the US health care system an estimated $5 billion annually [4].

Antibiotic resistance is a natural phenomenon with ancient origins and a product of evolution, occurring in the environment long before the Anthropocene [5]; however, it has accelerated since the commercial development of antibiotics in the twentieth century. In particular, the increased administration of antibiotics to humans and animals has sparked a tremendous selective pressure favoring bacteria that harbor genetic elements that confer antibiotic resistance [6]. Many of these elements are transposable and increase in prevalence as a direct result of human activity mediated by exposure to antibiotics, disinfectants, and heavy metals [7, 8].

The rapid escalation of next-generation sequencing technologies, particularly metagenomics, has produced a remarkable amount of data profiling the ‘resistome’ of environmental niches, with the human gut microbiome itself emerging as a major reservoir for ARGs among commensal organisms [9, 10]. Collectively, the human gut microbiome harbors more ARGs compared to other environments [6, 10,11,12]. Furthermore, commensal bacteria likely play a key role in the evolution and dissemination of ARGs, even if they are not the intended target of antibiotic therapies [13, 14], and recent evidence demonstrates that maternal gut microbes harboring ARGs are transferred to newborn infants during or shortly after birth [15, 16]. This represents a substantial risk to human health, as ARGs can be transferred from commensal bacteria obtained at birth to pathogens [17] and both can easily be disseminated between individuals. Moreover, resistance genes can persist for years in the human gut, even without long-term antibiotic exposure. For instance, ARGs have been identified up to four years after antibiotic treatment has ceased, indicating that antibiotic resistance can persist for an indeterminate amount of time in the absence of selective pressure [18].

Historically, Bifidobacterium are speculated to have been more abundant in healthy breastfed infants in the US and Europe than they are today, given their differential abundance among US and European infants relative to sub-Saharan Africa and southeast Asia, where traditional birth and infant feeding practices are predominant [19,20,21,22]. Nonetheless, when bifidobacteria are absent, other taxa including Proteobacteriaceae, Bacteroidaceae, Staphylococcaceae, Clostridiaceae predominate in the breastfed infant gut [23,24,25]. Although these taxa have been associated with potentially negative long-term health outcomes to varying degrees [26], they are also the primary reservoirs of clinically relevant ARGs [10]. We recently demonstrated that extensive and durable changes occur in the breastfed infant gut microbiome resulting from the colonization of Bifidobacterium longum subsp. infantis (B. infantis) EVC001 [25], while others have demonstrated that breastfed infants with higher levels of bifidobacteria have a reduced abundance and lower frequency of genes associated with antibiotic resistance [24, 27].

Here, we used shotgun metagenomics to characterize the effect of an intervention with B. infantis EVC001 on the abundance of ARGs in breastfed infants. We found that colonization by B. infantis EVC001 resulted in a significant reduction in ARGs and the bacteria that harbor them.

Materials and methods

Sample collection

The details of the clinical trial design have been previously reported [25, 28]. Briefly, mother-infant dyads were recruited in the Davis and Sacramento metropolitan region of Northern California (USA), and informed consent was obtained from the mothers prior to enrollment. Only one subject, in the control group, received antibiotics directly prior to sample collection and sequencing. Otherwise, 29% of control infants and 31% of EVC001-fed infants were born by cesarean section, 32% of control infants and 45% of EVC001-fed infants were born to mothers given antibiotics for labor and 16% of control infants and 31% of EVC001-fed infants were born to mothers who tested as Group B Streptococcus positive. None of these clinical variables were significantly different between the groups (Additional file 1: Table S1). Infants in the EVC001-fed group were fed 1.8 × 1010 CFU of B. infantis EVC001 (ATCC SD-7035) for 21 consecutive days from day 7 to day 27 postnatal. B. infantis EVC001 was delivered as 156 mg of live bacteria (> 1.8 × 1010 CFU) diluted in 469 mg of lactose as an excipient and packaged in single use sachets. Mothers were trained by lactation consultants to mix the contents of the sachet in 5 mL of expressed breast milk and feed this to the infant each day. The probiotic was stored at − 20 °C by the families during the study and stability at − 20 °C was confirmed by a plate count.

DNA isolation and shotgun metagenomics sequencing

DNA was previously extracted from approximately 100 mg of frozen stools collected from infants on day 21 postnatal [25]. Briefly, the DNA was subjected to bead beating prior to column purification using a Zymo Fecal DNA Miniprep kit, according to the manufacturer’s instructions. Samples with remaining raw material and sufficient extracted DNA meeting input quality standards, as defined by the manufacturer, were selected for sequencing. Sixty samples meeting these criteria were identified. Metagenomic shotgun sequencing was performed at the California Institute for Quantitative Biosciences (QB3) (University of California, Berkeley) on an Illumina HiSeq 4000 platform using a paired-end sequencing approach with a targeted read length of 150 bp and an insert size of 150 bp.

Quality filtering and removal of human sequences

Demultiplexed fastq sequences were quality filtered, including adaptor trimming using Trimmomatic v0.36 [29] with default parameters. Quality-filtered sequences were screened to remove human sequences using GenCoF v1.0 [30] against a non-redundant version of the Genome Reference Consortium Human Build 38, patch release 7 (GRCh38_p7; www.ncbi.nlm.nih.gov). Human sequence-filtered raw reads were deposited in the Sequence Read Archive (SRA; https://www.ncbi.nlm.nih.gov/sra) under the reference number, PRJNA390646.

Taxonomic and strain profiling

Taxonomic profiling of the metagenomic samples was performed using MetaPhlAn2 [31], which uses a library of clade-specific markers to provide panmicrobial (bacterial, archaeal, viral, and eukaryotic) profiling (http://huttenhower.sph.harvard.edu/metaphlan2). Strain characterization was performed using PanPhlan [32]. PanPhlan is used in combination with MetaPhlAn2 to characterize strain-level variants in marker genes for a selected organism. For PanPhlan analysis, the pangenomes from Bifidobacterium longum (https://bitbucket.org/CibioCM/panphlan) were used as a reference. Both MetaPhlAn2 and PanPhlan were used with their default settings as described in the updated global profiling of the Human Microbiome Project (2017) [33, 34].

Absolute quantification of Enterobacteriaceae by qRT-PC real-time PCR

Quantification of the total Enterobacteriaceae was performed by quantitative real-time PCR using the group-specific 16S rRNA gene primers En-lsu3F (TGCCGTAACTTCGGGAGAAGGCA) and En-lsu3 (TCAAGGCTCAATGTTCAGTGTC) [35]. Each reaction contained 10 μL of 2× Power SYBR Green PCR master mix (Applied Biosystems), 0.4 μm of each primer, and 1 μL of template DNA. Thermal cycling was performed on a QuantStudio 3 Real-Time PCR System and consisted of a denaturation step of 10 min at 95 °C followed by 40 cycles of 15 s at 95 °C and 1 min at 60 °C. Standard curves for absolute quantification were generated using genomic DNA extracted from a pure culture of E. coli isolated from the fecal sample of an infant (E. coli 7005–1) at concentrations ranging from 103 to 107 CFU/mL.

Antibiotic resistance gene analysis

Human-filtered reads were first subjected to the Humann2 [36] pipeline (http://huttenhower.sph.harvard.edu/humann2) and subsequently screened for ARGs using the DIAMOND v0.9.10 program [37] to conduct a BLASTX-type search against the Comprehensive Antibiotic Resistance Database (CARD) [38] vAugust 2017. Outputs from DIAMOND were parsed and filtered using custom scripts to collect the top hits for each sample considering positive-hits reads with (i) E-value ≤10− 5 [39]; (ii) amino acid identity ≥90%; and (iii) alignment length ≥ 25 amino acids [40,41,42]. Tabular outputs were collapsed into a biom format [43] and taxonomic annotation of reads was assigned according to the CARD database and NCBI GeneBank using custom scripts. The final taxonomic affiliation analysis of the ARGs was performed with a MEtaGenome Analyzer (MEGAN v6.9) [44], which was used to compute and interactively explore the taxonomic content of the data set using the lowest common ancestor (LCA) method on the NCBI Taxonomy phylogenetic Tree of Life (LCA parameters: mini-score 35, top percentage 10%). The portion of sequence types or subtypes identified as ARGs in the total metagenome were defined as the relative abundance, using cross-sample normalization in ppm (one read in one million mapped reads) [42, 45].

PCR amplification and cloning of antibiotic-resistant genes

Seven differentially abundant ARGs identified to be significantly higher in control infants compared to EVC001-colonized infants (P ≤ 0.001; Bonferroni) were selected for detection in samples by PCR and functional screening. Illumina reads annotated within the significant ARGs were retrieved and independently assembled in SPADES v3.11 [46]. The assembled sequences were then cross-checked via BLASTX and BLASTP against the NCBI non-redundant protein database (nr, https://www.ncbi.nlm.nih.gov/protein/). To validate the presence of these ARGs in the fecal DNA, primers specific to the open reading frames of each of the seven differentially abundant ARGs were designed using Primer3 (http://bioinfo.ut.ee/primer3-0.4.0/) (Additional file 2: Table S5). Target genes were then amplified separately from 1 μL of fecal DNA from each of six samples with the highest abundance of ARGs according to the metagenomics analysis. Amplified fragments were separated by agarose-gel electrophoresis, purified using the NucleoSpin® kit (Macherey-Nagel) and sequenced using Sanger technology (DNA Sequencing Facility UC Davis). The percent identity of the obtained nucleotide sequences to the corresponding open reading frame of the assembled ARGs was calculated using global pairwise alignment (Needelman-Wunsch). Sequences were also searched against non-redundant protein and nucleotide databases of the NCBI using BLASTX and BLASTN, respectively, and against protein domain databases (Pfam) and COG databases using RPS-BLAST.

Minimal inhibitory concentrations

Minimal inhibitory concentrations (MICs) were determined according to the Clinical and Laboratory Standards Institute guidelines for microdilution susceptibility testing [47, 48]. Strains grown in LB broth overnight were adjusted to 1 × 106 CFU/mL and inoculated into Mueller-Hinton Broth containing one of six different clinically-relevant antibiotics (ampicillin, tetracycline, cefotaxime, cefazolin, cefepime, and gentamicin) ranging from 0.5 to 512 μg/mL in 96-well polystyrene microtiter plates. The microtiter plates were incubated for 24 h at 37 °C. The optical density (OD) of each well was measured at 600 nm using an automated microtiter plate reader (BIO-TEK, Synergy HT). The MIC corresponded to the lowest antibiotic concentration at which no growth was detected. All tests were performed in triplicate.

Whole genome sequencing and assembly of bacterial isolates

Approximately 100 mg of fecal sample from day 21 postnatal (subjects 7005, 7084, 7122, and 7174) were serially diluted onto EMB agar and incubated overnight at 37 °C. Three colonies from each subject that were either dark in color and/or had a green metallic sheen were selected for subsequent analysis. Selected isolates were grown overnight in LB broth overnight at 37 °C. For each strain, a 1 mL aliquot was centrifuged at 10,000×g for 5 min and the supernatant was removed. The cell pellet was transferred into DNA/RNA Shield Microbe Lysis tubes (Zymo Research, Irvine CA) and high-molecular weight genomic DNA was extracted using a Quick-DNA Fecal/Soil Microbe Miniprep Kit (Zymo Research, Irvine, CA). DNA was extracted following the manufacturer’s protocol with a mechanical lysis in a FastPrep96 (MP Biomedicals, Santa Ana, CA) for 15 s at 1,800 rpm. gDNA integrity was assessed by gel electrophoresis using a high-molecular weight 1 Kb Extension ladder (Invitrogen, Carlsbad, CA). The presence of a gDNA band at 40 kb and no shearing revealed intact gDNA. The gDNA was quantified using a Quant-iT™ dsDNA Assay Kit, high sensitivity (Invitrogen). gDNA purity was assessed using the Take3 microwell UV-Vis system (BioTek, Winooski, VT). Individually barcoded libraries were prepared for each isolate using 400 ng of high-molecular weight gDNA with an Oxford Nanopore 1D Rapid Barcoding Kit (SQK-RBK004) (Oxford Nanopore Technologies, Oxford UK) according to the manufacturer’s protocol. Barcoded samples were pooled and a 1× HighPrep PCR bead clean-up (MagBio, Gaithersburg, MD) of the fragmented and barcoded libraries prior to Rapid adapter ligation was included at the recommendation of Oxford Nanopore. The final 12-plexed pool was loaded onto an R9.4 flow cell and run for 15 h. A secondary run was performed using the same protocol for the seven isolates whose initial coverage was below 6×. Reads were basecalled in real time using MinKnow (ONT, Oxford UK). Data for both runs were combined for subsequent processing. Basecalled reads were demultiplexed and adapters were trimmed using Porechop (version 0.2.3, https://github.com/rrwick/Porechop). Reads were assembled with Canu v1.5 [49] with default parameters. Assembled genomes were converted into local blast databases and the CARD database protein sequences were used as a query against the assembled genomes using TBLASTN with a min E-value set at 0.001. Genome assemblies were deposited into the NCBI Gene Bank (https://www.ncbi.nlm.nih.gov/genbank/) under accession number PRJNA472982.

Statistical analysis

A Mann-Whitney test was used for statistical comparisons between the two groups. Significantly different ARGs between the EVC001-fed infants and controls were estimated using a Kruskal-Wallis one-way analysis of variance, coupled with a Bonferroni correction for cross-sample normalization. A Fisher’s exact test was used to establish a significant presence/absence of gene families in a strain-level analysis. Rarefaction curves were computed to estimate the distribution of the identified ARGs across samples. A nonparametric two-sample t-test was used to compare rarefaction curves using Monte Carlo permutations (n = 999). A Bray-Curtis dissimilarity matrix was constructed to estimate the global resistome differences among the samples and visualized via a Principal Coordinate Analysis (PCoA). A Permutational Multivariate Analysis of Variance Using Distance Matrices (adonis) was used to assess global resistome differences between treatments and the effect-size (R2) of colonization by EVC001 on the resistome. P-values for the PCoA panel was computed using F-tests based on sequential sums of squares from permutations of the raw data. P-values throughout the manuscript are represented by asterisks (*, P < 0.05; **, P < 0.01; ***, P < 0.001; ****, P < 0.0001).

Results

B. infantis EVC001 reshapes the gut microbiome in breastfed infants

Using shotgun metagenome sequencing, we characterized the taxonomic and antibiotic resistance profiles within the gut microbiome of 60 healthy, term, breastfed infants in Northern California (USA) at 21 days postnatal. Details of the study design and subject characteristics have been previously reported [25, 28] and a summary of the demographic data of the subjects analyzed here is presented in Additional file 1: Table S1. After quality filtering, Illumina sequencing led to a total of 1.6 billion paired-end (PE) reads, of which approximately 3.6% were discarded as human genome sequence contaminants, resulting in an average of 27 million PE reads per sample (Table 1). High-quality, human-filtered reads were subjected to taxonomic profiling following the updated analysis pipeline of the Human Microbiome Project [34] (September, 2017; see Materials and Methods).

Table 1 Overview of recovered metagenomics sequencing results from infants fed EVC001 and controls

A total of 202 bacterial species belonging to 76 genera, 43 families, 21 orders, 13 classes, and 7 phyla were identified across the samples (Additional file 3: Table S2). There were differences in the taxonomic distribution between the infants fed EVC001 and those who were not. Among the infants fed EVC001 (n = 29), 10 bacterial genera comprised 99% of the community, with the Bifidobacterium genus representing 90% of the total relative abundance of any identified genus out of a total of 55 identified genera (P < 0.0001, Kruskal-Wallis test) (Fig. 1a).

Fig. 1
figure1

Taxonomic classification of metagenomic reads for EVC001-fed infants and controls. a Relative abundance of the top bacterial genera identified between the two groups of infants. b Relative abundance of bacterial species belonging to the Bifidobacterium genus identified among groups. c Hierarchical clustering based on a strain-level analysis of Bifidobacterium longum subspecies. Gene family profiles of a subgroup of reference genomes were selected from a global (n = 38) strain analysis. Each column represents the presence or absence of genes in a sample or a reference genome with respect to the total pangenome. All samples from EVC001-fed infants clustered together with B. longum subsp. infantis ATCC 15697 (B. infantis), whereas the samples from infants in the control group clustered separately with other B. longum subspecies (e.g., B. suis, B. longum DJ01A, and B. longum NCC2705). Functional analysis of gene families confirmed that the EVC001 samples were dominated by B. infantis due to the presence of unique genes (e.g., Blon_2348 in B. infantis), while genes present only in B. longum subsp. longum (e.g., araD; araA), were abundant in the communities from control infants. P-values were computed for each gene via Fisher’s exact test according to group

In the control group (n = 31), 68 genera were identified, of which Bifidobacterium comprised an average of 38% of the total microbiome, whereas there was a greater proportion of other genera than in the EVC001-fed group, particularly Clostridium (P = 0.01, Kruskal-Wallis test; Fig. 1a). Within the Bifidobacterium genus, eight species were identified, of which Bifidobacterium longum was the most abundant, representing 86% of the total identified bacterial species within the EVC001-colonized infants and 19% within the controls (P < 0.0001, Kruskal-Wallis test; Fig. 1b). Other detected bifidobacteria included B. breve and B. bifidum, which accounted for 9.4 and 7%, respectively, in the microbiome of the control infants and considerably less (1.4 and 0.4%, respectively) in the EVC001-colonized group.

To discriminate the B. longum species at the subspecies level and determine the abundance of B. longum subsp. infantis and B. longum subsp. longum to specifically relate changes in the microbiome composition to colonization with B. infantis, we performed a strain-level analysis within the B. longum species using the pangenome gene-families database provided by PanPhlan. This database includes genes from 38 strains of B. longum subspecies (e.g., B. longum subsp. longum, B. longum subsp. infantis, and B. longum subsp. suis). PanPhlan recovered an average of 98.8% of all genes present in Bifidobacterium longum subsp. infantis ATCC 15697 [50] from every sample in the EVC001-fed group, representing 2,449 pangenome gene families. In contrast, 19 infants in the control group lacked detectable reads that mapped to B. longum subspecies genes in their metagenomes. The remaining control samples (n = 12) reported 43% coverage of B. infantis genes, and Bifidobacterium longum subsp. longum NCC2705 displayed the highest gene coverage (79%) across 1,708 pangenome gene families.

Samples and four representative reference genomes were hierarchically clustered based on pair-wise similarities between strains calculated via the Jaccard distance between gene family profiles (Fig. 1c). The resulting heatmap showed that Bifidobacterium longum subsp. infantis was substantially more abundant than other Bifidobacterium longum subspecies in the EVC001 group. Gene loci unique to the B. infantis reference genome and samples from B. infantis EVC001-fed infants revealed key genes that were more abundant, including human milk oligosaccharide (HMO) utilization clusters [50, 51]. These genes were completely absent among the 29 of 31 infants who were not fed B. infantis EVC001, indicating that B. infantis was exceptionally rare (only 6% of infants) unless the infants were fed B. infantis EVC001. Genes unique to B. longum subsp. longum that enable characteristic arabinose consumption (araD and araA) were significantly enriched among infants harboring B. longum subsp. longum and rare among infants colonized by B. infantis EVC001. Together, these findings suggest that B. infantis was the dominant B. longum subspecies among infants fed B. longum subsp. infantis EVC001, and that B. infantis was exceptionally rare among these infants unless the infants were fed the strain during the clinical trial. The high level of B. infantis persistence previously reported in infants fed EVC001 in conjunction with the pangenome analysis here suggests stable colonization of the infants by this strain. This is in agreement with other models examining the efficient and durable colonization of host-associated gut microbes in their coevolved host [52, 53] and our prior longitudinal analysis of fecal samples by 16S rDNA sequencing [25]. In terms of microbial diversity, we have previously reported that there was no difference in community richness (alpha diversity) measured via Shannon index, but there are differences in the relative abundance of taxa and community stability (beta diversity) as assessed by UniFrac distance and Jaccard index [25].

Colonization by EVC001 is associated with a reduced ARG burden

We identified a total of 599,631 infant gut microbial genes from the shotgun sequencing data, of which 80,925 were unique to 29 infants fed B. infantis EVC001, and 313,683 microbial genes were unique to samples from 31 infants not fed B. infantis EVC001. Both groups shared a total of 205,023 microbial genes. The metagenomes were screened for ARGs using a BLASTX type search against the curated Comprehensive Antibiotic Resistance Database (CARD). After quality filtering the BLAST results (see Materials and Methods), a total of 652 ARGs were identified (Additional file 4: Table S3). The EVC001-fed group reported an average of 0.01% of ARGs among the total microbial genes (min = 0.001%; max = 0.18%; SEM = 0.006%), with 285 different ARGs (Fig. 2, a), of which 33 were found only in the EVC001 group at very low percentages (< 0.05%). Among the infants not fed B. infantis EVC001, these ARGs accounted, on average, for 0.09% of the total metagenomic reads (min = 0.004%; max = 0.24%; SEM = 0.01) with 612 of the different ARGs that were identified. Of these 612 ARGs, 360 uniquely belonged to this group. Thus, infants fed EVC001 had, on average, 90% fewer ARGs in their microbiome (P < 0.0001, Mann-Whitney test).

Fig. 2
figure2

Relative abundance of the total resistome profile in each metagenome sample. a Relative abundance of antibiotic resistance genes (ARGs) compared with the overall metagenome for each sample. Each point represents a sample resistome (control, n = 31; EVC001-fed, n = 29). Box plots denote the interquartile range (IQR), with horizontal lines representing the 25th percentile, median, and 75th percentiles. The whiskers represent the lowest and highest values within 1.5 times the IQR from the first and third quartiles, respectively. The asterisks on the top indicate significant P-values (Mann-Whitney test). b Relative abundance of ARGs according to their taxonomic identification. The shade of color represents genera belonging to the same bacterial class. The asterisks on the top indicate significant P-values (Kruskal-Wallis test)

To compare the microbial taxonomic affiliation of ARGs, the 652 ARGs identified among the optimal BLAST hits were assigned to different taxa according to the NCBI taxonomy guidelines coupled with the Lowest Common Ancestor (LCA) method in MEGAN [44]. A total of 41 bacterial genera were taxonomically assigned to the 652 ARGs, of which Escherichia, Staphylococcus, Bacteroides, and Clostridioides were associated with the majority of the ARGs (68.9, 5, 4, and 2.6%, respectively). Considering the taxonomic content within the resistome, the metagenomes from control infants had 17 bacterial genera with a relative abundance > 0.001%, with Escherichia-ARGs accounting for about 0.054% of the total metagenome (Fig. 2b). In the EVC001-fed group, only 12 bacterial genera had a relative abundance of associated ARGs > 0.001%. Escherichia was also the genus that carried the majority of ARGs but contributed significantly less to the overall metagenome (0.003%) compared with the controls (P = 0.001, Kruskal-Wallis test; Fig. 2b).

EVC001 significantly decreases the abundance of key antibiotic-resistant genes

Among the ARGs uniquely identified in the samples from infants in the control group, three were present in a relative abundance greater than 0.1% and were associated with the Clostridium genus. Specifically, tetA(P) and tetB(P), which are ARGs found on the same operon, were identified. tetA(P) is an inner membrane tetracycline efflux protein and tetB(P) is a ribosomal protection protein, both of which confer resistance to tetracycline [54, 55]. mprF was found only in the samples from infants in the control group, and acts by negatively charging phosphatidylglycerol on the bacterial membrane and confers resistance to antibiotic cationic peptides that disrupt the cell membrane, including host-produced defensins [56].

After cross-sample normalization, 38 ARGs were found to significantly differ between the two groups (P < 0.01, Kruskal-Wallis test). All 38 ARGs were significantly lower in the EVC001-fed group compared to the controls (Fig. 3a). Notably, none of the ARGs were significantly higher in the samples from the EVC001-fed group compared to the control group (P > 0.05, Kruskal-Wallis test). Genes enriched in the metagenome of infants in the control group were found to confer resistance to beta-lactams, fluoroquinolones, and macrolides, and 12 genes conferred resistance to multiple drug classes.

Fig. 3
figure3

Comparison of the most significant antibiotic resistance gene types. a Relative abundance of the most significantly different antibiotic resistance genes (ARGs) identified among EVC001-fed infants and controls (P < 0.02; Bonferroni). Percentages are relative to the overall metagenome. These ARGs confer resistance to different drug classes, including beta-lactams, fluoroquinolones, and macrolides. The ARGs are grouped by gene name, followed by CARD identification entry (ARO). The colored bars represent respective drug class to which the ARG is known to confer resistance to. b Heatmap showing a hierarchical cluster analysis of the total ARGs identified (n = 652). Two groups were identified, one characterized by a lower-ARG carriage, containing most of the samples from infants fed EVC001 and one characterized by a higher overall carriage, containing most samples from infants in the control group. Genes clustered based on similar biological mechanisms implicated in drug resistance (see Results). P-values on the bar were computed using a Kruskal-Wallis test normalized with a Bonferroni correction

Hierarchical clustering of the samples and total ARGs (n = 652) using the complete-linkage method generated two main clusters of samples (Fig. 3b). The majority of the samples from the EVC001-fed infants were clustered together within the lower-ARGs abundance panel. Row clustering by ARGs resulted in two groups. The most abundant genes clustered together and were annotated as being directly related to mechanisms of antimicrobial resistance. In particular, the proteins encoded by mdtB and mdtC form a heteromultimer complex that results in a multidrug transporter [57]. AcrD is an aminoglycoside efflux pump and its expression is regulated by baeR and cpxAR, which were also identified among the significant ARGs and best characterized in E. coli [58]. We also identified AcrB and TolC, which form the multidrug efflux complex, AcrA-AcrB-TolC, that confers multidrug resistance [59]. Moreover, RosA and RosB were significantly more abundant among infants not fed EVC001; these genes form an efflux pump/potassium antiporter system (RosAB) [60]. Three genes belonging to the multidrug efflux system, EmrA-EmrB-TolC, first identified in E. coli [61], were also significantly more abundant in the samples from control infants. In this complex, EmrB is an electrochemical-gradient powered transporter, whereas EmrA is the linker, and TolC is the outer membrane channel [62]. The complex confers resistance to fluoroquinolones, nalidixic acid, and thiolactomycin.

The majority (76%) of the significantly different ARGs were taxonomically assigned to bacteria belonging to the Enterobacteriaceae family (e.g., Escherichia coli) and its abundance was proportional to the presence of ARGs (R = 0.58; P < 0.00001, Pearson) (Fig. 3b). The absolute abundance (determined by qPCR) of Enterobacteriaceae was significantly lower (P < 0.0001) in EVC001-colonized infants (Fig. 4). Other ARGs reported multiple taxonomic assignments within the Proteobacteria phylum. According to NCBI’s taxonomic assignment and the CARD database, they could originate from any one of multiple, closely related species. These included: the efflux pump acrD; the MdtG protein, which appears to be a member of the major facilitator superfamily of transporters, that confers resistance to fosfomycin and deoxycholate [63]; BaeR, a response-regulator conferring multidrug resistance [64]; and marA, a global activator protein overexpressed in the presence of different antibiotic classes [65]. A global statistical analysis of ARGs by treatment group is reported in Additional file 5: Table S6.

Fig. 4
figure4

Quantification of Enterobacteriaceae family by group-specific qPCR. The data are represented as Log10 CFU per μg of genomic DNA extracted from stool samples. Data in boxplots show the median, first, and third quartiles (P < 0.0001, Mann-Whitney Test)

Colonization by EVC001 decreases the total abundance and composition of ARGs

To compare the overall impact of EVC001 colonization on the diversity of ARGs, the alpha-diversity (e.g., number of unique ARGs) within each sample was compared using rarefaction curves. Notably, the diversity of the ARGs was independent from the number of sequences per sample (Fig. 5a). Overall, the EVC001-fed infants had half as many unique ARGs as the control infants (P = 0.001, t-test). The global resistome differences among samples and the effect-size of EVC001 colonization on the overall diversity of the two study groups were assessed. A Bray-Curtis dissimilarity matrix transposed into a principal coordinate analysis (PCoA) showed that samples from the EVC001-colonized group were clustered closely together compared to the control, which had a wider distribution (Fig. 5b; P = 0.001, F-test). This indicates that samples from the EVC001-colonized infants had a less abundant and less diverse resistome compared with the control group samples. Colonization with EVC001 contributed to a greater than a 30% reduction in global AR diversity in the infant gut microbiome than in the gut of the controls (R2 = 0.31, P = 0.001, adonis). Finally, there was no statistically significant difference in the abundance of ARGs detected in the control group, whether babies were born via C-section or vaginally (P = 0.30; Mann–Whitney) or whether infants were exposed to intrapartum antibiotics (P = 0.5; Mann–Whitney). Conversely, infants fed B. infantis EVC001 had significantly lower ARG abundance, independently from delivery mode, antibiotic exposure or any other clinical variable.

Fig. 5
figure5

Diversity analyses of infant resistomes according to B. infantis EVC001 colonization. a Rarefaction curves showing the number of unique antibiotic resistance genes (ARGs) identified in relation to the increasing number of sequences. Both EVC001 and the control group presented similar curve trends, suggesting that the sequencing depth is not associated with the diversity of antibiotic resistance. P-values were computed with a nonparametric two-sample t-test using Monte Carlo permutations (n = 999). b Global resistome profiles computed via a principal coordinate analysis (PCoA) based on a Bray-Curtis dissimilarity matrix. The EVC001-fed samples clustered closely, indicating that they shared a similar resistome compared to the controls, which had a more dispersed distribution. The effect of B. infantis EVC001 colonization by itself accounted for 31% of the total explained variation (adonis). The P-value was computed using F-tests based on the sequential sums of squares from permutations of the raw data

In vitro validation of in silico predicted ARGs

To validate the ARG presence in fecal samples, PCR primer pairs were designed to target seven of the most abundant ARGs in the resistome of the control infants. Amplicons were obtained in at least half of the analyzed fecal samples, with the exception of the primer pair targeting the mfd gene, which did not amplify from any sample. A nucleotide sequence analysis of the generated amplicons confirmed the sequence identity, with a nucleotide identity of > 70% to the open reading frame (ORF) of the predicted target gene. The nucleotide sequence analysis revealed a high homology (85–99%) with the genomic regions annotated to encode the expected functions in gut bacteria, and the predicted amino acid sequences contained highly conserved structural and functional domains in the corresponding encoded proteins (Additional file 6: Table S4).

To confirm the presence of full length, functional ARGs and the relationship of these ARGs to individual resistance phenotypes at the strain level, bacteria isolated from EMB agar were obtained from the fecal samples of four representative control infants. EMB was used to select for lactose-fermenting coliforms (e.g., E. coli) which had been identified by shotgun metagenome sequencing to harbor the most ARGs across individuals in both groups. Whole-genome sequencing of 12 isolates was performed on a MinION sequencer and assembly led to an average coverage of 18× (min 5.4; max 40) (Table 2). Taxonomic identification was confirmed via BLASTN against the NCBI nucleotide database (https://github.com/rrwick/Porechop), revealing three isolates classified as Raoultella planticola and the remaining nine as Escherichia coli. CARD protein sequence collection was used as a query against the 12 assembled isolates via TBLASTN. The presence of the 40 most differentially abundant ARGs identified via shotgun metagenomics (P < 0.02; Bonferroni) was confirmed in all or some of the 12 genomes (average identity > 93%), except for rosB (ARO:3003049), and rob (ARO:3004108). The latter genes were absent from the E. coli and R. planticola genomes but were predicted in other taxa not isolated here (Table 2).

Table 2 Summary of assembly statistics and ARGs assignments for 12 bacterial isolates obtained from four control subjects

The minimum inhibitory concentration (MIC) to ampicillin, cefepime, cefotaxime, cefazolin, tetracycline, and gentamicin was determined for these isolates. With the exception of three isolates obtained from the same infant (7174), all of the isolates exhibited ampicillin resistance. Among the multidrug-resistance isolates, resistance to ampicillin, cefazolin, and tetracycline was the most common. No resistance to gentamicin was detected (Table 3).

Table 3 Minimum inhibitory concentrations (MIC) in μg/mL of antibiotics

Discussion

B. infantis is a well-known organism with a long history of evolutionary adaptation to the breastfed infant gut [50, 51]. We previously demonstrated that B. infantis EVC001 was able to remodel the infant gut microbiome and lower the abundance of common gut taxa, particularly genera belonging to the Proteobacteria and Firmicutes phyla (e.g., Escherichia and Clostridium) [23, 25]. In the present study, shotgun metagenomics enabled us to confirm the presence of B. infantis as the dominant B. longum subspecies among EVC001-fed infants using strain-level analytical tools [32]. Subspecies within Bifidobacterium longum present substantial differences in their genetic architecture that have a profound impact on carbohydrate utilization phenotypes [50, 51]. In particular, human milk oligosaccharides (HMOs), which account for the third largest component of human milk, are not digested by the infant and in the absence of B. infantis, are lost in the stool [25, 66,67,68]. Several functions of HMOs have been proposed, including immune signaling molecules, as well as a role as prebiotics [69]. Thus, changes in the HMO concentrations in milk might trigger shifts in the microbiome composition with an impact on health outcomes [70].

To confirm that the dominant B. longum subspecies in the EVC001-fed infant gut microbiome was indeed B. longum subsp. infantis, we collected a pangenome database, composed of 12,000 genes from 38 available subspecies of B. longum. A hierarchical cluster analysis based on gene presence or absence, showed that EVC001-fed infants clustered with the B. infantis reference genome (ATCC 15697), whereas the control group clustered separately with other B. longum subspecies. At the gene level, the HMO cluster-1, which is the largest cluster of HMO-utilization genes and is absent among the B. longum subsp. longum strains [50], was present only in the EVC001-fed group but was absent from the control group. In contrast, the L-arabinose isomerase (araA) and L-ribulose-5-phosphate 4-epimerase (araD) that confer the genetic capacity to ferment xylose and arabinose are absent in B. infantis but present in B. longum subspecies longum [50] and were only found in samples from the control group, which were colonized with Bifidobacterium longum subsp. longum. Notably, 60% of the control group did not have a detectable amount of any B. longum subspecies genes in their metagenomes, suggesting that these bacteria are missing from the microbiome of these infants, regardless of whether they were delivered vaginally or by cesarean section. This could reflect the generational loss of key gut bacterial symbionts in the general population as a result of increasing rates of cesarean section, antibiotic use, and formula feeding over the past century [21, 71, 72].

Given that the infant gut microbiota typically has low colonization resistance during the first two years of life, there is ample opportunity for the infant to acquire antibiotic-resistant populations of commensal organisms, including from the hospital environment [71, 73]. Once established, these resistant populations could contribute to the spread of ARGs. Antibiotic-resistant infections place a substantial health and economic burden on the US health care system and population [74], and there is a limited means by which the spread of ARGs, or the organisms that harbor them, can be restricted [75]. Recent studies have shown that the microorganisms transmitted at birth from the mother and the surrounding environment influence the initial colonization of the infant gut microbiome; however, few studies have investigated the ARG burdens in the gut microbiome of infants colonized by nosocomial bacteria. Recently, Taft et al. showed that an abundance of bifidobacteria is associated with a significant reduction in ARGs in Bangladesh [24]. In the present study, we demonstrated that even in the absence of an antibiotic selective pressure, the healthy term breastfed infant gut of a Northern California population harbors a variety of genes encoding resistance to several clinically relevant antibiotic classes. Furthermore, we demonstrate the possibility of an intervention via modulation of the gut microbiome to decrease the abundance of ARGs. Our results suggest that the effects of targeted probiotic intervention can recapitulate the decrease in ARG abundance observed in Bifidobacterium-dominated infants from developing nations [24]. In fact, feeding EVC001 resulted in a significant decrease of 90% in the ARG burden (0.01%) relative to the infants who were not fed EVC001 (0.1%). The computed relative abundance of ARGs in our study was lower compared to what has been found in adults (approximately 0.4–1%) [10, 76], but was similar to other studies in infants [77, 78].

After accounting for cross-sample normalization, we could not identify any ARG increase in the samples from infants fed B. infantis EVC001; however, 38 ARGs were significantly increased in the control samples relative to the infants fed EVC001. These genes are known to confer resistance to a variety of drug classes, including fluoroquinolones, beta-lactams, macrolides, and tetracycline antibiotics. A total of 12 significant ARGs were classified as multidrug-resistant genes and reported the widest differences between the two groups. Previous targeted or functional studies on healthy infants have found similar ARGs, including those conferring tetracycline resistance (tet), which we only identified in the control infants [79]; and several multidrug ARGs [80], including CRP, emrD, and CpxR (ARO:3004055), which were the most significantly different genes between the two groups (P < 0.0001). Similar ARG abundance and identities reported in our study have also been found in the gut microbiome of preterm infants [81].

These results suggest that the colonization of B. infantis EVC001 in infants not only remodeled the gut microbiome composition in breastfed infants, but also may provide an immediate translational impact through a reduction in ARG carriage. Diversity analyses indicated that there were significant reductions in the number of unique ARGs in the EVC001-fed group compared to the controls, and this trend was independent of the sequencing depth. This aspect might translate to clinical relevance, since infants whose infections originate within the gut microbiome could be less likely to exhibit resistance to a wider spectrum of drug classes. Further studies are needed to determine the durability of these ARG profiles following cessation of the probiotic.

Although the use of shotgun metagenomics coupled with clinical databases (e.g., CARD) offer a direct means of identifying ARGs on a global, high-throughput scale, short-read sequencing might represent a challenge when linked to the resistance of specific genes that may originate from multiple strains of the same species [82]. Moreover, the computational annotation of ARGs can be difficult to translate into actual drug resistance, with some evidence showing a disconnect between the predicted and actual phenotypes [83]. Additionally, alternative mechanisms independent from genetic variations, could confer drug resistance, including biofilm formation [84,85,86]. Similarly, the absence of ARGs may not necessarily translate into non-resistant phenotypes. Nevertheless, by sequencing the ORFs of the most abundant ARGs in the fecal DNA of control infants, we confirmed that these genes were present in the analyzed fecal samples. More importantly, homology searches against well-curated sequence databases (e.g., pfam and COG) confirmed that the obtained gene sequences contained functionally conserved domains. Genome sequencing and phenotypic confirmation of antibiotic-resistant phenotypes confirmed the association between the presence of many of these genes and antibiotic resistance using clinically defined breakpoints, highlighting the fact that these genes are present in predicted organisms and demonstrably confer clinically relevant levels of resistance.

Finally, it is important to note that this study is limited to ARGs found in public databases and thus contains only known genes involved in antibiotic resistance. Therefore, the identified ARGs may be underrepresented, as our analysis excluded chromosomal mutations and novel, resistance genes not found within the CARD database. However, the majority of ARGs are associated with Proteobacteria and Firmicutes, whose populations were significantly reduced in infants fed B. infantis EVC001 compared to those who were not. Natural, nontransferable antibiotic resistance phenotypes have long been established among bifidobacteria, and recent systematic analyses have differentiated between natural and acquired Bifidobacterium resistance phenotypes [87]; however, B. infantis lacks any known virulent genes or plasmids [50], and there are currently no studies describing evidence of horizontal gene transfer events.

Conclusions

Colonization of newborn infants with antibiotic-resistant bacteria represents a major risk to the continued dissemination of ARGs within human populations. Colonization of the gut microbiome with B. infantis EVC001 early in life reduces the abundance of these ARGs as well the taxa that harbor and potentially transmit them between individuals. In conjunction with appropriate drug stewardship practices by the medical community, this approach of modulation of the gut microbiome could help reduce the burden and diversity of ARGs and limit their transmission across individuals. Additional research is required to determine whether these changes have clinically relevant benefits, such as changing the prevalence of antibiotic-resistant nosocomial infections.

Availability of data and materials

The sequencing libraries generated in this study are publicly deposited at the NCBI sequence read archive (SRA) under accession number PRJNA390646. Whole-genome assemblies are deposited in NCBI under accession number PRJNA472982.

References

  1. 1.

    Zhang L, Kinkelaar D, Huang Y, Li Y, Li X, Wang HH. Acquired antibiotic resistance: are we born with it? Appl Environ Microbiol. 2011;77(20):7134–41.

  2. 2.

    Klevens RM, Edwards JR, Richards CL Jr, Horan TC, Gaynes RP, Pollock DA, et al. Estimating health care-associated infections and deaths in US hospitals, 2002. Public Health Rep. 2007;122(2):160–6.

  3. 3.

    Carmeli Y. Strategies for managing today’s infections. Clin Microbiol Infect. 2008;14(s3):22–31.

  4. 4.

    Dodds DR. Antibiotic resistance: A current epilogue. Biochem Pharmacol. 2017;134:139–46.

  5. 5.

    Hernández J, Stedt J, Bonnedahl J, Molin Y, Drobni M, Calisto-Ulloa N, et al. Human-associated extended-spectrum β-lactamase in the Antarctic. Appl Environ Microbiol. 2012;78(6):2056–8.

  6. 6.

    Forsberg KJ, Reyes A, Wang B, Selleck EM, Sommer MO, Dantas G. The shared antibiotic resistome of soil bacteria and human pathogens science 2012;337(6098):1107–1111.

  7. 7.

    Toussaint A, Chandler M. Prokaryote genome fluidity: toward a system approach of the mobilome. Bacterial Mol Netw. 2012:57–80.

  8. 8.

    Gillings MR. Evolutionary consequences of antibiotic use for the resistome, mobilome and microbial pangenome. Front Microbiol. 2013;4.

  9. 9.

    Salyers AA, Gupta A, Wang Y. Human intestinal bacteria as reservoirs for antibiotic resistance genes. Trends Microbiol. 2004;12(9):412–6.

  10. 10.

    Hu Y, Yang X, Qin J, Lu N, Cheng G, Wu N, et al. Metagenome-wide analysis of antibiotic resistance genes in a large cohort of human gut microbiota. Nat Commun. 2013;4:2151.

  11. 11.

    Penders J, Stobberingh EE, Savelkoul PH, Wolffs PF. The human microbiome as a reservoir of antimicrobial resistance. Front Microbiol. 2013;4.

  12. 12.

    Sommer MO, Dantas G, Church GM. Functional characterization of the antibiotic resistance reservoir in the human microflora science. 2009;325(5944):1128–31.

  13. 13.

    Andremont A. Commensal flora may play key role in spreading antibiotic resistance. ASM news. 2003;69(12):601–7.

  14. 14.

    Marshall BM, Ochieng DJ, Levy SB. Commensals: underappreciated reservoir of antibiotic resistance. Microbe. 2009;4(5):231–8.

  15. 15.

    Alicea-Serrano AM, Contreras M, Magris M, Hidalgo G, Dominguez-Bello MG. Tetracycline resistance genes acquired at birth. Arch Microbiol. 2013;195(6):447–51.

  16. 16.

    Yassour M, Jason E, Hogstrom LJ, Arthur TD, Tripathi S, Siljander H, et al. Strain-level analysis of mother-to-child bacterial transmission during the first few months of life. Cell Host Microbe. 2018;24(1):146–54 e4.

  17. 17.

    Jernberg C, Löfmark S, Edlund C, Jansson JK. Long-term impacts of antibiotic exposure on the human intestinal microbiota. Microbiology. 2010;156(11):3216–23.

  18. 18.

    Jakobsson HE, Jernberg C, Andersson AF, Sjölund-Karlsson M, Jansson JK, Engstrand L. Short-term antibiotic treatment has differing long-term impacts on the human throat and gut microbiome. PLoS One. 2010;5(3):e9836.

  19. 19.

    Young SL, Simon MA, Baird MA, Tannock GW, Bibiloni R, Spencely K, et al. Bifidobacterial species differentially affect expression of cell surface markers and cytokines of dendritic cells harvested from cord blood. Clin Diagn Lab Immunol. 2004;11(4):686–90.

  20. 20.

    Huda MN, Lewis Z, Kalanetra KM, Rashid M, Ahmad SM, Raqib R, et al. Stool microbiota and vaccine responses of infants. Pediatrics. 2014;134(2):e362–e72.

  21. 21.

    Tannock GW, Lee PS, Wong KH, Lawley B. Why Don't all infants have Bifidobacteria in their stool? Front Microbiol. 2016;7.

  22. 22.

    Henrick BM, Hutton AA, Palumbo MC, Casaburi G, Mitchell RD, Underwood MA, et al. Elevated Fecal pH Indicates a Profound Change in the Breastfed Infant Gut Microbiome Due to Reduction of Bifidobacterium over the Past Century. mSphere. 2018;3(2):e00041–18.

  23. 23.

    Casaburi G, Frese SA. Colonization of breastfed infants by Bifidobacterium longum subsp. infantis EVC001 reduces virulence gene abundance. Human Microbiome J. 2018.

  24. 24.

    Taft DH, Liu J, Maldonado-Gomez MX, Akre S, Huda MN, Ahmad S, et al. Bifidobacterial dominance of the gut in early life and acquisition of antimicrobial resistance. mSphere. 2018;3(5):e00441–18.

  25. 25.

    Frese SA, Hutton AA, Contreras LN, Shaw CA, Palumbo MC, Casaburi G, et al. Persistence of Supplemented Bifidobacterium longum subsp. infantis EVC001 in Breastfed Infants. mSphere. 2017;2(6).

  26. 26.

    Vatanen T, Kostic AD, d’Hennezel E, Siljander H, Franzosa EA, Yassour M, et al. Variation in microbiome LPS immunogenicity contributes to autoimmunity in humans. Cell. 2016;165(4):842–53.

  27. 27.

    Pärnänen K, Karkman A, Hultman J, Lyra C, Bengtsson-Palme J, Larsson DJ, et al. Maternal gut and breast milk microbiota affect infant gut antibiotic resistome and mobile genetic elements. Nat Commun. 2018;9(1):3891.

  28. 28.

    Smilowitz JT, Moya J, Breck MA, Cook C, Fineberg A, Angkustsiri K, et al. Safety and tolerability of Bifidobacterium longum subspecies infantis EVC001 supplementation in healthy term breastfed infants: a phase I clinical trial. BMC Pediatr. 2017;17(1):133.

  29. 29.

    Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30(15):2114–20.

  30. 30.

    Czajkowski MD, Vance DP, Frese SA, Casaburi G. GenCoF: A graphical user interface to rapidly remove human genome contaminants from metagenomic datasets. Bioinformatics. 2018;35(13):2318–9.

  31. 31.

    Truong DT, Franzosa EA, Tickle TL, Scholz M, Weingart G, Pasolli E, et al. MetaPhlAn2 for enhanced metagenomic taxonomic profiling. Nat Methods. 2015;12(10):902.

  32. 32.

    Scholz M, Ward DV, Pasolli E, Tolio T, Zolfo M, Asnicar F, et al. Strain-level microbial epidemiology and population genomics from shotgun metagenomics. Nat Methods. 2016;13(5):435–8.

  33. 33.

    Consortium HMP. Structure, function and diversity of the healthy human microbiome. Nature. 2012;486(7402):207.

  34. 34.

    Lloyd-Price J, Mahurkar A, Rahnavard G, Crabtree J, Orvis J, Hall AB, et al. Strains, functions and dynamics in the expanded human microbiome project. Nature. 2017;550(7674):61.

  35. 35.

    Matsuda K, Tsuji H, Asahara T, Kado Y, Nomoto K. Sensitive quantitative detection of commensal bacteria by rRNA-targeted reverse transcription-PCR. Appl Environ Microbiol. 2007;73(1):32–9.

  36. 36.

    Abubucker S, Segata N, Goll J, Schubert AM, Izard J, Cantarel BL, et al. Metabolic reconstruction for metagenomic data and its application to the human microbiome. PLoS Comput Biol. 2012;8(6):e1002358.

  37. 37.

    Buchfink B, Xie C, Huson DH. Fast and sensitive protein alignment using DIAMOND. Nat Methods. 2015;12(1):59–60.

  38. 38.

    McArthur AG, Waglechner N, Nizam F, Yan A, Azad MA, Baylay AJ, et al. The comprehensive antibiotic resistance database. Antimicrob Agents Chemother. 2013;57(7):3348–57.

  39. 39.

    Minot S, Sinha R, Chen J, Li H, Keilbaugh SA, Wu GD, et al. The human gut virome: inter-individual variation and dynamic response to diet. Genome Res. 2011;21(10):1616–25.

  40. 40.

    Zhang T, Zhang X-X, Ye L. Plasmid metagenome reveals high levels of antibiotic resistance genes and mobile genetic elements in activated sludge. PLoS One. 2011;6(10):e26041.

  41. 41.

    Kristiansson E, Fick J, Janzon A, Grabic R, Rutgersson C, Weijdegård B, et al. Pyrosequencing of antibiotic-contaminated river sediments reveals high levels of resistance and gene transfer elements. PLoS One. 2011;6(2):e17038.

  42. 42.

    Yang Y, Li B, Ju F, Zhang T. Exploring variation of antibiotic resistance genes in activated sludge over a four-year period through a metagenomic approach. Environ Sci Technol. 2013;47(18):10197–205.

  43. 43.

    McDonald D, Clemente JC, Kuczynski J, Rideout JR, Stombaugh J, Wendel D, et al. The biological observation matrix (BIOM) format or: how I learned to stop worrying and love the ome-ome. GigaScience. 2012;1(1):7.

  44. 44.

    Huson DH, Auch AF, Qi J, Schuster SC. MEGAN analysis of metagenomic data. Genome Res. 2007;17(3):377–86.

  45. 45.

    McCall CA, Bent E, Jørgensen TS, Dunfield KE, Habash MB. Metagenomic comparison of antibiotic resistance genes associated with liquid and dewatered biosolids. J Environ Qual. 2016;45(2):463–70.

  46. 46.

    Bankevich A, Nurk S, Antipov D, Gurevich AA, Dvorkin M, Kulikov AS, et al. SPAdes: a new genome assembly algorithm and its applications to single-cell sequencing. J Comput Biol. 2012;19(5):455–77.

  47. 47.

    Wikler MA. Methods for dilution antimicrobial susceptibility tests for bacteria that grow aerobically: approved standard. CLSI (NCCLS). 2006;26:M7-A.

  48. 48.

    CLSI. Performance Standards for Antimicrobial Susceptibility Testing: CLSI supplement. M100S. 26th ed. Wayne, PA: Clinical and Laboratory Standards Institute; 2016.

  49. 49.

    Koren S, Walenz BP, Berlin K, Miller JR, Bergman NH, Phillippy AM. Canu: scalable and accurate long-read assembly via adaptive k-mer weighting and repeat separation. Genome Res. 2017:gr. 215087.116.

  50. 50.

    Sela D, Chapman J, Adeuya A, Kim J, Chen F, Whitehead T, et al. The genome sequence of Bifidobacterium longum subsp. infantis reveals adaptations for milk utilization within the infant microbiome. Proc Natl Acad Sci. 2008;105(48):18964–9.

  51. 51.

    Underwood MA, German JB, Lebrilla CB, Mills DA. Bifidobacterium longum subspecies infantis: champion colonizer of the infant gut. Pediatr Res. 2015;77:229.

  52. 52.

    Maldonado-Gómez MX, Martínez I, Bottacini F, O’Callaghan A, Ventura M, van Sinderen D, et al. Stable engraftment of Bifidobacterium longum AH1206 in the human gut depends on individualized features of the resident microbiome. Cell Host Microbe. 2016;20(4):515–26.

  53. 53.

    Frese SA, Benson AK, Tannock GW, Loach DM, Kim J, Zhang M, et al. The evolution of host specialization in the vertebrate gut symbiont lactobacillus reuteri. PLoS Genet. 2011;7(2):e1001314.

  54. 54.

    Roberts MC. Update on acquired tetracycline resistance genes. FEMS Microbiol Lett. 2005;245(2):195–203.

  55. 55.

    Bannam TL, Johanesen PA, Salvado CL, Pidot SJ, Farrow KA, Rood JI. The Clostridium perfringens TetA (P) efflux protein contains a functional variant of the motif a region found in major facilitator superfamily transport proteins. Microbiology. 2004;150(1):127–34.

  56. 56.

    Peschel A, Jack RW, Otto M, Collins LV, Staubitz P, Nicholson G, et al. Staphylococcus aureus resistance to human defensins and evasion of neutrophil killing via the novel virulence factor MprF is based on modification of membrane lipids with l-lysine. J Exp Med. 2001;193(9):1067–76.

  57. 57.

    Kim H-S, Nikaido H. Different functions of MdtB and MdtC subunits in the heterotrimeric efflux transporter MdtB2C complex of Escherichia coli. Biochemistry. 2012;51(20):4188–97.

  58. 58.

    Rosenberg EY, Ma D, Nikaido H. AcrD of Escherichia coli is an aminoglycoside efflux pump. J Bacteriol. 2000;182(6):1754–6.

  59. 59.

    Eicher T, Brandstätter L, Pos KM. Structural and functional aspects of the multidrug efflux pump AcrB. Biol Chem. 2009;390(8):693–9.

  60. 60.

    Bengoechea JA, Zhang L, Toivanen P, Skurnik M. Regulatory network of lipopolysaccharide O-antigen biosynthesis in Yersinia enterocolitica includes cell envelope-dependent signals. Mol Microbiol. 2002;44(4):1045–62.

  61. 61.

    Lomovskaya O, Lewis K, Matin A. EmrR is a negative regulator of the Escherichia coli multidrug resistance pump EmrAB. J Bacteriol. 1995;177(9):2328–34.

  62. 62.

    Andersen JL, He G-X, Kakarla P, KC R, Kumar S, Lakra WS, et al. Multidrug efflux pumps from Enterobacteriaceae, Vibrio cholerae and Staphylococcus aureus bacterial food pathogens. Int J Environ Res Public Health. 2015;12(2):1487–547.

  63. 63.

    Fàbrega A, Martin RG, Rosner JL, Tavio MM, Vila J. Constitutive SoxS expression in a fluoroquinolone-resistant strain with a truncated SoxR protein and identification of a new member of the marA-soxS-rob regulon, mdtG. Antimicrob Agents Chemother. 2010;54(3):1218–25.

  64. 64.

    Nagakubo S, Nishino K, Hirata T, Yamaguchi A. The putative response regulator BaeR stimulates multidrug resistance of Escherichia coli via a novel multidrug exporter system, MdtABC. J Bacteriol. 2002;184(15):4161–7.

  65. 65.

    Randall L, Woodward MJ. The multiple antibiotic resistance (mar) locus and its significance. Res Vet Sci. 2002;72(2):87–93.

  66. 66.

    Newburg DS. Oligosaccharides in human milk and bacterial colonization. J Pediatr Gastroenterol Nutr. 2000;30:S8–S17.

  67. 67.

    McSweeney PL, Fox PF. Advanced dairy chemistry: volume 3: lactose, water, salts and minor constituents: springer; 2009.

  68. 68.

    Albrecht S, Schols HA, van den Heuvel EG, Voragen AG, Gruppen H. Occurrence of oligosaccharides in feces of breast-fed babies in their first six months of life and the corresponding breast milk. Carbohydr Res. 2011;346(16):2540–50.

  69. 69.

    Bode L. Recent advances on structure, metabolism, and function of human milk oligosaccharides. J Nutr. 2006;136(8):2127–30.

  70. 70.

    Davis JC, Lewis ZT, Krishnan S, Bernstein RM, Moore SE, Prentice AM, et al. Growth and morbidity of Gambian infants are influenced by maternal Milk oligosaccharides and infant gut microbiota. Sci Rep. 2017;7:40466.

  71. 71.

    Palmer C, Bik EM, DiGiulio DB, Relman DA, Brown PO. Development of the human infant intestinal microbiota. PLoS Biol. 2007;5(7):e177.

  72. 72.

    Tannock GW. Commentary: remembrance of microbes past. Int J Epidemiol. 2005;34(1):13–5.

  73. 73.

    Brooks B, Firek BA, Miller CS, Sharon I, Thomas BC, Baker R, et al. Microbes in the neonatal intensive care unit resemble those found in the gut of premature infants. Microbiome. 2014;2(1):1.

  74. 74.

    Golkar Z, Bagasra O, Pace DG. Bacteriophage therapy: a potential solution for the antibiotic resistance crisis. J Infect Dev Ctries. 2014;8(02):129–36.

  75. 75.

    Lee C-R, Cho IH, Jeong BC, Lee SH. Strategies to minimize antibiotic resistance. Int J Environ Res Public Health. 2013;10(9):4274–305.

  76. 76.

    Forslund K, Sunagawa S, Kultima JR, Mende DR, Arumugam M, Typas A, et al. Country-specific antibiotic use practices impact the human gut resistome. Genome Res. 2013;23(7):1163–9.

  77. 77.

    Versluis D, D’Andrea MM, Garcia JR, Leimena MM, Hugenholtz F, Zhang J, et al. Mining microbial metatranscriptomes for expression of antibiotic resistance genes under natural conditions. Sci Rep. 2015;5:11981.

  78. 78.

    Gibson MK, Wang B, Ahmadi S, Burnham C-AD, Tarr PI, Warner BB, et al. Developmental dynamics of the preterm infant gut microbiota and antibiotic resistome. Nat Microbiol. 2016;1:16024.

  79. 79.

    Francino M. Antibiotics and the human gut microbiome: dysbioses and accumulation of resistances. Front Microbiol. 2016;6:1543.

  80. 80.

    Fouhy F, Ogilvie LA, Jones BV, Ross RP, Ryan AC, Dempsey EM, et al. Identification of aminoglycoside and β-lactam resistance genes from within an infant gut functional metagenomic library. PLoS One. 2014;9(9):e108016.

  81. 81.

    Rose G, Shaw AG, Sim K, Wooldridge DJ, Li M-S, Gharbia S, et al. Antibiotic resistance potential of the healthy preterm infant gut microbiome. PeerJ. 2017;5:e2928.

  82. 82.

    Willmann M, Peter S. Translational metagenomics and the human resistome: confronting the menace of the new millennium. J Mol Med (Berl). 2017;95(1):41–51.

  83. 83.

    Davis MA, Besser TE, Orfe LH, Baker KN, Lanier AS, Broschat SL, et al. Genotypic-phenotypic discrepancies between antibiotic resistance characteristics of Escherichia coli from calves in high and low antibiotic use management settings. Appl Environ Microbiol. 2011.

  84. 84.

    Corona F, Martinez JL. Phenotypic resistance to antibiotics. Antibiotics. 2013;2(2):237–55.

  85. 85.

    Drenkard E, Ausubel FM. Pseudomonas biofilm formation and antibiotic resistance are linked to phenotypic variation. Nature. 2002;416(6882):740–3.

  86. 86.

    Balaban NQ, Merrin J, Chait R, Kowalik L, Leibler S. Bacterial persistence as a phenotypic switch. Science. 2004;305(5690):1622–5.

  87. 87.

    Duranti S, Lugli GA, Mancabelli L, Turroni F, Milani C, Mangifesta M, et al. Prevalence of antibiotic resistance genes among human gut-derived bifidobacteria. Appl Environ Microbiol. 2017;83(3):e02894–16.

Download references

Acknowledgements

Sequencing in this study was conducted at the Vincent J. Coates Genomics Sequencing Laboratory at UC Berkeley, which is supported by the NIH S10 OD018174 Instrumentation Grant.

Funding

This study was funded by Evolve Biosystems, Inc.

Author information

GC, SAF, JTS, and MAU contributed to the study or experimental design. RMD, SAF, LC, and RM conducted the experiments. GC and DV analyzed the data. JTS organized the clinical trial design and its operation. All authors contributed to writing, editing, and approving the manuscript.

Correspondence to Jennifer T. Smilowitz or Mark A. Underwood.

Ethics declarations

Ethics approval and consent to participate

The study was conducted under the approval of the University of California, Davis Institutional Review Board and registered at ClinicalTrials.gov under study NCT02457338. All experiments were performed in accordance with relevant guidelines and regulations.

Consent for publication

Consent for the participation and publication of the study results was obtained from subjects or their legal guardians.

Competing interests

GC, RMD, SAF, LC, RM, and DV are employed by Evolve BioSystems, Inc.

Additional information

Publisher’s Note

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

Additional files

Additional file 1:

Table S1. Study population demographics. (DOC 50 kb)

Additional file 2:

Table S5. Primers sequences used to amplify seven differentially expressed ARGs from the fecal DNA (DOC 30 kb)

Additional file 3:

Table S2. Global taxonomic profile by sample. Relative abundance (%) of individual taxa identified across samples. Stratification refers to different taxonomic rank (k = kingdom; p = phylum; c = class; o = order; f = family; g = genus; s = species; t = strain). Within the same taxonomic rank (e.g., Family, “f__”) the total relative abundance sum up to 100%. (XLS 427 kb)

Additional file 4:

Table S3. Global ARG profile by sample. Number of quality filtered hits for every sample against the CARD database. First column reports the annotation in the following order: protein ID based on NCBI reference sequence database; antibiotic resistance gene ID according to the CARD database and protein name. (XLS 312 kb)

Additional file 5:

Table S6. Kruskal-Wallis test with post-hoc corrections for every ARG identified according to treatment status. (XLSX 51 kb)

Additional file 6:

Table S4. Nucleotide sequence analysis of ORFs. Global alignment results of selected open reading frames (ORFs) from Illumina reads assembly. BLASTX and BLASTN scores of selected ORFs identified among ARGs are reported. (XLSX 11 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Keywords

  • Antibiotic resistance gene
  • Metagenomics
  • Probiotics
  • Host-microbe interactions
  • Microbiology