|Home | About | Journals | Submit | Contact Us | Français|
The mammalian gut ecosystem has significant influence on host physiology1–4, but the mechanisms that sustain this complex environment in the face of different stresses remain obscure. Perturbations to this ecosystem, such as through antibiotic treatment or diet, are currently interpreted at the level of bacterial phylogeny5–7. Less is known about the contributions of the abundant population of phage to this ecological network. Here, we explore the phageome as a potential genetic reservoir for bacterial adaptation by sequencing murine fecal phage populations following antibiotic perturbation. We show that antibiotic treatment leads to the enrichment of phage-encoded genes that confer resistance via disparate mechanisms to the administered drug as well as genes that confer resistance to antibiotics unrelated to the administered drug, and we demonstrate experimentally that phage from treated mice afford aerobically cultured naïve microbiota increased resistance. Systems-wide analyses uncover post-treatment phage-encoded processes related to host colonization and growth adaptation, indicating that the phageome broadly enriches for functionally beneficial genes under stress-related conditions. We also show that antibiotic treatment expands the interactions between phage and bacterial species, leading to a more highly connected phage-bacterial network for gene exchange. Our work implicates the phageome in the emergence of multidrug resistance and indicates that the adaptive capacity of the phageome may represent a community-based mechanism for protecting the gut microflora, preserving its functional robustness during antibiotic stress.
Antibiotic treatment, an important and often necessary therapeutic intervention, can negatively impact the mammalian gut environment, potentially giving rise to immune2 and metabolic deficiencies8. Studies on the disruption of intestinal homeostasis have focused on the resulting alterations in microbial composition6,7. However, investigation of the gut ecosystem has evinced a myriad of resident phage9, and it remains unclear how perturbation of the gut environment affects these symbionts. Phage can contribute genes that are advantageous to their microbial hosts10,11, in turn promoting their own survival and propagation12. This gene flow suggests that phage may play an important role in the adaptation of the microbiome to stressful environments. In this work, we used a comparative metagenomic approach to explore the effects of antibiotic perturbation on functions encoded in the phageome as well as to examine how antibiotic treatment alters the phage-bacterial ecological network.
We treated groups of young adult mice (n = 5) orally with physiologically relevant concentrations of ciprofloxacin (a quinolone that inhibits DNA synthesis) or ampicillin (a β-lactam that inhibits cell-wall synthesis), each with a respective control. We obtained collective fecal samples from each group after eight weeks of treatment and purified phage as previously described9,13. DNA was extracted from phage and whole-genome amplified before performing shotgun 454 GS FLX+ pyrosequencing. We obtained a total of 440,792 quality reads, with a median read length of 477 nt (210 Mb total; Fig. S1). Evaluation of contamination by quantification of bacterial 16S ribosomal RNA genes indicated that contaminating bacterial sequences represented less than 0.1% of our data, which was subsequently accounted for in all statistical analyses (see Supplementary Discussion and Fig. S2).
Phage DNA sequences were compared to the non-redundant NCBI protein and environmental protein databases (BLASTX; E-value < 10−5). Approximately 70% of reads were not assigned to previously sequenced genes (Fig. S3), suggesting that the mouse phageome, like many other viral communities14,15, harbors uncharacterized genetic material. We used the most significant BLAST alignment of a sequence, when available, to determine its phylogenetic origin. Most of the identifiable phage in our mouse phageomes (Fig. S4a) were from the Caudovirales order, comprising the tailed phage families Siphoviridae, Podoviridae, and Myoviridae, many of which are known to have a temperate life cycle. Since phage genomes incorporate bacterial genes, we also identified bacterial taxa; we found that 97% of phage-encoded bacterial genes were attributed to the four phyla known to dominate the gut (the Firmicutes, Bacteroidetes, Proteobacteria, and Actinobacteria; Fig. S4b), consistent with the known hosts of the phage we detected.
We hypothesized that antibiotic treatment leads to increases in phage-encoded genes for drug resistance. To investigate this, we compared DNA sequences in the phageome to an assembled database of antibiotic resistance proteins (BLASTX; E-value < 10−3, see Methods). We found that reads annotated as antibiotic resistance genes were highly enriched in phage metagenomes from mice treated with ciprofloxacin or ampicillin compared to those from control mice (Z = 7.3 and Z = 7.0, respectively; Fig. S5; read annotations enumerated in Table S1). We catalogued the resistance reservoir by annotating phage-encoded genes based on the drug class to which they confer resistance (Fig. 1a). Our analysis revealed that resistance to the administered drug class was enriched in phage metagenomes from antibiotic-treated mice, such that resistance to DNA synthesis inhibitors was enriched in ciprofloxacin treatment (Z = 2.6), and resistance to cell-wall synthesis inhibitors was enriched in ampicillin treatment (Z = 5.0). Additionally, upon drug treatment, new resistance genes were found in the phageome. For example, phage from ciprofloxacin-treated mice carried genes encoding numerous quinolone efflux pumps (e.g., norM, mexD, mexF), and phage from ampicillin-treated mice carried genes encoding sensor and response regulators of cell-wall synthesis inhibitors (e.g., vanRS).
Of note, resistance to other, orthogonal drug classes was also overrepresented in the antibiotic-perturbed phageomes. Both treatments led to significant enrichment of resistance to antibiotics that target protein synthesis, and ciprofloxacin treatment also led to significant enrichment of resistance to cell-wall synthesis inhibitors (Fig. 1a). This cross-resistance was mediated by drug-specific inactivators (e.g., chloramphenicol acetyltransferases), as well as multidrug resistance exporters (e.g., mdtK; Table S1). Together, these findings implicate the phage metagenome as a potential source of multidrug resistance during antibiotic treatment of the host.
We aimed to understand the specific mechanisms represented in phage-encoded genes conferring resistance to the most significantly enriched classes of antibiotics: inhibitors of cell-wall synthesis, DNA synthesis, and protein synthesis. We categorized resistance genes according to primary resistance mechanisms, which include modification or protection of the drug target (target modification), enzymatic inactivation of the drug (drug inactivation), and transport of the drug out of the cell (efflux)16. Using this framework to classify phage-encoded resistance genes, our analysis revealed that antibiotic treatment led to disparate resistance mechanism profiles for each drug class (Fig. 1b). Analysis of resistance to cell-wall synthesis inhibitors showed that all types of resistance mechanisms were significantly enriched under both ciprofloxacin treatment and ampicillin treatment. In contrast, resistance to DNA synthesis and protein synthesis inhibitors occurred predominantly by efflux. Resistance to protein synthesis inhibitors occurred by target modification and drug inactivation mechanisms at low levels, and, in accord with its rarity, resistance to DNA synthesis inhibitors via drug inactivation was not detected. These data likely reflect resistance mechanisms that are both environmentally available in the gut ecosystem and impose sustainable in vivo fitness costs. As continued treatment with an antibiotic invariably leads to its own resistance, mechanisms that enable cross-resistance encoded by the phageome may be an important consideration when selecting subsequent therapeutics.
We next sought to demonstrate that phage from antibiotic-treated mice confer increased drug resistance to the host-associated bacterial community. We assessed the frequency of resistant isolates from aerobically cultured naïve microbiota that were infected ex vivo with phage from antibiotic-treated or control (untreated) mice. Our results show that this fraction of microbiota infected with phage from mice administered ciprofloxacin or ampicillin yielded two to three times more colonies resistant to the respective drug than the aerobically cultured fraction of microbiota infected with phage from control mice (Fig. 1c). These data indicate that phage from antibiotic-treated mice can contribute relevant functional advantages to their microbial hosts.
We next took a systems-level approach and classified other phage-encoded genes into functional pathways described by the Kyoto Encyclopedia of Genes and Genomes (KEGG) database. We depicted enriched functional changes in phage metagenomes following antibiotic treatment as a network diagram (Fig. 2; abundances shown in Fig. S6). Among the most significantly enriched pathways were functional properties related to the mode of action of the administered drug (Table S2). Phageomes from ampicillin-treated mice were enriched for the amino sugar and nucleotide metabolism pathway (part of the broader carbohydrate metabolism process; Z = 5.6) indicating overrepresentation of genes related to synthesis of cell wall constituents, and increases in these components have been shown to be requisite for drug resistance in clinical isolates17. Additionally, we found that phageomes from ciprofloxacin-treated mice were enriched for replication- and repair-related pathways, including base excision repair (Z = 6.1), nucleotide excision repair (Z = 7.4), and homologous recombination (Z = 11.2). Included in these pathways are members of the GO system for the repair of DNA oxidative lesions, which has been demonstrated to reduce cytotoxicity due to a range of antibiotic classes18. Also represented are members of the DNA damage-inducible SOS system, known to provide protection against antibiotic-mediated cell death and induce the development of resistance-conferring mutations19. Furthermore, hyper-recombination has been shown to promote multidrug resistance phenotypes20. These results show that under drug treatment, the phageome encodes diverse mechanisms for modulating antibiotic susceptibility.
We also observed that phage metagenomes from antibiotic-treated mice were enriched for microbial functions that contribute to host metabolism (Fig. 2 and Table S2). Phageomes from ciprofloxacin-treated mice were uniquely enriched for pathways relevant to the metabolism of cofactors and vitamins, including thiamine, an essential nutrient provided by the microbiome. Microbiota ferment polysaccharides indigestible by the host alone; metabolism of these sugars enables bacterial survival in and colonization of the gut environment21 and, as a beneficial consequence, provides energy to the host4. We found that polysaccharide degradation genes, specifically related to metabolism of starch, cellulose, lactose, and fructans (plant-derived fructose polymers), were enriched with antibiotic treatment (Fig. 3a). Genes coding for carbohydrate active enzymes (CAZymes), which enable bacteria to ferment a variety of dietary and host-sourced glycans, were also enriched under antibiotic treatment and represented by a range of glycoside hydrolase and glycosyltransferase families (Fig. S7 and Table S3). Because many gut microbes express only a specific array of carbohydrate-degrading enzymes22, bacteria which acquire these genes from the phage reservoir may gain additional foraging capacity and, consequently, a selective growth advantage. These results suggest that the phageome may be an adaptive repository for functions important for the host-commensal relationship that may otherwise be depleted by antibiotic perturbation. (See Supplementary Information for additional discussion.)
We next aimed to elucidate the phylogenetic basis of these phage-encoded bacterial functions. In Fig. 3b, we illustrate the taxonomic composition of all sequences of bacterial origin and sequences annotated with enriched functions following drug perturbation. Examining bacterial phylotypes that contribute antibiotic resistance genes to the phage metagenome, we found a comparatively high representation of the Clostridia class and a low representation of the Bacilli class. Notably, a large fraction of CAZyme-annotated sequences originated from the Bacteroidetes class, which comprises members found to have diverse capabilities for carbohydrate metabolism22. Investigation of thiamine metabolism reveals that the Bacilli and Verrucomicrobia classes constitute a large proportion of these annotations. Since the phageome reflects emergent properties of its environment15, phylogenetic analyses of phage-encoded elements may more broadly enable the identification of bacteria actively contributing to specific functions in the gut environment.
Our results show that the phageome harbors a diversity of potentially beneficial functional elements in the face of antibiotic perturbation. However, the extent to which this genetic reservoir is accessible to members of the microbiota remains unclear23. To investigate this, we sought to elucidate the phage-bacterial ecological network and how it changes under antibiotic treatment. We approximated the network of phage-microbe interactions with relationships identified through the reconstruction and analysis of individual viral genomes. De novo assembly was accomplished using stringent parameters (see Methods). Reconstructed viral genomes are composed of a mosaic of bacterial genes, and we utilized the phylogenetic origins of these sequences to determine putative phage-bacterial associations (Fig. 4). Our resulting network recapitulated known interactions, including the lysogenic relationships of foodborne pathogens, such as bacteriophage phi3626 infection of C. perfringens24 and Siphoviridae Listeria phage A500 infection of L. monocytogenes25. Importantly, antibiotic treatment leads to widespread restructuring of the phage-bacterial ecological network (Fig. 4). These data show that new links between phage and bacteria are formed with drug treatment, giving rise to significantly greater network connectivity (Fig. S8). This increased connectivity is reflective of phage broadly, as more bacterial species are associated with a given phage (Fig. S9). These results suggest that antibiotic treatment increases the frequency of phage integration and stimulates broad host range, which promotes a functional reservoir that is both genetically diverse and highly accessible by gut bacteria.
While the phageome is a highly connected network for gene exchange, the functional consequence of acquiring genetic material from this reservoir depends on the molecular context of the host bacterium. Acquisition of a single gene may enhance an existing function by gene dosage or enable a novel phenotype. As some proteins rely on additional machinery, subsequent horizontal gene transfer events may be required to produce a phenotype. Moreover, redundantly acquired genes may gain new functionality through paralogous evolution.
Here we demonstrate that antibiotic treatment enriches the phage metagenome for stress-specific and niche-specific functions, while mediating changes in the topology of the phage-bacterial ecological network to potentiate accessibility of these genetic elements. Functional resilience of the microbiome following environmental perturbation has been empirically documented and has engendered interest in the restorative forces that return the commensal flora to its pre-perturbed state6,7,26. Our results implicate phage encapsulation of adaptive signatures as a community-based mechanism for functional robustness in the gut environment during stress. Of note, antibiotic treatment can also prime the gut environment for pathogen invasion2, and our findings pose potential implications for the emergence of drug resistance and evasion strategies in pathogenic populations. Cohabitation of phage and bacteria in the gut ecosystem is likely governed by complex and dynamic interactions, particularly during stress perturbation. Additional work is needed to discern the selective pressures imposed on each member of this community and the resulting mechanisms that influence the encoding and progressive enrichment of functionally beneficial genes in the phageome. Phage-mediated gene flow may be an important phenotypic buffer for bacterial communities, and further investigation of the adaptive reservoir of the phageome and the dynamic nature of the phage-bacterial ecological network may prove critical to understanding the influence of the gut ecosystem on host physiology.
All experiments involving animals were pre-reviewed and approved by the Boston Children’s Hospital Institutional Animal Care and Use Committee (IACUC). Groups of separately housed 6-week-old female FVB/NJ mice (n = 5; Jackson Laboratory) were treated with ampicillin (142.5 mg/L) or ciprofloxacin (62.5 mg/L) in their drinking water. The corresponding dosage was 28.5 mg/kg/day ampicillin and 12.5 mg/kg/day ciprofloxacin, based on the average mouse weight of 20 g and an approximate intake of 4 ml per day. Control mouse groups were supplied with standard drinking water (ampicillin) or alkaline water (ciprofloxacin). Treatments were refreshed twice per week. Mice were housed in sterile conditions and received autoclaved chow during the course of this study. We harvested fresh collective fecal samples from each group to obtain ample material (3–4 g) for purification. Samples were stored at −80 °C prior to use.
Viral purification from fecal samples of mice after eight weeks was performed as previously described9,13. An aliquot of the viral preparation (1.5 g ml−1 layer from ultracentrifugation) was stained with SYBR gold and visualized with epifluorescence microscopy to verify the absence of bacterial contamination. Viral particles were concentrated and desalted using an Ultra-4 Centrifugal Filter Unit (Ultracel-30K MWCO; Millipore) to a volume of ~200 μl. Concentrated viral samples were treated with DNase (0.2 mg/ml) and samples were passed through a 0.22 μM filter to ensure no procedural contamination was introduced. Genomic DNA was extracted using the QiaAMP DNA mini kit (Qiagen) using the protocol for viral DNA detailed in the manual. Genomic DNA was amplified using the Illustra Genomiphi v2 kit (GE) according to manufacturer’s instructions. For each sample, we pooled amplified DNA from three separate reactions to minimize bias.
Viral DNA was submitted to the Tufts Genomic Core for library preparation and shotgun sequencing. Equimolar concentrations of multiplexed (Rapid Library MID) samples were pooled on a single plate and pyrosequenced using the 454 GS FLX+ platform. Resulting sequences were filtered by removing duplicates using the tool available at http://microbiomes.msu.edu/replicates/28 with the following parameters: sequence identity cutoff = 97%; length difference requirement = 0; number of beginning base pairs to check = 20.
To facilitate annotation of antibiotic resistance genes, we assembled a database consisting of the Antibiotic Resistance Genes Database27 and UniProt proteins annotated with the GO term “antibiotic breakdown” (GO: 0017001) to achieve a total of 12,687 protein sequences. The functional annotations of the proteins in this database are supported by either experimental validation or high-quality computational prediction. Viral DNA sequences were compared to this database using BLASTX, and sequences with an E-value < 10−3 were deemed significant. This cutoff was selected to maintain a consistent stringency, by accounting for database size, with functional annotation to the NCBI databases.
Phage DNA sequences were compared to the non-redundant NCBI protein (nr) and environmental protein (env_nr) databases (BLASTX; E-value < 10−5). Sequences were annotated with KEGG29 (v61.0) using custom Perl scripts interfaced with the KEGG API. The most significant BLAST hit that resulted in a KEGG annotation was used, and we included all annotations in the event that a hit had multiple annotations. KEGG ortholog annotations were compared to glycoside hydrolases and glycosyltransferases found on www.cazy.org30 to identify CAZyme-encoding genes.
To compare the functional annotations of two metagenomic datasets, A and B, we generated a distribution for each A and B reflecting the number of annotations from 10,000 trials of random sampling with replacement, sampled at the number of reads in the comparison dataset. To account for the effects of contaminating bacterial DNA in our comparative metagenomic analyses, we assumed that contamination would be uniformly distributed and therefore we randomly discarded bacterial reads in each sampling trial according to the amount of contamination detected by quantitative PCR. We compared the number of annotations identified in a given sample to the comparison distribution and determined a Z-score, calculated as (x−μ)/σ, where, for example, x is the raw number of annotations in sample A, μ is the mean number of annotations in the distribution for sample B, and σ is the standard deviation of B’s distribution. In essence, this results in two Z-scores, one comparing A to B and one comparing B to A, and consequently the minimum |Z| was used. According to the Central Limit Theorem, random sampling with replacement results in a normal distribution, so Z-scores > 1.65 (P < 0.05) were considered enriched. In KEGG pathway analysis, Bonferroni was used to correct for multiple hypotheses, where the P-value cutoff 0.05/n (n was the number of total third-level pathways identified in our phage metagenomes) was converted to a Z-score.
Microbiota were isolated from fecal samples of naïve mice as previously described31, except PBS + 0.1% cysteine supernatants were plated on four separate LB agar plates. Colonies were grown aerobically for 24 hours before plates were scraped with 2 ml LB, amassed, and cultured at 37 °C and 300 RPM for two hours. This mixture was stored as 150 μl aliquots in 15% glycerol at −80 °C. For each experiment, an entire aliquot was used as the inoculum to minimize growth biases. Phage were isolated as described above. Due to limited sample availability, we harvested phage from mice that had been treated for five weeks for ciprofloxacin experiments and phage from mice that had been treated for three weeks for ampicillin experiments along with phage from control mice, respectively. Phage preparations from drug-treated mice and control mice were diluted to the same volume, split equivalently, and incubated with 0.25 ml microbiota (cultured to exponential phase in 0.2% maltose) with 5 mM CaCl2 and 10 mM MgSO4. Phage-microbiota mixtures were allowed to adsorb for one hour at 37 °C (no shaking). Phage-infected microbiota were then pelleted, resuspended in fresh LB, and plated on LB agar plates with ciprofloxacin (1 μg/ml) or LB agar plates with ampicillin (4 μg/ml). A 10 μl aliquot was serially diluted and plated onto no-drug LB agar plates. Frequency was calculated as: # of colonies on drug plate divided by # of colonies on no-drug plate. Additionally, the basal frequency of resistant isolates from microbiota was measured as above in the absence of phage, and this frequency was confirmed to be lower than that from microbiota infected with phage from either treated or untreated mice.
Quantitative PCR was used to measure levels of the 16S rRNA gene in our viral preparations. We used the universal primers 8F and 338R9, and qPCR was performed using the SYBR Green I Master kit and the LightCycler 480 (Roche) according to manufacturer’s instructions.
Contigs were assembled using the Roche 454 GS De Novo Assembler with default parameters, except for a minimum overlap of 100 bp and a minimum identity of 100% to minimize erroneous alignments. Phage-bacterial associations were determined by computing the combination of phage phylogenetic annotations and bacterial phylogenetic annotations on a given contig, and non-redundant phage-bacterial associations were amassed for all contigs. To evenly compare the number of phage-bacterial associations across samples, we performed 50 assemblies using 60,000 randomly selected sequences with replacement from each sample. The data presented in Fig. 4 are representative of these 50 assemblies such that an association was illustrated if it was present in at least one assembly analysis. We computed the mean number of phage-bacterial associations for a given sample (Fig. S8) from the number of non-redundant associations found in each assembly trial. Significance was determined by comparing sample values using the Mann-Whitney U-test. We computed the mean bacteria-to-phage ratio for a given sample (Fig. S9) from the number of bacterial species associated with a given phage for all phage using the union of associations from 50 assemblies shown in Fig. 4. Significance was determined by comparing sample values using the Mann-Whitney U-test.
We thank K. Bodi and J. Schiemer at the Tufts Genomic Core for their sequencing assistance and A. Green and K. Pardee for ultracentrifugation help. This work was supported by the Howard Hughes Medical Institute and the NIH Director’s Pioneer Award Program.
Author ContributionsAll authors designed the study. C.S.S. oversaw the mouse work. S.R.M. and H.H.L. performed and analyzed the experiments with conceptual insight provided by J.J.C. S.R.M., H.H.L., and J.J.C. prepared the manuscript.
Virome data sets reported here are accessible in the NCBI Short Read Archive under accession number SRA073587.
Reprints and permissions information is available at www.nature.com/reprints.
The authors declare no competing financial interests.