Mohammad Mazharul Islam
Samodha C. Fernando
Rajib Saha- 1Department of Chemical and Biomolecular Engineering, University of Nebraska-Lincoln, Lincoln, NE, United States
- 2Department of Animal Science, University of Nebraska-Lincoln, Lincoln, NE, United States
The complex microbial ecosystem within the bovine rumen plays a crucial role in host nutrition, health, and environmental impact. However, little is known about the interactions between the functional entities within the system, which dictates the community structure and functional dynamics and host physiology. With the advancements in high-throughput sequencing and mathematical modeling, in silico genome-scale metabolic analysis promises to expand our understanding of the metabolic interplay in the community. In an attempt to understand the interactions between microbial species and the phages inside rumen, a genome-scale metabolic modeling approach was utilized by using key members in the rumen microbiome (a bacteroidete, a firmicute, and an archaeon) and the viral phages associated with them. Individual microbial host models were integrated into a community model using multi-level mathematical frameworks. An elaborate and heuristics-based computational procedure was employed to predict previously unknown interactions involving the transfer of fatty acids, vitamins, coenzymes, amino acids, and sugars among the community members. While some of these interactions could be inferred by the available multi-omic datasets, our proposed method provides a systemic understanding of why the interactions occur and how these affect the dynamics in a complex microbial ecosystem. To elucidate the functional role of the virome on the microbiome, local alignment search was used to identify the metabolic functions of the viruses associated with the hosts. The incorporation of these functions demonstrated the role of viral auxiliary metabolic genes in relaxing the metabolic bottlenecks in the microbial hosts and complementing the inter-species interactions. Finally, a comparative statistical analysis of different biologically significant community fitness criteria identified the variation in flux space and robustness of metabolic capacities of the community members. Our elucidation of metabolite exchange among the three members of the rumen microbiome shows how their genomic differences and interactions with the viral strains shape up a highly sophisticated metabolic interplay and explains how such interactions across kingdoms can cause metabolic and compositional shifts in the community and affect the health, nutrition, and pathophysiology of the ruminant animal.
Introduction
Within the ruminant species, the microbial community has co-evolved with its host and has helped the host animal obtain energy from low quality fiber rich diets (Hungate, 1975; Hobson, 1988; Henderson et al., 2015). The feed ingested by the ruminant animal undergoes extensive microbial digestion and fermentation in the cattle rumen, producing a range of short chain fatty acids (SCFAs) and energy for the host, and also releases methane, hydrogen, and carbon di-oxide to the atmosphere (Nocek and Russell, 1988). This anaerobic environment is densely (>1010−1011 cells/g of contents) populated by diverse and interdependent species of Bacteria, Protozoa, Fungi, Archaea and Viruses, which are involved in breakdown of complex carbohydrates and polymers from plants, hydrogen transfer, and inter-species transaction of fermentation products and of oligomers and monomers (Bryant and Burkey, 1953; Hungate, 1966; Russell and Rychlik, 2001; Flint et al., 2008). Among the major bacterial players are the phyla Bacteroidetes and Firmicutes, with an abundance of 21% and 12% of the bacterial population in the rumen, respectively (Henderson et al., 2015). Majority of the Bacteroidetes (specifically P. ruminicola and P. bryantii) have a significant role in breakdown of starch and xylan polysaccharides and in the metabolism of proteins and peptides (Wallace et al., 1997; Flint et al., 2008; Thomas et al., 2011). Additionally, Firmicutes including Butyrivibrio fibrisolvens, Ruminococcus flavefaciens, Ruminococcus albas, Eubacterium cellulosolvens etc. play an essential role in the metabolism of cellulose (Flint et al., 2008; Henderson et al., 2015). Archaea are also abundant (approximately 40%) within the rumen and many of them (including Methanobrevibacter gottschalkii, Methanobrevibacter ruminantium, Methanosarcina barkeri, Methanosarcina mazeii etc.) are involved in methanogenesis (St-Pierre and Wright, 2013; Henderson et al., 2015; Seedorf et al., 2015; Danielsson et al., 2017). These methanogens are responsible for the release of methane and other gases in the atmosphere while consuming SCFAs, carbon di-oxide, and hydrogen from other microbes including Bacteroidetes or Firmicutes. The stability and diversity of the rumen microbiome is critical for the animals' health, nutrition, immunity, and survival. Prior studies concluded that disruption in the community composition by sudden administration of grain or glucose to a ruminant animal previously on a dried forage ration often leads to the damage of the rumen tissues and death of the animal within a short period (approximately 18 hours) due to an explosive growth of the Firmicute Streptococcus bovis and Lactobacillus spp., and accumulation of lactic acid (Hungate et al., 1952; Owens et al., 1998; Russell and Rychlik, 2001), thus causing lactic acidosis. In addition, methane production in cattle correlates with rumen methanogenic Archaea and bacterial community structure and dietary composition (Moss et al., 2000; Danielsson et al., 2017). Hence, rumen microbial community remains one of the most interesting and poorly explored natural ecosystems.
The complex interactions between bacterial host and viral phages associated with them drive the ecological dynamics and behavior in many natural systems (Lenski, 1988). Viral modifications of microbial and cyanobacterial (Thompson et al., 2011; Hurwitz and U'Ren, 2016) metabolism was identified in a substantial number of natural systems including marine ecosystem, infectious human diseases, aquifer sediments, and animal gut ecosystems (Sullivan et al., 2006; Suttle, 2007; Minot et al., 2011; Hurwitz et al., 2013; Pan et al., 2014; Weitz et al., 2015; Crummett et al., 2016; De Smet et al., 2016; Howe et al., 2016). Viruses affect intestinal and ruminal microbial ecosystems in cattle through a myriad of processes including cell lysis, energy production, reproduction, and reprogramming of microbial metabolism via Auxiliary Metabolic Genes or AMGs (Berg Miller et al., 2012; Parmar et al., 2016). Studies in recent years have demonstrated that viral AMGs augment the breakdown of complex plant carbohydrates and boost energy production and harvest, while accelerating viral replication inside the host (Anderson et al., 2017). However, a complete understanding of the complex virome-microbiome interaction and the roles of AMGs in the metabolic reprogramming of the host is still in its infancy.
Although advancements in high-throughput sequencing provide access to the vast diversity and makeup of the complex rumen microbial ecosystem, our understanding of the factors that shape rumen microbial communities and interactions among them is lacking. Little is known about the processes shaping the distribution of rumen viruses or the modulation of microbe-driven processes in the rumen. In addition, the study of the rumen ecosystem suffers due to the lack of truly selective media or unique end products, which renders rapid identification of the community composition and metabolic states nearly impossible (Hungate, 1966; Krause and Russell, 1996a,b). The study of ruminal ecology is further complicated by the observation that approximately 75% of the ruminal bacteria are tightly attached to feed particles or are found in biofilms (Hungate, 1966) and cannot be analyzed by only studying fecal contents. On the other hand, many of the virome-microbiome interactions studies appear to be focused on interactions between infectious human viruses and bacteria in an effort to understand respiratory infectious diseases (Levin et al., 1977; Pettigrew et al., 2011; Bosch et al., 2013; Opatowski et al., 2013; Shrestha et al., 2013). A recent study of virus-bacterial interactions in the rumen identified approximately 28,000 viral sequences present in the rumen, which found that the majority of viruses belong to a diversity of viral families including Siphoviridae, Myoviridae, and Podoviridae interacting with host bacterial phyla such as Firmicutes and proteobacteria (Berg Miller et al., 2012).
While a clear understanding of complex biological systems is often challenging, explicit mathematical relation-based modeling promises in silico evaluation and analysis of the biological phenomena. With the gradual increase in computational capacity and the abundance of in silico genome-scale metabolic reconstruction tools, metabolic network models combined with constraint-based analysis provide a host of methods to explore, make discovery and hypotheses, and redesign biological systems at a genome-level (Varma and Palsson, 1994; Mahadevan et al., 2002; Papin et al., 2002; Burgard and Maranas, 2003; Price et al., 2003; Palsson, 2006; Feist et al., 2009; Oberhardt et al., 2009b; Feist and Palsson, 2010; Henry et al., 2010; Orth et al., 2010; Schellenberger et al., 2010; Thiele and Palsson, 2010; Zomorrodi et al., 2012; Maranas and Zomorrodi, 2016; Islam and Saha, 2018). Genome-scale metabolic modeling offers the opportunity not only to map the metabolic landscape of single organisms but also to explore microbe-microbe and phage-microbe interactions. A number of computational tools were developed to model the interactions and dynamics in multi-species microbial communities in the past years (Stolyar et al., 2007; Salimi et al., 2010; Bordbar et al., 2011; Tzamali et al., 2011; Zhuang et al., 2011, 2012; Zomorrodi and Maranas, 2012; Zomorrodi et al., 2014; Henry et al., 2016; Mendes-Soares et al., 2016; Chan et al., 2017; Hanemaaijer et al., 2017; Mendes-Soares and Chia, 2017; Zomorrodi and Segre, 2017; Zuniga et al., 2017; Zengler and Zaramela, 2018). For modeling of multi-species communities, extensive experimental data on different omics' level and a significant knowledge of the inter-species interactions are needed. However, the knowledge-base for microbe-microbe as well as host-microbe interactions is still very poor, and carrying out experiments for studying inter-species interactions appears to be challenging (Fritz et al., 2013). Thus, a combined computational-experimental approach can accelerate new discoveries in the realm of microbe-microbe and microbe-host metabolic interactions. Despite the limitations of current in silico reconstructed host-microbe interaction models, such efforts are of utmost importance because they allow for a detailed metabolic resolution of the complex relationships within microbial communities and with their host. Recently, Heinken et al. (2013) reconstructed and analyzed the first integrated stoichiometric model of murine and Bacteroides thetaiotaomicron metabolism and demonstrated the beneficial interaction of the host and the commensal microbes in the gut. There have been several efforts to reconstruct individual and integrated community models of human gut microbes in recent years, using representative species from dominant classes of microorganisms (Heinken et al., 2013; Shoaie et al., 2013). However, a genome-scale metabolic analysis of the complex community in the rumen was never attempted before.
We developed a simplified and representative rumen community metabolic model with Ruminococcus flavefaciens, Prevotella ruminicola, and Methanobrevibacter gottschalkii as representative organisms from the three major functional guilds in the rumen ecosystem (i.e., Firmicutes, Bacteroidetes, and Archaea, respectively) mentioned above. These three organisms are responsible for fiber digestion, starch and protein digestion, and methane production, respectively. These species were chosen because of their high abundance levels in the rumen, co-abundance reported in previous studies, and the availability of genome annotation (Hungate, 1966; Krause and Russell, 1996a; Wallace et al., 1997; Weimer et al., 1999; Flint et al., 2008; Janssen and Kirs, 2008; Tymensen and McAllister, 2012; Henderson et al., 2015; Thomas et al., 2017; Zhu et al., 2017). We reconstructed the draft models for these species by using the ModelSEED database (Henry et al., 2010), and then performed extensive manual curation, including chemical and charge-balancing, eliminating thermodynamically infeasible cycles, and ensuring network connectivity. The curated models of R. flavefaciens (467 genes, 1033 metabolites, 1,015 reactions), P. ruminicola (546 genes, 1069 metabolites, 1,088 reactions), and M. gottschalkii (319 genes, 900 metabolites, 847 reactions) were integrated into a community model using a multi-level optimization framework (Zomorrodi and Maranas, 2012; Zomorrodi et al., 2014). The community model was used to estimate metabolite secretion profiles and community compositions. To enrich our understanding of the inter-species interactions in the ecosystem, we employed a detailed and comprehensive heuristic procedure that utilized existing GapFind-GapFill tools (Satish Kumar et al., 2007) and a subsequent series of knowledgebase-driven validations. We identified 22 novel interactions involving the transfer of fatty acids, vitamins, coenzymes, amino acids, and sugars among the community members. In addition, we bridged the network gaps in the pentose phosphate pathway, amino acid synthesis and utilization, nucleotide synthesis and degradation, purine metabolism, glycerophospholipid metabolism, and starch metabolism in the metabolic models of these organisms. To elucidate the functional role of the virome on the microbial ecosystem, we used local alignment search and identified metabolic functions of the viruses associated with the community members that drive nucleotide synthesis, reducing power generation, the reprogramming of the bacterial carbon metabolism to pentose phosphate pathway and folate biosynthesis, and viral replication. The identified functions of viral AMGs were incorporated into the model as additional metabolic functions. The addition of viral functionalities resulted in significant changes in bacterial metabolism, including relaxing metabolic bottleneck in the models, complementing microbe-microbe interactions, utilizing nutrients more efficiently and energy harvest by the host. We validated our results based on meta-transcriptomics, meta-proteomics and metabolomics studies on the rumen published recently (Saleem et al., 2012, 2013; Li and Guan, 2017; Wang et al., 2019).
Overall, these findings support the hypothesis that viral AMGs play a crucial role in enhancing host fitness and robustness. We also studied the effect of using different community-level objective functions (i.e., growth, short-chain fatty acids production, plant feed utilization, greenhouse gas release, and small sugar molecule production) on the metabolic capacity of the community members. We found that the flux ranges of the microbial species are robust irrespective of the choice of a community objective. Hence, maximizing the community biomass is a rational choice since a stable community in the rumen needs to survive and grow at a reasonable rate to perform its necessary role in host nutrition and pathophysiology despite constant washout events like fecal secretion.
Results
Genome-Scale Models of R. flavefaciens, P. ruminicola, and M. gottschalkii
Three genome-scale metabolic models of representative organisms of the rumen microbiome were reconstructed for P. ruminicola, R. flavefaciens, and M. gottschalkii using the ModelSEED database (Henry et al., 2010). This was followed by manual editing, refinements and further curation of the models. To this end, the draft reconstructions contained network gaps and 3-25% of reactions with chemical imbalances. One of them (R. flavefaciens) had thermodynamically inconsistent production of redox cofactors, which warranted an extensive manual curation to be performed on the models (see methods section for details). The manual curation steps ensured that there is no chemical or charge imbalance present in the models. The numbers of reactions that carries unrealistically high fluxes without any nutrient uptake (thus defined as reactions with unbounded fluxes) also decreased substantially (79% for R. flavefaciens, 79% for P. ruminicola, and 50% for M. gottschalkii). The rest of the unbounded reactions from the nucleotide degradation pathways were not fixed since they do not critically affect the biological significance of the models and are usually common in existing metabolic models (Fritzemeier et al., 2017). The initial manual curation process reconnected a number of blocked reactions in the models, which will be addressed in the subsequent sections. The reconciliation of cofactor inconsistency in R. flavefaciens model was performed by removing spurious duplicate reactions that caused unlimited ATP generation, thereby ensuring that all models are energetically consistent. The summary of model statistics before and after the manual curation are shown in Table 1. Supplementary Materials 1–3 contain the details of the microorganism models.
Table 1. Statistics for the draft and curated models of R. flavefaciens, P. ruminicola, and M. gottschalkii.
Metabolic Interactions in the Microbial Community
Individual metabolic models of the three species were integrated into an interacting community framework. This in silico community model was simulated at the growth condition estimated for a standard-sized domestic cow (see details in Methods section). The simulated community biomass flux (0.513 h−1) was comparable to the experimentally observed values of passage/dilution rates of 0.043–1.0 h−1 in the rumen depending on dietary regimen and other factors (Goetsch and Galyean, 1982; Stokes et al., 1985; Tellier et al., 2004). Dilution rate or passage rate of ruminal content is defined as the rate at which the digesta leaves a compartment in the gut (in this case, rumen) (Goetsch and Galyean, 1982; Tellier et al., 2004). It is important that the rumen microbiome growth rate (the community biomass flux) is within the margin of experimentally observed dilution rate because a stable community needs to survive and grow at a reasonable rate to perform its necessary role in host nutrition and pathophysiology despite constant washout events like fecal secretion. The extent and directionalities of the major metabolic exchanges in the community are shown in Figure 1, which shows that the model was able to capture the known rumen microbiome behavior in terms of metabolite production, consumption, and inter-species transactions (Bryant and Burkey, 1953; Hungate, 1966; Russell and Rychlik, 2001; Flint et al., 2008; Henderson et al., 2015). Complex plant carbohydrates and proteins are utilized by R. flavefaciens and P. ruminicola, and the short-chain fatty acids (SCFAs) are absorbed by the rumen epithelium. M. gottschalkii accepts hydrogen, formate, carbon di-oxide from the other two members, and releases methane to the atmosphere, similar to what was described in previous reports (St-Pierre and Wright, 2013; Henderson et al., 2015; Seedorf et al., 2015; Danielsson et al., 2017). It should be noted that all the metabolic transaction may not be active in all physiological condition of the community. When optimizing for the overall community biomass, the flux distribution is shaped in a way that optimizes the production of biomass and not for other carbon molecules like short-chain fatty acids or small sugars. At the same time, when optimizing for other community objectives like production of SCFAs, some of the inactive fluxes including these metabolic transactions can become active. A comparative analysis of the effect of different community objective functions is presented in subsequent sections.
Figure 1. Initial community simulation results showing the interactions between the bacterial and archaeal members. Rfl, Prm, and Mgk represent R. flavefaciens, P. ruminicola, and M. gottschalkii, respectively. The numbers inside the circles for each microbe represent the biomass flux (growth rate) of the respective microbe (hr−1). The arrows represent metabolic fluxes in mmol/gDCW.hr (dashes for inter-species/shared metabolites and solids for transfer to and from the rumen epithelium). The numbers along the arrows represent the minimum and maximum flux values.
De novo Interactions
The broad range of metabolic capabilities of different microorganisms in the rumen ecosystem is indicative of numerous inter-species interactions at metabolic, signaling, and regulatory levels. While the major functions of Bacteroides, Firmicutes and Archaea, and the interplay between them is commonly known (Hungate, 1975; Wolin, 1979; Flint, 1997; Raizada et al., 2003; Henderson et al., 2015; Nagaraja, 2016; Anderson et al., 2017; Liu et al., 2017), the current knowledge is only partial. As described in the methods section, a detailed step-by-step procedure was employed to identify 22 potential metabolic transactions. The overarching idea of the procedure developed is that complementation of the metabolic functions among community members gives rise to inter-species interactions. The identified de novo interactions are shown in Figure 2 and listed in detail in Supplementary Material 4. In addition to the transfer of sugar monomers from R. flavefaciens to P. ruminicola, a number of co-enzymes and vitamins were found to be exchanged between M. gottschalkii and other members. Model statistics after filling the gaps using the GapFind-GapFill algorithms are shown in Table 2.
Figure 2. Identified de novo interactions in the community. Metabolites in black text were previously known to be exchanged, metabolites in color text are identified in this work. A cartoon inside the circles shows the main pathway map for each organism.
Table 2. Model statistics after filling the gaps using GapFind-GapFill.
The GapFill optimization procedure was repeated after the metabolic functions from the viruses were incorporated into the metabolic model of each of the microbes. Among the 22 identified de novo interactions as mentioned earlier, 20 were retained upon gap filling after adding these virome functionalities. The transfer of D-mannose and fructose from R. flavefaciens to M. gottschalkii did not appear in these gap filling results.
Viral Auxiliary Metabolic Genes and Shifts in Flux Distributions in the Metabolic Models
Upon filtering more than 3,000 candidate proteins from BLAST search results for high bit score and expectancy values (<10−34), the models of R. flavefaciens, P. ruminicola, and M. gottschalkii were amended with 29, 26, and 18 metabolic reactions, respectively. For a detailed result of the alignment search (see Supplementary Material 5). The reactions included additional metabolic functions and also novel metabolic capacities in the existing models. The addition of viral reactions resulted in significant changes in flux distributions. Up to 11% of reaction in all three of the models (see Supplementary Material 6 for details) had increased their flux ranges significantly, either by decreasing the minimum flux or by increasing the maximum flux or both. Reactions with a standard deviation of their minimum or maximum flux values >1 across the simulations under different community objective functions (described in the following section) were considered at this stage. Overall, all the reactions that had a changed flux distribution were relaxed in AMG-amended models compared to the post-gapfilled models. Figure 3 illustrates the change of flux space upon addition of viral AMGs.
Figure 3. Changes in flux space after viral AMGs were added to the metabolic models. The lighter shade for each of the colors represent the relaxed flux space after addition of AMGs in each of the microbial metabolic models.
Variability of the Metabolic Fluxes Under Different Community Objective Functions
The metabolic spaces of the community members were assessed under different community objective functions after viral AMGs were added to the host models. To simulate different functions of the rumen microbial community, i.e., (i) growth, (ii) SCFA production, (iii) feed utilization, (iv) methane and carbon-di-oxide release, and (v) small sugar molecule production, the appropriate community-level objective functions were chosen and optimized (see Figure 4A for details). Statistical analysis was performed to investigate the deviation of flux ranges under these five objective functions. Approximately 44, 18, and 6% of the reactions from R. flavefaciens, P. ruminicola, and M. gottschalkii models, respectively, had non-zero standard deviation among their flux ranges in these conditions (see Supplementary Material 6 for details). However, under these five objective functions, the variation of flux distributions across the three community members was very small (with average standard deviation of 0.0001 and the range between 0 and 0.68). As representatives of the corresponding flux ranges, Figures 4B–D show the density functions of the standard deviations of the exchange fluxes (i.e., both uptakes and imports) from all of the three community members. From this analysis, it is evident that while viral AMGs play an important role in relaxing the flux ranges, the overall flux spaces of the intracellular and extracellular reactions in a community are very robust against the choice of a community level objective function.
Figure 4. Variability in metabolic fluxes under different community objective functions. (A) Visual representation showing the choice of different community-level objective functions. The density functions (right) show the insignificance of the variations in exchange flux space upon optimizing for different objective functions, i.e., maximizing total community biomass, maximizing total Short-chain Fatty Acids (SCFA) production, maximizing total complex carbohydrate uptake, minimizing total Methane and carbon-di-oxide production, and maximizing total sugar production by the community: (B) R. flavefaciens, (C) P. ruminicola, and (D) M. gottschalkii. The bandwidth is the standard deviation of the smoothing kernel of the density function.
Discussion
Individual Model Performance and Community Interactions
The reconstructed and curated genome-scale metabolic models of P. ruminicola, R. flavefaciens, and M. gottschalkii were evaluated for fitness and species-specific metabolic capacities. Each of the metabolic models were simulated for growth and energy production at the growth condition estimated for a standard-sized domestic cow (see methods). The simulated growth rates were 0.086, 0.065, and 0.362 hr−1 for R. flavefaciens, P. ruminicola, and M. gottschalkii, respectively. While an experimental growth rate observation for these microbes are not available, the growth rates are comparable to the reported dilution rates in the rumen (Huntington et al., 1981; Estell and Galyean, 1985). The degradation of plant cellulose, biosynthesis of branched chain amino acid and short-chain fatty acid (SCFA) and production of hydrogen were predicted by R. flavefaciens model and were in agreement with experimental observations (Helaszek and White, 1991; Flint et al., 2008; Zheng et al., 2014). Protein digestion by P. ruminicola and methane production from carbon di-oxide and hydrogen consumption by M. gottschalkii were also validated (Wallace et al., 1997; Qiao et al., 2014; Henderson et al., 2015; Seedorf et al., 2015).
The representative community successfully captured most of the known metabolic interactions with the three functional guilds of microbes, as was evident from the flux values of the shared metabolites (Figure 1). While each of the community members were tested for known metabolic contribution in the community, a community simulation with total biomass as the biological objective rendered some of the metabolic transactions inactive. For example, even though the production of succinate and propionate by R. flavefaciens and P. ruminicola were validated, the maximization of community biomass did not drive the production of those SCFAs for growth. It should be noted that despite both P. ruminicola and M. gottschalkii being known and tested in silico for carbon di-oxide production in rumen, the community simulated showed that only P. ruminicola produced all the carbon di-oxide that was released from the rumen. This shows that metabolic exchanges can have separate dynamics based on community fitness criteria set during the optimization process.
Insight Into de novo Community Interactions
The predicted interactions between R. flavefaciens and P. ruminicola include the transfer of various small sugar monomers and fatty acids (as shown in Figure 2). These interactions are highly warranted given the metabolic functions each of these organisms perform (Wallace, 1992; Nagaraja, 2016). Both of these cellulolytic organisms contribute toward breakdown of complex plant material and the production of SCFAs and small sugar molecules for the host and other members in the community. Several NMR and GC-MS based metabolomic analyses of the ruminal fluid showed the secretion of glucose, mannose, acetate, formate, and maltoheptaose in the rumen (Saleem et al., 2012, 2013), which match our observations. These interactions build mutualistic and commensal relationships in the rumen ecosystem. Prevotella is also known for its proteolytic functions and high dipeptidyl peptidase activity (Nagaraja, 2016), which was demonstrated by the consumption and degradation of amino acids in the community simulation. This observation agrees with the meta-transcriptomic analysis by Li and Guan (Li and Guan, 2017), where they found downregulated amino acid synthesis and upregulated alanine, aspartate and glutamate catabolism in Prevotella, and the metabolomic study by Wang et al. (2019), where a co-occurrence network analysis among the microbiota and metabolites showed positive correlation between Prevotella and several amino acids. Li and Guan (Li and Guan, 2017) also observed downregulated fatty acids synthesis pathway in Prevotella, which can be complemented by the uptake of fatty acids like octadecanoate produced by R. flavefaciens and M. gottschalkii. Prevotella also lacks the capability of glycan degradation (Li and Guan, 2017), which provides the explanation for our in silico observation that glycans are degraded by R. flavefaciens and then the degradation products are taken up by Prevotella. Methanobrevibacter acts as a sulfate reducer and a hydrogenotroph (Samuel and Gordon, 2006; Hansen et al., 2011; Zheng et al., 2014; Nagaraja, 2016), which can be attributed to its consumption of thiosulfate, formate, and hydrogen. The consumption of hydrogen by Methanobrevibacter may potentially facilitate the extraction of energy from nutrients and increase digestion efficiency by redirecting the ruminal fermentation toward more oxidized end products (Armougom et al., 2009). Previous studies have inferred an absence of Glycosaminoglycan degradation functions in M. gottschalkii (Li and Guan, 2017), which explains our observation of acetyl glucosamine secretion by M. gottschalkii and subsequent consumption by P. ruminicola. At the same time, the role of Methanobrevibacter in supplying a number of important amino acids, vitamins and co-enzymes to the other organisms was observed, which was not reported before the current study. These roles suggest that Methanobrevibacter is a major player in the rumen ecosystem and warrants further studies on its role in specificity and efficiency of bacterial digestion in rumen.
Relaxation of Metabolic Bottlenecks in the Presence of Viral Auxiliary Metabolic Genes
The addition of Auxiliary Metabolic Genes (AMGs) from the viruses associated with each of the host community members had profound impact on shifting and relaxing the flux ranges of key rate-limiting steps in the host organisms, termed as “metabolic bottlenecks.” As seen in Figure 3, in R. flavefaciens, 127 reactions were observed to have wider flux ranges after AMGs were added, compared to the gap-filled model. Viral AMGs complemented important function in the host model involving nucleotidyltransferases and DNA repair (details in Supplementary Material 5). The pathways that were relaxed due to viral metabolic genes are Calvin-Benson cycle (CBB), amino acid biosynthesis, Sugar utilization, nucleoside catabolism, fermentation, glycolysis/gluconeogenesis, cofactor biosynthesis, single carbon metabolism (tetrahydropterines), and various transport mechanisms. In addition to that, Pentose Phosphate Pathway (PPP) carried higher flux compared to the gapfilled model. Therefore, viral AMGs not only add some additional metabolic functions to host but also complement and relax some key pathways that drive the fitness of the host. In P. ruminicola, 25 reactions were relaxed after AMGs were added, compared to the gap-filled model. In addition to the nucleotide polymerization reactions, viral AMGs also complemented complex sugar breakdown pathways, as listed in Supplementary Material 5. The pathways that were relaxed due to viral metabolic genes are folate biosynthesis, energy production, cofactor synthesis, and glycogen synthesis, which are important for boosting the energy production in the microbes and the generation of reducing power. It has also been previously hypothesized that these phenomena aid in viral replication (Anderson et al., 2017). In M. gottschalkii, 11 reactions had wider flux ranges after viral AMGs responsible for propanoate and amino acids metabolism were added. The relaxed pathways include amino acid biosynthesis and degradation, and nucleotide metabolism. In addition, coenzyme biosynthesis, glycine and serine degradation, diaminopimelic acid (DAP) pathway for lysine synthesis, methanogenesis, methionine biosynthesis, peptidoglycan biosynthesis, pentose phosphate pathway, SCFA production, and purine and pyrimidine conversion pathways carried higher fluxes compared to the post-gapfilled model. These show the important role of viral AMGs in cell growth and replication for the host (M. gottschalkii).
The inclusion of viral AMGs as metabolic functions resulted in noticeable changes in the inter-species transports, as manifested by the change in the flux ranges of the shared metabolite transactions (see Supplementary Material 6 for details). The overall changes in metabolic interactions among R. flavefaciens, P. ruminicola, and M. gottschalkii are shown in Figure 5. While M. gottschalkii could produce biotin, nitrite, and octadecanoate at a high rate, P. ruminicola became less reliant on these metabolites after addition of viral AMGs, possibly due to the complementation in the galacturonan and pectin digestion ability and alanine, aspartate, and glutamate metabolism. Galacturonate and Deoxyadenosine were not consumed by R. flavefaciens at the same rate they were being produced by P. ruminicola because R. flavefaciens was able to digest the complex carbohydrates in the plant feed more efficiently after addition of viral AMGs. The transfer of mannose and fructose sugars to M. gottschalkii from R. flavefaciens was not identified as interactions, which ultimately resulted in reduction of the number of de novo interactions in the community to 20 when viral AMGs were present. Upon inspection of the functionalities coming from the phages associated with M. gottschalkii, it became apparent that the D-glyceraldehyde-3-phosphate glycolaldehyde transferase and D-mannose-6-phosphate aldose-ketose-isomerase system that drive the conversion between fructose and mannose became unnecessary in the presence of virome functionalities in the models. At the same time, the phosphofructokinase/sedoheptulokinase (EC 2.7.1.14) enzymes, including Sedoheptulose 1,7-bisphosphate D-glyceraldehyde-3-phosphate-lyase and ATP:Sedoheptulose 7-phosphate 1-phosphotransferase could shift the carbon flux from glycolysis (via fructose-6-phosphate) to Pentose Phosphate pathway through tetrose sugars and therefore facilitate the adaptation of M. gottschalkii even without the support of sugars from R. flavefaciens. The ability to maintain species fitness without these interactions leads to better community fitness in situations when these interactions are threatened (e.g., antibiotic treatment).
Figure 5. Shifts in metabolism and inter-species interactions after the inclusion of viral auxiliary metabolic genes. Inside the circles for each organism, increase in pathway fluxes are shown in thicker green lines and decrease in pathway fluxes are shown in purple lines. Decreased metabolic transactions are shown in thinner dashed lines.
Robust Community Behavior Under Different Community Objectives
While a large number of reactions had varying flux ranges under different community-level objective functions, the change was very negligible, with an average standard deviation of 0.0001. The variation in the number of reactions with non-zero standard deviation between different metabolic models could possibly be correlated to their metabolic roles in the community and the choice of the community-level objective function. For example, R. flavefaciens is a degrader of complex carbohydrate and producer of small sugar molecules, which thus drives the growth of other community members as well as the SCFA production. Therefore, most of the community-level objectives we studied, i.e., growth, SCFA production, feed utilization, and small sugar molecule production, were directly related to this species. It seems logical that switching between these objective functions will largely affect R. flavefaciens' flux space. On the other hand, one of the community-level objectives (methane and CO2 production) was directly related to the metabolic profile of M. gottschalkii, and a very low percentage (6%) of reactions had varied flux ranges under all the five different community-level objectives in this species.
Although living organisms have evolved to maximize their survival in varying environments, it is difficult to assess the primary driver of the flux distributions of primary metabolism (Burgard and Maranas, 2003). Disagreement is prevalent among the scientific community as to whether metabolic flux in a network is distributed to satisfy optimal biomass production, maximum energy generation, or most efficient utilization of substrates. It is especially true for a naturally occurring microbial community such as the rumen ecosystem in which thousands of microbial, viral, fungal, and archaeal species co-evolved with different metabolic functions and interplay among one another. Therefore, maximizing community biomass may not always be the best choice for an objective function. However, choosing the best community objective is not as straightforward as it seems. With that in mind, the variation of flux ranges under different community objective functions was studied. From Figure 4 (also Supplementary Material 6), it is observed that the choice of a community objective has negligible effect while simulating the possible metabolic capacities (flux space) of the community members. While there have been attempts at optimizing a number of different objective functions for prokaryotic microorganisms, optimizing for growth often seems realistic for both eukaryotic and prokaryotic organisms (Burgard and Maranas, 2003). Due to the lack of knowledge about the overall goal of the rumen community, maximizing community biomass is a logical choice, given that a stable community in the rumen needs to survive and grow at a reasonable rate to perform its necessary role in host nutrition and pathophysiology despite constant washout events like fecal secretion.
Moving Ahead: Challenges and Future Prospects
In summary, we developed a workflow to elucidate novel interactions between participating members of a complex community and its phages. Our in silico predictions of metabolite exchange between three members of rumen microbiome agrees with the current knowledgebase about their metabolic functionalities and roles in the rumen ecosystem as well as recently published multi-omics datasets. We also identified possible viral auxiliary metabolic genes associated with the three members in our community that reinforced their metabolic capacities and helped in relaxing several bottlenecks in the metabolic networks. This was manifested by the enhancement of reaction fluxes in important metabolic pathways in the models and the metabolic robustness achieved by the microbes. Our community metabolic model serves to discover unidentified metabolite transactions and answer key ecological questions of ruminant nutrition through virome-microbiome interactions, while promising to address important biological aspects of ruminant nutrition and greenhouse emission. The development of additional bioinformatics tools, advancement in high-throughput sequencing technologies and cultivation-independent “omics” approaches may drive further development of new mathematical frameworks for analyzing rumen ecosystems. In vitro cell culture systems and in situ experimentation of simplified microbial, viral, and fungal communities in gnotobiotic bovine rumen to capture spatiotemporal community dynamics will enhance our analyzing power to a greater extent. Therefore, what we need is a combination of computational and experimental efforts to enrich the current knowledgebase regarding in situ metabolic and taxonomic profiling, species identification and characterization, annotation, and advanced tools to accommodate for large-scale data analysis and integration, which will accelerate further efforts in deciphering the complexity of this ecosystem.
Methods
Model Reconstructions and Refinements
The initial draft genome-scale metabolic reconstructions of P. ruminicola, R. flavefaciens, and M. gottschalkii were created and downloaded from the ModelSEED biochemical database (in September 2017). The models included reactions for central carbon metabolism, secondary biosynthesis pathway, energy and cofactor metabolism, lipid synthesis, elongation and degradation, nucleotide metabolism, amino acid biosynthesis and degradation. Flux Balance Analysis (FBA) was employed during model testing, validation, and analyzing flux distributions at different stages of our work (Varma and Palsson, 1993, 1994; Oberhardt et al., 2009a). The optimization formulation, in its most common form, is given below.
Here, I and J are the sets of metabolites and reactions in the metabolic model, respectively. Sij is the stoichiometric coefficient of metabolite i in reaction j and vj is the flux value of reaction j. Parameters LBj and UBj denote the minimum and maximum allowable fluxes for reaction j, respectively. vbiomass is the flux of the biomass reaction which mimics the cellular growth yield.
Curation of Metabolic Models
Correcting Reaction Imbalances
For balancing the reactions imbalanced in protons, we checked for the protonation state consistent with the reaction set in the draft model and performed addition/deletion of one or multiple protons on either the reactant or the product side. For the remaining imbalanced reactions, we corrected the reaction stoichiometry in order to ensure that the atoms on both sides of the reactions balance out.
Identifying and Eliminating Thermodynamically Infeasible Cycles
One of the limitations of constraint-based genome-scale models is that the mass balance constraints only describe the net accumulation or consumption of metabolites, without restricting the individual reaction fluxes. While biochemical conversion cycles like TCA cycle or urea cycle are ubiquitous in a metabolic network model, there can be cycles which do not consume or produce any metabolite. Therefore, the overall thermodynamic driving force of these cycles are zero, implying that no net flux can flow around this cycle (Schellenberger et al., 2011). To identify Thermodynamically Infeasible Cycles in our model, we turned off all the nutrient uptakes to the cell and used an optimization formulation called Flux Variability Analysis (FVA) which maximizes and minimizes each of the reaction fluxes subject to mass balance constraints (Mahadevan and Schilling, 2003). The reaction fluxes which hit either the lower bound or upper bound are defined as unbounded reactions and were grouped together as a linear combination of the null basis of their stoichiometric matrix. To eliminate the cycles, we either removed duplicate reactions, turned off lumped reaction or selectively turned reactions on/off based on available cofactor specificity information (see Supplementary Material 7 for details). The mathematical formulation of Flux Variability Analysis (FVA) is given below. vapp−obj, threshold (from constraint 3, which is optional) is a predetermined threshold value of the appropriate objective flux vapp−obj to ensure that the feasible flux space satisfy the targeted value.
Ensuring Known Metabolic Functions
Each of the three metabolic models were checked for the capacity to produce biomass and metabolites they were known to produce (Helaszek and White, 1991; Wallace et al., 1997; Flint et al., 2008; Qiao et al., 2014; Zheng et al., 2014; Henderson et al., 2015; Seedorf et al., 2015). To ensure these metabolic functionalities, the reactions missing in metabolic pathways were systematically identified and added manually from biochemical databases (Kanehisa and Goto, 2000; Apweiler et al., 2004; Henry et al., 2010) after an extensive search for each of the missing enzyme activities in related organisms such as Ruminococcus albus, Bacteroides thetaiotaomicron, and Methanobrevibacter smithii for R. flavefaciens, P. ruminicola, and M. gottschalkii, respectively. A missing metabolic function was only added if the genes between any pair of organisms were found to be orthologous and amending the models did not result in an increase in the number of thermodynamically infeasible cycles.
Community Formation and Simulation
Once the individual microbe models were curated, they were integrated to form a community model using existing optimization framework, namely OptCom (Zomorrodi and Maranas, 2012). In this bi-level multi-objective optimization framework, the individual flux balance analysis problems for each community member are treated as inner-level optimization problem, and the community objective is optimized for in the outer-level problem. The mathematical description of the OptCom procedure is given below.
The inner problem(s) represents the steady-state FBA problem for each of the community members k and limits on uptake or export flux of a shared metabolite using the parameters and , respectively, which are imposed by the outer problem. Constraint (6) in the outer problem describes a mass balance for each shared metabolite in the extra-cellular environment of shared metabolite pool, where the terms and represent the total uptake and export of the shared metabolite i by community members, respectively. These constraints are the key equations for modeling the known metabolic interactions among participants of the community. When further information on de novo interaction became available, they were also incorporated using these outer-level constraints. Incorporating the known interaction in the community, a rumen ecosystem for a 1,000 lb cow using the three-member simplified community was simulated. The community nutrient uptakes were designed based on the data obtained from previous experiments (Anderson et al., 2017) and details are given in Supplementary Material 8.
Identification of Viral Auxiliary Metabolic Genes
In addition to the syntrophic, mutualistic, and competitive microbial interactions, viruses impact microbial populations through cell lysis and reprogramming of host metabolism by Auxiliary Metabolic genes (AMGs). Most of the bacterial species in the cattle rumen have their own associated viruses (Berg Miller et al., 2012; Ross et al., 2013; Parmar et al., 2016). To identify the viral AMGs, the viruses/phages associated with each of the community members were searched for. Berg Miller et al. (2012) suggested 13 different phages associated with R. flavefaciens, P. ruminicola and M. gottschalkii. A local alignment search (BLAST) of the viral proteomes downloaded from several databases (Apweiler et al., 2004; Uniprot, 2007; Hubisz et al., 2011; Zhou et al., 2011; Brister et al., 2015; Arndt et al., 2016) to The National Center for Biotechnology Information (NCBI) non-redundant proteins sequence database (O'Leary et al., 2016) was performed. The search yielded more than 3,000 candidate proteins, which were filtered for expectancy values (<10−34).
Identification of Unknown Interactions and Bridging of Network Gaps
One of the key limitations of any genome-scale reconstruction is the presence of gaps in the resultant metabolic network. These gaps occur when there are dead-end metabolites in the model, meaning that the consumption and/or production reactions for a particular metabolite are absent in the model. Any such gaps from a specific model could be reconciled by borrowing the required metabolic functionalities from the models of neighboring organism. For each of the metabolic models in the rumen community, a set of Mixed-integer Linear Programming (MILP) optimization procedures (named GapFind and GapFill) were used to identify and eliminate network gaps in these reconstructions (Satish Kumar et al., 2007). To identify the unknown inter-species interactions in the community, a protocol was developed that takes each suggestion from the GapFill algorithm and performs a series of tests (shown in Figure 6) to categorize the results as unacceptable, possible inter-species interactions, or acceptable solutions to bridge network gaps. The metabolites that cannot be produced or consumed in a network are called the problem metabolites. For the problem metabolites in each of the member species in the community, we performed three separate sets of gap filling procedures. Two of them used the two other community members as the source database, and one used the Modelseed biochemical database (Henry et al., 2010) as the source database (current as of February 2018). For each of the gap filling solutions, if it was from the Modelseed database, it was tested for possible transport mechanisms (Ren et al., 2007), likelihood based on presence in taxonomically related organisms (with prioritizing organisms in the same lower taxonomic levels), and whether the solution created thermodynamically infeasible cycles, and finally that solution was either accepted or rejected. Similarly, for a gap filling suggestion coming from another member in the community, the same sets of checks were performed, and the solution was either accepted or rejected. However, for these gap filling suggestions, if the solution was a transport mechanism and the target organism was known to have a transport mechanism as per the current knowledgebases (Ren et al., 2007), it was concluded to be a de novo identification of a potential interaction that exists in the community. It should be noted that for each transport reaction, strong bioinformatic evidence or clear ortholog of experimentally characterized transporter in a closely related organism was ensured before it was accepted as a gapfill solution. When the individual metabolic model for each of the microbes was augmented with viral AMGs, the interaction identification procedure (shown in Figure 6) was repeated to compare the shifts in inter-species interactions. Therefore, the gapfiling procedure and the protocol described above served two very important steps in the model curation process by bridging network gaps as well as identifying possible metabolic interactions between organisms.
Figure 6. Workflow for identifying possible interspecies interactions and filtering GapFill suggestions.
Computational Resources
The General Algebraic Modeling System (GAMS) version 24.7.4 with IBM CPLEX solver was used to run all the optimization algorithms in this work. For each of the algorithms, the required optimization framework was scripted in GAMS and then run on a Linux-based high-performance cluster computing system at the University of Nebraska-Lincoln. The models were parsed from Systems-Biology Markup Language (SBML) level 2 document using standard programming languages (i.e., Python) to generate the input files required by GAMS.
Data Availability Statement
All datasets generated or analyzed for this study are included in the manuscript and the Supplementary Files.
Author Contributions
RS and SF conceived the study. RS supervised the study. MI performed all the experiments and analyzed the results. MI, SF, and RS wrote the manuscript. All authors read and approved the manuscript.
Funding
This study is based upon work supported by University of Nebraska-Lincoln Faculty Startup Grant (21-1106-4308) to Rajib Saha and Animal Nutrition, Growth and Lactation Grant no. 2018-67015-27496, Effective Mitigation Strategies for Antimicrobial Resistance Grant no. 2018-68003-27545, and Multi-state research project accession no. 1000579 from the USDA National Institute of Food and Agriculture awarded to SF.
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Acknowledgments
The authors would like to thank Mr. Wheaton Schroeder and Mr. Matthew Van Beek from Chemical and Biomolecular Engineering Department at the University of Nebraska-Lincoln for their help during the model curation step.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2019.02412/full#supplementary-material
Supplementary Material 1. Model files for R. flavefaciens.
Supplementary Material 2. Model files for P. ruminicola.
Supplementary Material 3. Model files for M. gottschalkii.
Supplementary Material 4. Metabolic interactions in the community identified in the current study.
Supplementary Material 5. Basic Local Alignment (BLAST) search results for viral AMGs to be incorporated in the microbial host models.
Supplementary Material 6. Flux distributions in the models at different stages of gap filling and addition of viral AMGs.
Supplementary Material 7. Strategies for breaking infeasible cycles.
Supplementary Material 8. Details of the community nutrient uptake calculations.
References
Anderson, C. L., Sullivan, M. B., and Fernando, S. C. (2017). Dietary energy drives the dynamic response of bovine rumen viral communities. Microbiome 5:155. doi: 10.1186/s40168-017-0374-3
Apweiler, R., Bairoch, A., Wu, C. H., Barker, W. C., Boeckmann, B., Ferro, S., et al. (2004). UniProt: the Universal Protein knowledgebase. Nucleic Acids Res. 32, D115–D119. doi: 10.1093/nar/gkh131
Armougom, F., Henry, M., Vialettes, B., Raccah, D., and Raoult, D. (2009). Monitoring bacterial community of human gut microbiota reveals an increase in Lactobacillus in obese patients and Methanogens in anorexic patients. PLoS ONE 4:e7125. doi: 10.1371/journal.pone.0007125
Arndt, D., Grant, J. R., Marcu, A., Sajed, T., Pon, A., Liang, Y., et al. (2016). PHASTER: a better, faster version of the PHAST phage search tool. Nucleic Acids Res. 44, W16–W21. doi: 10.1093/nar/gkw387
Berg Miller, M. E., Yeoman, C. J., Chia, N., Tringe, S. G., Angly, F. E., Edwards, R. A., et al. (2012). Phage-bacteria relationships and CRISPR elements revealed by a metagenomic survey of the rumen microbiome. Environ. Microbiol. 14, 207–227. doi: 10.1111/j.1462-2920.2011.02593.x
Bordbar, A., Feist, A. M., Usaite-Black, R., Woodcock, J., Palsson, B. O., and Famili, I. (2011). A multi-tissue type genome-scale metabolic network for analysis of whole-body systems physiology. BMC Syst. Biol. 5:180. doi: 10.1186/1752-0509-5-180
Bosch, A. A., Biesbroek, G., Trzcinski, K., Sanders, E. A., and Bogaert, D. (2013). Viral and bacterial interactions in the upper respiratory tract. PLoS Pathog. 9:e1003057. doi: 10.1371/journal.ppat.1003057
Brister, J. R., Ako-Adjei, D., Bao, Y., and Blinkova, O. (2015). NCBI viral genomes resource. Nucleic Acids Res. 43, D571–D577. doi: 10.1093/nar/gku1207
Bryant, M. P., and Burkey, L. A. (1953). Cultural methods and some characteristics of some of the more numerous groups of bacteria in the bovine rumen. J. Dairy Sci. 36, 205–217. doi: 10.3168/jds.S0022-0302(53)91482-9
Burgard, A. P., and Maranas, C. D. (2003). Optimization-based framework for inferring and testing hypothesized metabolic objective functions. Biotechnol. Bioeng. 82, 670–677. doi: 10.1002/bit.10617
Chan, S. H. J., Simons, M. N., and Maranas, C. D. (2017). SteadyCom: predicting microbial abundances while ensuring community stability. PLoS Comput. Biol. 13:e1005539. doi: 10.1371/journal.pcbi.1005539
Crummett, L. T., Puxty, R. J., Weihe, C., Marston, M. F., and Martiny, J. B. H. (2016). The genomic content and context of auxiliary metabolic genes in marine cyanomyoviruses. Virology 499, 219–229. doi: 10.1016/j.virol.2016.09.016
Danielsson, R., Dicksved, J., Sun, L., Gonda, H., Muller, B., Schnurer, A., et al. (2017). Methane production in dairy cows correlates with rumen methanogenic and bacterial community structure. Front. Microbiol. 8:226. doi: 10.3389/fmicb.2017.00226
De Smet, J., Zimmermann, M., Kogadeeva, M., Ceyssens, P. J., Vermaelen, W., Blasdel, B., et al. (2016). High coverage metabolomics analysis reveals phage-specific alterations to Pseudomonas aeruginosa physiology during infection. ISME J. 10, 1823–1835. doi: 10.1038/ismej.2016.3
Estell, R. E., and Galyean, M. L. (1985). Relationship of rumen fluid dilution rate to rumen fermentation and dietary characteristics of beef steers. J. Anim. Sci. 60, 1061–1071. doi: 10.2527/jas1985.6041061x
Feist, A. M., Herrgard, M. J., Thiele, I., Reed, J. L., and Palsson, B. O. (2009). Reconstruction of biochemical networks in microorganisms. Nat. Rev. Microbiol. 7, 129–143. doi: 10.1038/nrmicro1949
Feist, A. M., and Palsson, B. O. (2010). The biomass objective function. Curr. Opin. Microbiol. 13, 344–349. doi: 10.1016/j.mib.2010.03.003
Flint, H. J. (1997). The rumen microbial ecosystem–some recent developments. Trends Microbiol. 5, 483–488. doi: 10.1016/S0966-842X(97)01159-1
Flint, H. J., Bayer, E. A., Rincon, M. T., Lamed, R., and White, B. A. (2008). Polysaccharide utilization by gut bacteria: potential for new insights from genomic analysis. Nat. Rev. Microbiol. 6, 121–131. doi: 10.1038/nrmicro1817
Fritz, J. V., Desai, M. S., Shah, P., Schneider, J. G., and Wilmes, P. (2013). From meta-omics to causality: experimental models for human microbiome research. Microbiome 1:14. doi: 10.1186/2049-2618-1-14
Fritzemeier, C. J., Hartleb, D., Szappanos, B., Papp, B., and Lercher, M. J. (2017). Erroneous energy-generating cycles in published genome scale metabolic networks: identification and removal. PLoS Comput. Biol. 13:e1005494. doi: 10.1371/journal.pcbi.1005494
Goetsch, A. L., and Galyean, M. L. (1982). Effect of dietary concentrate level on rumen fluid dilution rate. Can. J. Anim. Sci. 62, 649–652. doi: 10.4141/cjas82-076
Hanemaaijer, M., Olivier, B. G., Roling, W. F., Bruggeman, F. J., and Teusink, B. (2017). Model-based quantification of metabolic interactions from dynamic microbial-community data. PLoS ONE 12:e0173183. doi: 10.1371/journal.pone.0173183
Hansen, E. E., Lozupone, C. A., Rey, F. E., Wu, M., Guruge, J. L., Narra, A., et al. (2011). Pan-genome of the dominant human gut-associated archaeon, Methanobrevibacter smithii, studied in twins. Proc. Natl. Acad. Sci. U.S.A. 108, 4599–4606. doi: 10.1073/pnas.1000071108
Heinken, A., Sahoo, S., Fleming, R. M., and Thiele, I. (2013). Systems-level characterization of a host-microbe metabolic symbiosis in the mammalian gut. Gut Microbes 4, 28–40. doi: 10.4161/gmic.22370
Helaszek, C. T., and White, B. A. (1991). Cellobiose uptake and metabolism by Ruminococcus flavefaciens. Appl. Environ. Microbiol. 57, 64–68.
Henderson, G., Cox, F., Ganesh, S., Jonker, A., Young, W., Global Rumen Census, C., et al. (2015). Rumen microbial community composition varies with diet and host, but a core microbiome is found across a wide geographical range. Sci. Rep. 5:14567. doi: 10.1038/srep14567
Henry, C. S., Bernstein, H. C., Weisenhorn, P., Taylor, R. C., Lee, J. Y., Zucker, J., et al. (2016). Microbial community metabolic modeling: a community data-driven network reconstruction. J. Cell Physiol. 231, 2339–2345. doi: 10.1002/jcp.25428
Henry, C. S., Dejongh, M., Best, A. A., Frybarger, P. M., Linsay, B., and Stevens, R. L. (2010). High-throughput generation, optimization and analysis of genome-scale metabolic models. Nat. Biotechnol. 28, 977–982. doi: 10.1038/nbt.1672
Howe, A., Ringus, D. L., Williams, R. J., Choo, Z. N., Greenwald, S. M., Owens, S. M., et al. (2016). Divergent responses of viral and bacterial communities in the gut microbiome to dietary disturbances in mice. ISME J. 10, 1217–1227. doi: 10.1038/ismej.2015.183
Hubisz, M. J., Pollard, K. S., and Siepel, A. (2011). PHAST and RPHAST: phylogenetic analysis with space/time models. Brief Bioinform. 12, 41–51. doi: 10.1093/bib/bbq072
Hungate, R. E. (1975). The rumen microbial ecosystem. Annu. Rev. Ecol. Syst. 6, 39–66. doi: 10.1146/annurev.es.06.110175.000351
Hungate, R. E., Dougherty, R. W., Bryant, M. P., and Cello, R. M. (1952). Microbiological and physiological changes associated with acute indigestion in sheep. Cornell Vet. 42, 423–449.
Huntington, G. B., Britton, R. A., and Prior, R. L. (1981). Feed intake, rumen fluid volume and turnover, nitrogen and mineral balance and acid-base status of wethers changed from low to high concentrate diets. J. Anim. Sci. 52, 1376–1387. doi: 10.2527/jas1981.5261376x
Hurwitz, B. L., Hallam, S. J., and Sullivan, M. B. (2013). Metabolic reprogramming by viruses in the sunlit and dark ocean. Genome Biol. 14:R123. doi: 10.1186/gb-2013-14-11-r123
Hurwitz, B. L., and U'Ren, J. M. (2016). Viral metabolic reprogramming in marine ecosystems. Curr. Opin. Microbiol. 31, 161–168. doi: 10.1016/j.mib.2016.04.002
Islam, M. M., and Saha, R. (2018). Computational approaches on stoichiometric and kinetic modeling for efficient strain design. Methods Mol. Biol. 1671, 63–82. doi: 10.1007/978-1-4939-7295-1_5
Janssen, P. H., and Kirs, M. (2008). Structure of the archaeal community of the rumen. Appl. Environ. Microbiol. 74, 3619–3625. doi: 10.1128/AEM.02812-07
Kanehisa, M., and Goto, S. (2000). KEGG: Kyoto encyclopedia of genes and genomes. Nucleic Acids Res. 28, 27–30. doi: 10.1093/nar/28.1.27
Krause, D. O., and Russell, J. B. (1996a). How many ruminal bacteria are there? J. Dairy Sci. 79, 1467–1475. doi: 10.3168/jds.S0022-0302(96)76506-2
Krause, D. O., and Russell, J. B. (1996b). An rRNA approach for assessing the role of obligate amino acid-fermenting bacteria in ruminal amino acid deamination. Appl. Environ. Microbiol. 62, 815–821.
Lenski, R. E. (1988). “Dynamics of interactions between bacteria and virulent bacteriophage,” in Advances in Microbial Ecology, ed K. C. Marshall (Boston, MA: Springer US), 1–44. doi: 10.1007/978-1-4684-5409-3_1
Levin, B. R., Stewart, F. M., and Chao, L. (1977). Resource-limited growth, competition, and predation—a model and experimental studies with bacteria and bacteriophage. Am. Nat. 111, 3–24. doi: 10.1086/283134
Li, F., and Guan, L. L. (2017). Metatranscriptomic profiling reveals linkages between the active rumen microbiome and feed efficiency in beef cattle. Appl. Environ. Microbiol. 83:e00061–17. doi: 10.1128/AEM.00061-17
Liu, K., Xu, Q., Wang, L., Wang, J., Guo, W., and Zhou, M. (2017). The impact of diet on the composition and relative abundance of rumen microbes in goat. Asian-Aust. J. Anim. Sci. 30, 531–537. doi: 10.5713/ajas.16.0353
Mahadevan, R., Edwards, J. S., and Doyle, F. J. (2002). Dynamic flux balance analysis of diauxic growth in Escherichia coli. Biophys. J. 83, 1331–1340. doi: 10.1016/S0006-3495(02)73903-9
Mahadevan, R., and Schilling, C. H. (2003). The effects of alternate optimal solutions in constraint-based genome-scale metabolic models. Metab. Eng. 5, 264–276. doi: 10.1016/j.ymben.2003.09.002
Maranas, C., and Zomorrodi, A. (2016). Optimization Methods in Metabolic Networks. Hoboken, NJ: Wiley.
Mendes-Soares, H., and Chia, N. (2017). Community metabolic modeling approaches to understanding the gut microbiome: bridging biochemistry and ecology. Free Radic. Biol. Med. 105, 102–109. doi: 10.1016/j.freeradbiomed.2016.12.017
Mendes-Soares, H., Mundy, M., Soares, L. M., and Chia, N. (2016). MMinte: an application for predicting metabolic interactions among the microbial species in a community. BMC Bioinformatics 17:343. doi: 10.1186/s12859-016-1230-3
Minot, S., Sinha, R., Chen, J., Li, H., Keilbaugh, S. A., Wu, G. D., et al. (2011). The human gut virome: inter-individual variation and dynamic response to diet. Genome Res. 21, 1616–1625. doi: 10.1101/gr.122705.111
Moss, A. R., Jouany, J.-P., and Newbold, J. (2000). Methane production by ruminants: its contribution to global warming. Ann. Zootech. 49, 231–253. doi: 10.1051/animres:2000119
Nagaraja, T. G. (2016). Microbiology of the rumen. Rumenology 2016, 39–61. doi: 10.1007/978-3-319-30533-2_2
Nocek, J. E., and Russell, J. B. (1988). Protein and energy as an integrated system. relationship of ruminal protein and carbohydrate availability to microbial synthesis and milk production. J. Dairy Sci. 71, 2070–2107. doi: 10.3168/jds.S0022-0302(88)79782-9
Oberhardt, M. A., Chavali, A. K., and Papin, J. A. (2009a). Flux balance analysis: interrogating genome-scale metabolic networks. Methods Mol. Biol. 500, 61–80. doi: 10.1007/978-1-59745-525-1_3
Oberhardt, M. A., Palsson, B. Ø., and Papin, J. A. (2009b). Applications of genome-scale metabolic reconstructions. Mol. Syst. Biol. 5, 320–320. doi: 10.1038/msb.2009.77
O'Leary, N. A., Wright, M. W., Brister, J. R., Ciufo, S., Haddad, D., Mcveigh, R., et al. (2016). Reference sequence (RefSeq) database at NCBI: current status, taxonomic expansion, and functional annotation. Nucleic Acids Res. 44, D733–D745. doi: 10.1093/nar/gkv1189
Opatowski, L., Varon, E., Dupont, C., Temime, L., Van Der Werf, S., Gutmann, L., et al. (2013). Assessing pneumococcal meningitis association with viral respiratory infections and antibiotics: insights from statistical and mathematical models. Proc. Biol. Sci. 280:20130519. doi: 10.1098/rspb.2013.0519
Orth, J. D., Thiele, I., and Palsson, B. O. (2010). What is flux balance analysis? Nat. Biotech. 28, 245–248. doi: 10.1038/nbt.1614
Owens, F. N., Secrist, D. S., Hill, W. J., and Gill, D. R. (1998). Acidosis in cattle: a review. J. Anim. Sci. 76, 275–286. doi: 10.2527/1998.761275x
Palsson, B. (2006). Systems Biology: Properties of Reconstructed Networks. Cambridge: Cambridge University Press.
Pan, D., Watson, R., Wang, D., Tan, Z. H., Snow, D. D., and Weber, K. A. (2014). Correlation between viral production and carbon mineralization under nitrate-reducing conditions in aquifer sediment. ISME J. 8, 1691–1703. doi: 10.1038/ismej.2014.38
Papin, J. A., Price, N. D., and Palsson, B. O. (2002). Extreme pathway lengths and reaction participation in genome-scale metabolic networks. Genome Res. 12, 1889–1900. doi: 10.1101/gr.327702
Parmar, N. R., Jakhesara, S. J., Mohapatra, A., and Joshi, C. G. (2016). Rumen virome: an assessment of viral communities and their functions in the rumen of an Indian buffalo. Curr. Sci. 111, 919–925. doi: 10.18520/cs/v111/i5/919-925
Pettigrew, M. M., Gent, J. F., Pyles, R. B., Miller, A. L., Nokso-Koivisto, J., and Chonmaitree, T. (2011). Viral-bacterial interactions and risk of acute otitis media complicating upper respiratory tract infection. J. Clin. Microbiol. 49, 3750–3755. doi: 10.1128/JCM.01186-11
Price, N. D., Papin, J. A., Schilling, C. H., and Palsson, B. (2003). Genome-scale microbial in silico models: the constraints-based approach. Trends Biotechnol. 21, 162–169. doi: 10.1016/S0167-7799(03)00030-1
Qiao, J. Y., Tan, Z. L., and Wang, M. (2014). Potential and existing mechanisms of enteric methane production in ruminants. Sci. Agricola 71, 430–440. doi: 10.1590/0103-9016-2013-0423
Raizada, N., Sonakya, V., Dalhoff, R., Hausner, M., and Wilderer, P. A. (2003). Population dynamics of rumen microbes using modern techniques in rumen enhanced solid incubation. Water Sci. Technol. 48, 113–119. doi: 10.2166/wst.2003.0234
Ren, Q., Chen, K., and Paulsen, I. T. (2007). TransportDB: a comprehensive database resource for cytoplasmic membrane transport systems and outer membrane channels. Nucleic Acids Res. 35, D274–D279. doi: 10.1093/nar/gkl925
Ross, E. M., Petrovski, S., Moate, P. J., and Hayes, B. J. (2013). Metagenomics of rumen bacteriophage from thirteen lactating dairy cattle. BMC Microbiol. 13:242. doi: 10.1186/1471-2180-13-242
Russell, J. B., and Rychlik, J. L. (2001). Factors that alter rumen microbial ecology. Science 292, 1119–1122. doi: 10.1126/science.1058830
Saleem, F., Ametaj, B. N., Bouatra, S., Mandal, R., Zebeli, Q., Dunn, S. M., et al. (2012). A metabolomics approach to uncover the effects of grain diets on rumen health in dairy cows. J. Dairy Sci. 95, 6606–6623. doi: 10.3168/jds.2012-5403
Saleem, F., Bouatra, S., Guo, A. C., Psychogios, N., Mandal, R., Dunn, S. M., et al. (2013). The bovine ruminal fluid metabolome. Metabolomics 9, 360–378. doi: 10.1007/s11306-012-0458-9
Salimi, F., Zhuang, K., and Mahadevan, R. (2010). Genome-scale metabolic modeling of a clostridial co-culture for consolidated bioprocessing. Biotechnol. J. 5, 726–738. doi: 10.1002/biot.201000159
Samuel, B. S., and Gordon, J. I. (2006). A humanized gnotobiotic mouse model of host-archaeal-bacterial mutualism. Proc. Natl. Acad. Sci. U.S.A. 103, 10011–10016. doi: 10.1073/pnas.0602187103
Satish Kumar, V., Dasika, M. S., and Maranas, C. D. (2007). Optimization based automated curation of metabolic reconstructions. BMC Bioinformatics 8:212. doi: 10.1186/1471-2105-8-212
Schellenberger, J., Lewis, N. E., and Palsson, B. O. (2011). Elimination of thermodynamically infeasible loops in steady-state metabolic models. Biophys. J. 100, 544–553. doi: 10.1016/j.bpj.2010.12.3707
Schellenberger, J., Park, J. O., Conrad, T. M., and Palsson, B. O. (2010). BiGG: a Biochemical Genetic and Genomic knowledgebase of large scale metabolic reconstructions. BMC Bioinformatics 11:213. doi: 10.1186/1471-2105-11-213
Seedorf, H., Kittelmann, S., and Janssen, P. H. (2015). Few highly abundant operational taxonomic units dominate within rumen methanogenic archaeal species in New Zealand sheep and cattle. Appl. Environ. Microbiol. 81, 986–995. doi: 10.1128/AEM.03018-14
Shoaie, S., Karlsson, F., Mardinoglu, A., Nookaew, I., Bordel, S., and Nielsen, J. (2013). Understanding the interactions between bacteria in the human gut through metabolic modeling. Sci. Rep. 3:2532. doi: 10.1038/srep02532
Shrestha, S., Foxman, B., Dawid, S., Aiello, A. E., Davis, B. M., Berus, J., et al. (2013). Time and dose-dependent risk of pneumococcal pneumonia following influenza: a model for within-host interaction between influenza and Streptococcus pneumoniae. J. R. Soc. Interface 10:20130233. doi: 10.1098/rsif.2013.0233
Stokes, M. R., Bull, L. S., and Halteman, W. A. (1985). Rumen liquid dilution rate in dairy cows fed once daily: effects of diet and sodium bicarbonate supplementation. J. Dairy Sci. 68, 1171–1180. doi: 10.3168/jds.S0022-0302(85)80944-9
Stolyar, S., Van Dien, S., Hillesland, K. L., Pinel, N., Lie, T. J., Leigh, J. A., et al. (2007). Metabolic modeling of a mutualistic microbial community. Mol. Syst. Biol. 3:92. doi: 10.1038/msb4100131
St-Pierre, B., and Wright, A. D. (2013). Diversity of gut methanogens in herbivorous animals. Animal 7, 49–56. doi: 10.1017/S1751731112000912
Sullivan, M. B., Lindell, D., Lee, J. A., Thompson, L. R., Bielawski, J. P., and Chisholm, S. W. (2006). Prevalence and evolution of core photosystem II genes in marine cyanobacterial viruses and their hosts. PLoS Biol. 4:e234. doi: 10.1371/journal.pbio.0040234
Suttle, C. A. (2007). Marine viruses–major players in the global ecosystem. Nat. Rev. Microbiol. 5, 801–812. doi: 10.1038/nrmicro1750
Tellier, R. C., Mathison, G. W., Okine, E. K., Mccartney, D., and Soofi-Siawash, R. (2004). Frequency of concentrate supplementation for cattle fed barley straw. 2. Ruminal dilution rates, pH and metabolite concentrations. Can. J. Anim. Sci. 84, 467–479. doi: 10.4141/A03-075
Thiele, I., and Palsson, B. O. (2010). A protocol for generating a high-quality genome-scale metabolic reconstruction. Nat. Protoc. 5, 93–121. doi: 10.1038/nprot.2009.203
Thomas, F., Hehemann, J. H., Rebuffet, E., Czjzek, M., and Michel, G. (2011). Environmental and gut bacteroidetes: the food connection. Front. Microbiol. 2:93. doi: 10.3389/fmicb.2011.00093
Thomas, M., Webb, M., Ghimire, S., Blair, A., Olson, K., Fenske, G. J., et al. (2017). Metagenomic characterization of the effect of feed additives on the gut microbiome and antibiotic resistome of feedlot cattle. Sci. Rep. 7:12257. doi: 10.1038/s41598-017-12481-6
Thompson, L. R., Zeng, Q., Kelly, L., Huang, K. H., Singer, A. U., Stubbe, J., et al. (2011). Phage auxiliary metabolic genes and the redirection of cyanobacterial host carbon metabolism. Proc. Natl. Acad. Sci. U.S.A. 108, E757–E764. doi: 10.1073/pnas.1102164108
Tymensen, L. D., and McAllister, T. A. (2012). Community structure analysis of methanogens associated with rumen protozoa reveals bias in universal archaeal primers. Appl. Environ. Microbiol. 78, 4051–4056. doi: 10.1128/AEM.07994-11
Tzamali, E., Poirazi, P., Tollis, I. G., and Reczko, M. (2011). A computational exploration of bacterial metabolic diversity identifying metabolic interactions and growth-efficient strain communities. BMC Syst. Biol. 5:167. doi: 10.1186/1752-0509-5-167
Uniprot, C. (2007). The universal protein resource (UniProt). Nucleic Acids Res. 35, D193–D197. doi: 10.1093/nar/gkl929
Varma, A., and Palsson, B. O. (1993). Metabolic capabilities of Escherichia coli: I. Synthesis of biosynthetic precursors and cofactors. J. Theor. Biol. 165, 477–502. doi: 10.1006/jtbi.1993.1202
Varma, A., and Palsson, B. O. (1994). Stoichiometric flux balance models quantitatively predict growth and metabolic by-product secretion in wild-type Escherichia coli W3110. Appl. Environ. Microbiol. 60, 3724–3731.
Wallace, R. J. (1992). Rumen microbiology, biotechnology and ruminant nutrition—the application of research findings to a complex microbial ecosystem. FEMS Microbiol. Lett. 100, 529–534. doi: 10.1111/j.1574-6968.1992.tb05751.x
Wallace, R. J., Mckain, N., Broderick, G. A., Rode, L. M., Walker, N. D., Newbold, C. J., et al. (1997). Peptidases of the rumen bacterium, Prevotella ruminicola. Anaerobe 3, 35–42. doi: 10.1006/anae.1996.0065
Wang, B., Ma, M. P., Diao, Q. Y., and Tu, Y. (2019). Saponin-induced shifts in the rumen microbiome and metabolome of young cattle. Front. Microbiol. 10:356. doi: 10.3389/fmicb.2019.00356
Weimer, P. J., Waghorn, G. C., Odt, C. L., and Mertens, D. R. (1999). Effect of diet on populations of three species of ruminal cellulolytic bacteria in lactating dairy cows. J. Dairy Sci. 82, 122–134. doi: 10.3168/jds.S0022-0302(99)75216-1
Weitz, J. S., Stock, C. A., Wilhelm, S. W., Bourouiba, L., Coleman, M. L., Buchan, A., et al. (2015). A multitrophic model to quantify the effects of marine viruses on microbial food webs and ecosystem processes. ISME J. 9, 1352–1364. doi: 10.1038/ismej.2014.220
Wolin, M. J. (1979). “The rumen fermentation: a model for microbial interactions in anaerobic ecosystems,” in Advances in Microbial Ecology. Vol. 3, ed M. Alexander (Boston, MA: Springer US), 49–77. doi: 10.1007/978-1-4615-8279-3_2
Zengler, K., and Zaramela, L. S. (2018). The social network of microorganisms - how auxotrophies shape complex communities. Nat. Rev. Microbiol. 16, 383–390. doi: 10.1038/s41579-018-0004-5
Zheng, Y., Kahnt, J., Kwon, I. H., Mackie, R. I., and Thauer, R. K. (2014). Hydrogen formation and its regulation in Ruminococcus albus: involvement of an electron-bifurcating [FeFe]-hydrogenase, of a non-electron-bifurcating [FeFe]-hydrogenase, and of a putative hydrogen-sensing [FeFe]-hydrogenase. J. Bacteriol. 196, 3840–3852. doi: 10.1128/JB.02070-14
Zhou, Y., Liang, Y., Lynch, K. H., Dennis, J. J., and Wishart, D. S. (2011). PHAST: a fast phage search tool. Nucleic Acids Res. 39, W347–W352. doi: 10.1093/nar/gkr485
Zhu, Z., Noel, S. J., Difford, G. F., Al-Soud, W. A., Brejnrod, A., Sorensen, S. J., et al. (2017). Community structure of the metabolically active rumen bacterial and archaeal communities of dairy cows over the transition period. PLoS ONE 12:e0187858. doi: 10.1371/journal.pone.0187858
Zhuang, K., Izallalen, M., Mouser, P., Richter, H., Risso, C., Mahadevan, R., et al. (2011). Genome-scale dynamic modeling of the competition between Rhodoferax and Geobacter in anoxic subsurface environments. ISME J. 5, 305–316. doi: 10.1038/ismej.2010.117
Zhuang, K., Ma, E., Lovley, D. R., and Mahadevan, R. (2012). The design of long-term effective uranium bioremediation strategy using a community metabolic model. Biotechnol. Bioeng. 109, 2475–2483. doi: 10.1002/bit.24528
Zomorrodi, A. R., Islam, M. M., and Maranas, C. D. (2014). d-OptCom: dynamic multi-level and multi-objective metabolic modeling of microbial communities. ACS Synth. Biol. 3, 247–257. doi: 10.1021/sb4001307
Zomorrodi, A. R., and Maranas, C. D. (2012). OptCom: a multi-level optimization framework for the metabolic modeling and analysis of microbial communities. PLoS Comput. Biol. 8:e1002363. doi: 10.1371/journal.pcbi.1002363
Zomorrodi, A. R., and Segre, D. (2017). Genome-driven evolutionary game theory helps understand the rise of metabolic interdependencies in microbial communities. Nat. Commun. 8:1563. doi: 10.1038/s41467-017-01407-5
Zomorrodi, A. R., Suthers, P. F., Ranganathan, S., and Maranas, C. D. (2012). Mathematical optimization applications in metabolic networks. Metab. Eng. 14, 672–686. doi: 10.1016/j.ymben.2012.09.005
Keywords: microbial community, microbiome-virome interaction, rumen, genome-scale metabolic modeling, viral auxiliary metabolic genes
Citation: Islam MM, Fernando SC and Saha R (2019) Metabolic Modeling Elucidates the Transactions in the Rumen Microbiome and the Shifts Upon Virome Interactions. Front. Microbiol. 10:2412. doi: 10.3389/fmicb.2019.02412
Received: 17 May 2019; Accepted: 07 October 2019;
Published: 22 October 2019.
Edited by:
Steve Lindemann, Purdue University, United StatesReviewed by:
Joshua Chan, Colorado State University, United StatesAmit Ghosh, Indian Institute of Technology Kharagpur, India
Sangram Bagh, Saha Institute of Nuclear Physics (SINP), India
Copyright © 2019 Islam, Fernando and Saha. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Rajib Saha, rsaha2@unl.edu