- Department of Chemical and Bioprocess Engineering, School of Engineering, Pontificia Universidad Católica de Chile, Santiago, Chile
The gut microbiome is a complex microbial community that has a significant influence on the host. Microbial interactions in the gut are mediated by dietary substrates, especially complex polysaccharides. In this environment, breakdown products from larger carbohydrates and short chain fatty acids are commonly shared among gut microbes. Understanding the forces that guide microbiome development and composition is important to determine its role in health and in the intervention of the gut microbiome as a therapeutic tool. Recently, modeling approaches such as genome-scale models and time-series analyses have been useful to predict microbial interactions. In this study, a bottom-up approach was followed to develop a mathematical model based on microbial growth equations that incorporate metabolic sharing and inhibition. The model was developed using experimental in vitro data from a system comprising four microorganisms of the infant gut microbiome (Bifidobacterium longum subsp. infantis, Lactobacillus acidophilus, Escherichia coli, and Bacteroides vulgatus), one substrate (fructooligosaccharides, FOS), and evaluating two metabolic products (acetate and lactate). After parameter optimization, the model accurately predicted bacterial abundance in co-cultures from mono-culture data. In addition, a good correlation was observed between the experimental data with predicted FOS consumption and acid production. B. infantis and L. acidophilus were dominant under these conditions. Further model validation included cultures with the four-species in a bioreactor using FOS. The model was able to predict the predominance of the two aforementioned species, as well as depletion of acetate and lactate. Finally, the model was tested for parameter identifiability and sensitivity. These results suggest that variations in microbial abundance and activities in the infant gut were mainly explained by metabolic interactions, and could be properly modeled using Monod kinetics with metabolic interactions. The model could be scaled to include data from larger consortia, or be applied to microbial communities where sharing metabolic resources is important in shaping bacterial abundance. Moreover, the model could be useful in designing microbial consortia with desired properties such as higher acid production.
Introduction
The human colonic microbiome is a complex microbial community that has a significant impact on host health. This is a diverse community that reaches high cell densities and includes four dominant phyla (Bacteroidetes, Firmicutes, Actinobacteria and Proteobacteria) (Lozupone et al., 2012; Qin et al., 2013). The gut microbiome coexists with the host and deploys important functions that impact host metabolism and gut physiology (Rajilić-Stojanović and de Vos, 2014). Even though its composition is variable among people (Goodrich et al., 2014), the functions these microorganisms performed are basically conserved (Lozupone et al., 2012; Qin et al., 2013). In certain cases, imbalances in the composition of the microbiome are a contributing factor to the onset of inflammatory bowel diseases as well as autoimmune and metabolic (Greenblum et al., 2012; Sevelsted et al., 2015; Marchesi et al., 2016; Tamburini et al., 2016; Cox et al., 2017). How the microbiome assembles in the first months of life appears to be important later in life (Tamburini et al., 2016). One important factor shaping the early microbiome is the type of feeding (Qin et al., 2013; Rajilić-Stojanović and de Vos, 2014). Human breast milk contains large amounts of oligosaccharides (HMO), which are selectively utilized by beneficial gut microbes. Bifidobacterium species such as B. longum subsp. infantis display multiple adaptations to utilize these substrates (Thomson et al., 2017). Lactobacilli are also abundant in the infant gut microbiome (Bäckhed et al., 2015). In contrast, formula-fed infants have a distinct microbiome composition, not dominated by Bifidobacterium and with a higher representation of members of Bacteroides (B. fragilis, B. vulgatus) and Enterobacteria (Escherichia coli, Klebsiella spp.) (Bäckhed et al., 2015). The activity of these microbes results in high amounts of acetate and lactate in infant feces, resulting in an acidic pH (Cinquin et al., 2004; Tamburini et al., 2016).
Microbial interactions are important for the assembly and functioning of the gut microbiome. Dominant ecological interactions found in the gut microbiome are competition and cooperation (Faust and Raes, 2012). These interactions broadly represent the sum of all physical, chemical and microbiological activities that microorganisms exert upon others (Roume et al., 2015; Vogt et al., 2015; Hecht et al., 2016; Rakoff-Nahoum et al., 2016). Considering that diet is a major driver guiding gut microbiome composition, microbial interactions are influenced by dietary compounds (Cameron et al., 2014; Medina et al., 2017; Tuncil et al., 2017). Cross-feeding of fermentation breakdown products of the microbiome appears to be common among gut species (Rogowski et al., 2015). This has been shown for example in the utilization of mucin and sialylated milk oligosaccharides between B. bifidum and B. breve (Egan et al., 2014a,b), or during fructan consumption between bifidobacteria and butyrate- producing bacteria (Moens et al., 2016). Cross-feeding is also observed when metabolic end products from one microorganism, such as amino acids or short chain fatty acids (SCFA), are used by another microorganism (Egan et al., 2014a; Moens et al., 2016). For example, lactate and acetate are end products of lactic acid bacteria, which could be utilized by butyrate-producing bacteria such as Faecalibacterium prausnitzii and Eubacterium rectale (Louis and Flint, 2017).
Modeling-based approaches have been recently developed to study and predict the composition and interactions in the gut microbiome (Magnúsdóttir et al., 2016). These include ecological-statistic models, genome-scale metabolic reconstructions (GSM) and ordinary differential equation (ODE)-based kinetic models (Trosvik et al., 2010a; Kettle et al., 2015). A Generalized Additive Model (GAM) (Hastie and Tibshirani, 1990) consists of a statistic regression technique that has been used in time-series analysis of ecological data to characterize and estimate cross-feeding and competition between microorganisms. GAMs do not need any assumption about functional relationships in the group for its formulation. However, they could be affected by overfitting when many parameters are needed for matching the data (Wood, 2008; Trosvik et al., 2010b). GAMs usually require a post cross-validation process to curate the model (Ward, 2014). After proper calibration and validation, these models provide accurate predictions by interpolation (Trosvik et al., 2010b).
Lately, GSMs have been successfully applied to explore microbial interactions among gut microbes (Magnúsdóttir et al., 2016). They require an extensive database for reconstruction, editing and gap-filling of full metabolic pathways (Thiele et al., 2014). Several techniques based on orthology, topology and stoichiometry of biological reactions facilitate the draft design and curation process (Thiele and Palsson, 2010). Characteristic features of the species to be reconstructed must be first identified (Kanehisa, 2006). After curation and defining specific environments and constraints, microbial interactions can be obtained for a few species (Thiele et al., 2013).
Recently, a kinetic model constructed from experimental data of gut microbes in a bioreactor was presented, aimed to model the dynamic behavior of the gut microbiome (Kettle et al., 2015). The analysis required a metabolic pathway input and a matrix describing the compounds produced during the fermentation, to generate an ODE system for simulation of microbiome abundance (Walker et al., 2011). Here, microbiome complexity was simplified assigning gut microbes to ten bacterial functional groups (BFGs), based on metabolic properties such as similar breakdown of complex substrates or similar SCFA production or consumption patterns (Kettle et al., 2015). The model showed a good fit with experimental data, which corresponded to a continuous flow bioreactor inoculated with human fecal microbiota.
In order to help understanding the forces dominating gut microbiome structure and composition, here we developed and assessed a mathematical model based on microbial growth equations, taking into account metabolic interactions among bacteria. We focused on the interactions of four gut microbes, Bifidobacterium longum subsp. infantis, Lactobacillus acidophilus, Bacteroides vulgatus and Escherichia coli, during their growth in vitro using fructooligosaccharides (FOS) as substrate. FOS is a well studied prebiotic with degree of polymerization of fructose of 3–6 units (Roberfroid et al., 2010). Experimental data was obtained from co-culture experiments, which were used later to construct and calibrate the model, including the impact of metabolic inhibition or stimulation on bacterial growth. The model was finally validated using additional experimental data of the consortium of the four species on FOS using a biological reactor.
Materials and Methods
Microorganisms and Media
Microorganisms used in this study were obtained from the UC Davis, Department of Viticulture and Enology Culture Collection (L. acidophilus ATCC 4356, B. infantis ATCC 15697, Escherichia coli K12), and the American Type Culture Collection (Bacteroides vulgatus ATCC 8482; Manassas, VA, United States). Bacteria were, respectively, cultured at 37°C for 24 h in de Man–Rogosa–Sharp (MRS), MRS supplemented with 0.05% L-cysteine-HCl (Loba Chemie, India), LB broth, or Reinforced Clostridium Medium (Becton-Dickinson) supplemented with 1 g/L L-cysteine. All bacteria excepting E. coli were routinely grown under anaerobic conditions in an anaerobic jar (Anaerocult, Merck, Germany) with anaerobic packs (Gaspak EM, Becton Dickinson). All media were pre-reduced in an anaerobic jar overnight before inoculation, and prior to each assay bacteria were sub-cultured twice.
Co-culture Batch Experiments
Combinations of L. acidophilus (La), E. coli (Ec), B. vulgatus (Bv) and B. infantis (Bi) were prepared in co-culture experiments. Culture media used was a modified version of previously described ZMB (Zhang et al., 2009), which was supplemented with hemin (0.01 g/L, Sigma–Aldrich, St. Louis, MO, United States) and L-cysteine-HCl (0.5 g/L, Sigma–Aldrich, St. Louis, MO, United States). Single amino acid groups in ZMB were replaced by Bacto-Tryptone (at 28 g/L). Carbon sources used were either lactose (10 g/L; Lyngby, Denmark) or FOS (10 g/L; Raftilose Synergy 1, Orafti, Malvern, PA, United States) as carbon source. Single cultures of B. infantis (Bi), B. vulgatus (Bv), E. coli (Ec) and L. acidophilus (La); and co-cultures BiBv, BiEc, BiLa, BvEc, BvLa and EcLa were prepared. An experiment with all bacteria (All) and a negative control with no bacteria were included. Fresh overnight cultures of each microorganism were washed in sterile mZMB, and 1 mL of each overnight culture was used to inoculate 10 mL of mZMB containing FOS. This experiment was performed in duplicate. Volumes of 200 μL of inoculated mZMB were placed in 96 well sterile microplates, covered with 30 μL of sterile mineral oil, and incubated in anaerobic jars at 37°C for either 24, 48, or 72 h. In parallel, growth was monitored every 12 h in a microplate reader (Tecan Infinite M200 PRO, Switzerland). Samples were recovered from each microplate and centrifuged at 12000 ×g for 2 min. Pellets and supernatants were stored at -20°C until use.
Quantification of Bacterial Abundance by qPCR
Total DNA from each sample was purified using the UltraClean® Microbial DNA Isolation Kit (Mo Bio Laboratories, Carlsbad, CA, United States), following manufacturer instructions and using a Disruptor Genie (Scientific Industries, Inc., Bohemia, NY, United States). Extracted DNA was quantified using a NanoQuant Plate in the Tecan Infinite M200 PRO plate reader, and diluted to 1 ng/μL to be used in qPCR reactions. For qPCR we used 0.2 μM of the following primers: for Bv, Bacteroidetes primer F (5′-GGTGTCGGCTTAAGTGCCAT-3′) and Bacteroidetes primer R (5′-CGGACGTAAGGGCCGTGC-3′); for Bi, Blon_0883F (5′-AGTTCGGCTCCAAAGACCTG-3′) and Blon_0883R (5′-CATGCCTCGATACGGTCGAA), targeting an ABC solute binding protein; for Ec, Eco1457F (5′-CATTGACGTTACCCGCAGAAGAAG) and Eco1652R (5′-CTCTACGAGACTCAAGCTTGC-3′) (Kassinen et al., 2004); and for La, LACTO_F (5′-TGGAAACAGRTGCTAATACCG-3′) and LACTO_R (5′-GTCCATTGTGGAAGATTCCC-3′) (Bartosch et al., 2004). qPCR reactions were performed using the qPCR PowerUp SYBR Green Master Mix in MicroAmp Fast Optical plates (Applied Biosystems, United States), and using a StepOnePlus Real-Time PCR System (Applied Biosystems, United States). Reactions were carried out for 2 min at 50°C, 2 min at 95°C and 40 cycles of 3 s at 95°C and 30 s at 62°C. Absolute quantification was performed including a standard curve using DNA from a pure culture of each species, with dilutions starting from 1 ng/μL to 0.1 pg/μL. To convert bacterial DNA concentrations into cell genome numbers, the following equation was used (equation 1).
Batch Bioreactor Culturing
Four independent batch co-culture experiments were performed in a 250 mL bioreactor (Mini-bio Applikon Biotechnology, Netherlands), using mZMB as culture media supplemented with FOS at 1%. In these experiments, the four microorganisms (Bi-La-Ec-Bv) were inoculated at an initial OD630 of 0.05. The bioreactor has two six-bladed Rushton turbines and operated at 100 rpm. The temperature was set at 37°C and the pH was maintained at 5.5 with automatic injection of 3N HCl and 3N NaOH. The dissolved oxygen concentration was set at 1 ppm by purging N2 (99.99% grade) before inoculation and during the lag phase. The foam level was controlled adding 100 μL antifoam in the inoculum (Polydimethylsiloxane base, Winkler, Chile). Two milliliter from the bioreactor were obtained every 2 h and centrifuged at 4000 ×g for 5 min. Supernatants were stored at -20°C for carbohydrate and SCFA quantification. Pellets were stored for DNA extraction, quantified and diluted to 10 ng/μL for qPCR assays as described above in an AriaMx Realtime PCR System (Agilent Technologies, Santa Clara, CA, United States).
Sample Analysis
Total carbohydrate quantification was performed using the phenol-sulfuric acid method (Tuomivaara et al., 2015). Acetate and lactate were quantified by HPLC using an Aminex HPX-87H ion exchange carbohydrate-organic acid column (Bio-Rad, United States) at 35°C with a flow rate of 0.450 mL/min (H2SO4 5 mM, mobile phase) on a LaChrom L-700 HPLC system (Hitachi, Japan), equipped with a Diode Array and a Refractive Index detectors as described previously (Mendoza et al., 2017).
Model Development
The equations used in the model are described in the Model development section in Supplementary Material. The model, the parameter identifiability and sensitivity analysis codes are also presented in Supplementary Material. As input for the determination of the parameters, mono-culture and paired co-culture abundance data are required, in addition to an estimation of acetate and lactate produced and carbohydrate consumed under these conditions. To simplify the analysis, some assumptions were taken into account: (a) an inhibition term was added to Monod kinetics (Model development, Supplementary Material); (b) a microorganism will prefer the consumption of the main carbon sources (glucose, lactose), over other intermediates produced during the fermentation; (c) the ability of a microorganism to produce or consume an intermediate was determined from its metabolic pathway and the literature, and later confirmed experimentally in mono-cultures.
Results
Model Description
In this work a kinetic black-box model was developed, aimed to predict the abundance of a bacterial population, substrate consumption and SCFA production, based on mono and co-culture data (Figure 1). The model is based on microbial growth equations, but it also considers the metabolic influence of one microorganism on another. This could be considered as a feedback control mechanism (Figure 1).
FIGURE 1. Model general representation. Initial substrate and product concentrations and lag phase are used as input (black bars). Microbial growth, consumption, and acid production are considered to interact with other bacteria. Final outputs are observed substrate, acids, and biomass.
Parameter Settings in Mono-culture
For single microorganisms, the general model consisted of 5 ODEs (Equations 2, 4, 5, and 6 in Supplementary Material), 17 parameters and constitutive Monod-like inhibition equations (Sacher et al., 2011). Mono-culture parameters (Table 1) were set as described in the Parameter fitting section in the Supplementary Material. 96 well-plates mono-cultures of Bi, Bv, Ec and La were prepared, in a semi-synthetic media (mZMB) and using FOS as the sole carbon source. Bacterial abundance, FOS consumption and acetate and lactate produced were measured to fit model parameters. An average of eight parameters were set for each bacterium (Table 1), which were found by the optimization task. The calculated error in the assay is shown in Supplementary Table S1. For any microorganism and under all conditions, parameter Ks (half-velocity constant) appeared insensitive.
TABLE 1. Parameters found via scatter search in mono-culture and then used in co-culture optimization.
Paired Co-culture and Parameter Fitting
The model was later expanded to include the metabolic interaction between two microorganisms. This model consists of 7 ODEs, 17 parameters per bacteria and two interaction parameters per co-culture. Every parameter not calibrated in mono-culture was set in this step. In order to fit the co-culture parameters, all paired combinations of microorganisms were cultured in FOS and analyzed as described above. Figure 2A shows the percentage of change in abundance for all six paired combinations, determined experimentally. As a comparison, Figure 2B shows these percentage changes according to the fitted models. Most of the times, the model was able to predict well the changes in abundance in all co-cultures. Experimentally, initial Ec cell numbers were higher than the other microorganisms. However, during growth Bi and La recovered in part their levels compared to Ec (Figures 2A,B). Co-culture data allowed the prediction of Bv predominance over La and Bi during growth on FOS, which was also observed experimentally.
FIGURE 2. Changes in bacterial population during growth on FOS, expressed as percentage of the co-culture in time. (A) experimental data of co-cultures; (B) model estimation of abundance in co-cultures; (C) abundance of the four-species co-culture in microplates, experimental (Left) and estimated by the model (Right); (D) abundance in the four-species co-culture in the bioreactor during growth on FOS, experimental (Left) and estimated by the model (Right).
Figures 3A,B compares the experimental consumption of FOS by the co-cultures with the values simulated with the fitted model. Most experimental and simulated combinations showed total carbohydrate depletion between 24 and 48 h. In general the model indicated a faster consumption compared to experimental data. One important exception was the BvLa paired co-culture, in which not all of the carbohydrate was consumed. This behavior was not captured by the model, which assumed that since both bacteria reached 100% consumption in single culture, the same rule should apply to their combination.
FIGURE 3. Heat map representing FOS concentration in co-cultures and its prediction by the model. (A) FOS concentration in paired co-cultures. (B) Model estimation of FOS concentration in paired co-cultures; (C) experimental and predicted FOS concentration in the four-species co-culture in microplates; (D) experimental and predicted FOS concentration in the four-species co-culture in the bioreactor.
Figure 4A shows the concentration of acetate produced over time. In certain cases the model predicted the experimental behavior of acetate production. Bi combinations displayed larger acetate amounts compared to other co-cultures, and in certain cases the model predicted higher values than what was observed. Interestingly, the model predicted that acetate production in co-culture BiEc will have a peak and later decrease. This was also observed experimentally, but at a different time and different intensity (Figure 4A). These results indicate that Bi growth is an important parameter for sensitivity assays.
FIGURE 4. Acetate production and estimation by the model. (A) Experimental data in paired co-cultures, compared to values predicted by the model; (B) experimental and predicted acetate values of the consortium in microplate assay; (C) experimental and predicted acetate values of the consortium in the bioreactor.
Figure 5A displays the concentration of lactate in co-cultures. A good agreement between observed and predicted data was obtained in co-cultures BiLa, BiEc and LaEc. Combinations BiBv and BvLa were predicted to produce lactate because of Bi and La activities; however, lactate amounts were negligible and not reproduced well by the model. In addition, BvEc co-culture showed production of lactate, but the model assumptions and structure did not consider this situation. The error calculated (equation 10 in Supplementary Material) for the parameter fitting process is shown in Supplementary Table S1.
FIGURE 5. Lactate production and estimation by the model. (A) Experimental data in paired co-cultures, compared to values predicted by the model; (B) experimental and predicted lactate values of the consortium in microplate assay; (C) experimental and predicted lactate values of the consortium in the bioreactor.
The parameters determined in paired co-cultures are shown in Table 1. The interaction parameters in Table 2 indicate the influence of one microorganism on another’s growth rate. A negative value indicates that one microorganism favors another’s growth, while a positive term indicates inhibition. Values near 0 suggest a greater interaction effect, while values near the limit indicate there is no effect on the other bacteria. A strong inhibition was found from Ec to Bv, and in general the effects observed were positive or neutral.
TABLE 2. Interaction parameters (efji in equation 8, Supplementary Material) found in co-cultures by the model.
Model Validation Using Bacterial Consortia
Finally, the model was validated using independent experimental data from co-culture of the four microorganisms using FOS as the sole carbon source. The experiment was set in microplates and analyzed as discussed above. To test the validity of the model in another set-up, the consortium was additionally cultured on FOS in a 250 mL pH/oxygen controlled stirred bioreactor. This batch system offers a much more controlled and reproducible anaerobic environment, which also provides much faster growth compared to microplates.
Figure 2C shows percentage abundance data obtained for each member of the consortium in microplate assays. The initial levels of Bv were much lower compared to the other three microorganisms. Interestingly, the amounts of La, Ec and Bi in the well-plates cultures were closely predicted by the model. Under these conditions, Bi dominated the co-culture using FOS, followed by La. A good prediction was also observed for the total carbohydrate concentration in spent media (Figure 3C). Finally, the amounts of acetate and lactate appeared overestimated by the model (Figures 4B, 5B).
As expected, growth of the consortium in the bioreactor resolved in a shorter time compared to the assays above (Figure 2D). Therefore, time was linearly adjusted for comparison and integration in the model. As in microplates, we observed a predominance of Bi and La. This observation was sustained during the course of the fermentation. Interestingly, the model was also able to predict this predominance (Figure 2D). In addition, both the model and data showed a full consumption of FOS at 12 h (Figure 3D). Finally, a good agreement of acetate and lactate amounts between the experimental evidence and the model was obtained (Figures 4C, 5C). Since La was a good competitor during growth on FOS in the bioreactor, lactate concentrations appeared higher compared to previous experiments (Figures 5A–C). The parameters that define the production of lactate and acetate in Bi appear to be important in the four-bacterium co-culture, considering the predominance of Bi.
Finally, we performed a simple additional simulation to test the prediction capabilities of the model where a bacteriostatic agent is used against each member of the consortium (Figure 6). In every co-culture where Bi was able to grow, it predominated over the others (Figures 6B–D). On the other hand, if Bi was inhibited, Ec predominated in the co-culture (Figure 6A).
FIGURE 6. Simulation of the effect of a bacteriostatic agents on the consortium. These agents are simulated to be directed and inhibit the growth of each member of the consortium. (A) Bifidobacterium infantis is unable to grow; (B) Bacteroides vulgatus is unable to grow; (C) Escherichia coli is unable to grow; (D) Lactobacillus acidophilus is unable to grow.
Parameter Identifiability Analysis
Parameter identifiability was used to find correlations between parameters (Parameter identifiability in Supplementary Material). This analysis is important for further reducing the number of fitted parameters by setting one of them and defining the other as a function. Inspection of the parameter covariance matrix is one way to find which parameters allow the model to be identifiable. As shown in Figure 7, highlighted cells display a high correlation (positive or negative). Usually parameters inside a cluster have a high correlation. In this case, this could be observed for all parameters from the same microorganism. For example, production of acetate and lactate in Bi are directly correlated, while some correlations between microorganisms were found. La’s parameters (Ysx – biomass yield, μmax – Maximum growth rate, Ia – Acetate inhibition constant, Il – Lactate inhibition constant) are inversely correlated to Ec bacterial parameters such as growth and inhibition constants. This suggests that the higher the La growth, the lower the E. coli biomass yield and higher inhibition. Several parameters associated to Bv growth were mostly directly correlated to Ec growth, indicating a more neutral or cooperative interaction.
FIGURE 7. Parameter model identifiability. Correlation values between each parameter in the model was calculated (for each microorganism including interaction). Only >|0.95| values are highlighted; red values are inversely correlated, while blue values are directly correlated. Parameters on both axes are indicated in Table 1.
Parameter Sensitivity Analysis
This analysis allows the determination of the influence of every parameter in each differential equation of the model. As shown in Figure 8, the effects of the parameters initially set are important in every ODE, due to the fact that Bi appears as the dominant microorganism in the consortium (Figures 2C,D). Specifically, the second parameter of the model (Bi’s μmax) has the highest influence on every other microorganism and their metabolic equations. Parameters K3 and K4 (Bi’s inhibition constants of acetate and lactate) also display a large influence on other microorganisms. In order to analyze the effects of the sensitive parameters found in the previous assay, Figure 9 shows the average and standard deviation after 5000 iterations of randomly changing a parameter by 5% in its amount. The strongest effect of changing the value of Bi’s μmax is on Ec cell numbers (Figure 9A), variable that can vary around 4% the value. On the other hand, a change in a parameter could also imply an advance or delay in the kinetics. Figure 9B shows the change in the FOS consumption kinetics due to effects of higher or lower values of Bi’s lactate inhibition constant. Here we observed that changing the parameter only altered the dynamics of the ODE. Finally, Figure 9C shows the last case found in the sensitivity analysis, a parameter that is not sensitive to any differential equation. For example, measured Bv was not affected even after changing 50% parameter 25 (Ec substrate yield Ysx).
FIGURE 8. Model parameter average sensitivity. Sensitivity (Y-axis) of each parameter (X-axis) for every ODE is shown. xi represents every ODE described in the model (x1 = dS/dt, x2 = dA/dt, x3 = dL/dt, x4 = dX1/dt, x5 = dX1m/dt, x6 = dX2/dt, x7 = dX2m/dt, x8 = dX3/dt, x9 = dX3m/dt, x10 = dX4/dt, x11 = dX4m/dt), where S: substrate; A: acetate; L: lactate. X: live biomass; Xm: total biomass. Parameters are in the same order in Figure 6.
FIGURE 9. Variation of the ODEs values (g/L) over time due a 5% change in the parameters in 5000 iterations. (A) worst case scenario, with parameter 2 (Bi’s μmax) affecting measured Ec ODE; (B) a change in kinetics scenario, with parameter 4 (Bi’s lactate inhibition constant), affecting FOS consumption; (C) a non-sensitive parameter (measured Bv ODE), for example to Ec substrate yield.
Discussion
The gut microbiome is a complex microbial community that modulates several host responses. This connection to host health makes it important to understand what forces guide microbiome composition and cause it to drift to an altered or dysbiotic microbiome (Cox et al., 2017). The interest in determining and predicting key factors in the establishment and maintenance of the gut microbiome is the major goal of several works (Trosvik et al., 2010a; Greenblum et al., 2012; Kettle et al., 2015; Shashkova et al., 2016).
Diet is a major modulator of the composition of the gut microbiome, and the nature of these substrates probably dictates which species predominate. In this study we evaluated if a mathematical model capturing metabolic interactions is able to recapitulate the composition and functions of a consortium of species of the gut microbiome. For this, we chose four representative bacteria of the infant gut microbiome, and using experimental data from mono and co-culture, a model was developed, calibrated and validated. Using a bioreactor, the developed model was assessed in a more controlled environment.
The system was studied during growth on FOS, a major prebiotic present in infant formula (Roberfroid et al., 2010). All members of the consortium display the ability to use this substrate (Roberfroid et al., 2010), including E. coli which could use small amounts of mono or disaccharides found in FOS. Moreover, different molecular mechanisms for FOS consumption have been described (Barrangou et al., 2003). In general the predictions by the model followed the in vitro behavior of the consortium, either in paired co-cultures, and growing the four-species consortium either in microplates or in a more controlled environment such as a biological reactor. This indicates that the model is able to predict changes in the bacterial abundance using only co-culture data for calibration.
It is very possible that interactions and parameters determined in this study are dependent on which prebiotic is used. FOS are commonly added to infant formula, but in combination with galactooligosaccharides (GOS), another important prebiotic (Garrido et al., 2013). Breast milk contains large concentrations of HMO, which are also a large catalog of oligosaccharides derived from lactose (Thomson et al., 2017). Moreover, the gut epithelium is covered with a mucin layer, containing oligosaccharides that could be used as carbon source by infant gut bacteria (Tailford et al., 2015). In a more realistic situation probably all these carbohydrates contribute to shape microbial interactions in different ways, since their chemical structure selects for specific microbial strains endowed with the cognate molecular machinery for utilization. However, if metabolic interactions are key in shaping microbiome composition, we could hypothesize that a mathematical model including these interactions could predict microbiome composition when other substrates are used.
We observed a good fit between experimental data and modeling results. This suggests that inhibitions observed in certain cases could be due to acetate and lactate production, variables that were quantified and included in the mathematical model. Both the reactor and the microplates had an initial pH of 5.5, however, pH was not regulated in the latter system. Considering this, similar results in both systems could also indicate that results obtained are independent of the pH.
A general good agreement was also observed for acid production and carbohydrate consumption. For Bi in mono-cultures and co-cultures where it predominates, the amounts of acetate and lactate produced are near a 3:2 ratio (Garrido et al., 2013). This was also observed during the growth of the consortium in the bioreactor. Acetate production by Ec was overestimated by the model (0.21 g of acetate per 1 g of FOS consumed). In general Ec was thought to benefit from other microorganism activities in that it uses mono or disaccharides released to the media (Table 2) (Ravcheev et al., 2013; Vuoristo et al., 2015). Another possibility might be protein fermentation by Ec (Lulit and Strohl, 1990). Lactate production of La determined by the model was around 0.63 g per 1 g of substrate, a similar yield in lactose reported (Fu and Mathews, 1999).
In some co-cultures, the concentration of either acetate or lactate was overestimated. This was evident in Bi co-cultures whenever it predominated. Parameters of the model could be much better estimated in experiments with improved resolution and more frequent measurements. Since the model in co-cultures defines the intervals where the parameters are most sensitive, it is possible that an increase in the number of samples would reduce the variation of underestimated parameters. The time points where the substrate is being fully consumed are critical, and microorganisms could find another substrate for growth (for secondary fermenters) or entering to a stationary phase. Also, for Bacteroides and Escherichia cultures, the microbial concentration could be overestimated by some intrinsic pathways of these genera (Neis et al., 2015; Vuoristo et al., 2015).
Moreover, while acetate and lactate are major metabolic products in this system, a more complete picture could be obtained if the model included other metabolites. Adding more equations of utilization and inhibition by metabolites such as ethanol, propionate, butyrate and amino acids could be important. Amino acid cross-feeding between Bacteroides and Lactobacillus supports bacterial growth in vitro and in silico (Magnúsdóttir et al., 2016).
The analysis of bacteriostatic agent effects on the culture suggested that Bi should be predominant if other bacteria are inhibited. However, when Bi is inhibited, La or Bv should grow more than Ec, because of their glycolytic properties (Ravcheev et al., 2013). This is a limitation of the model, probably due to missing functions that describe the breakdown of complex carbohydrates by Bv, or the protein fermentation as a carbon source of bacteria. In addition, further work could corroborate these hypotheses by adding the respective antibiotic and measuring the same variables used in this work.
A possible application of this initial ODE-based model is that it could be used to predict microbial composition in the gut based on diet, at least in simpler microbiome communities. This work indicates that it is possible to have a good approach to this goal if metabolic interactions are included. Moreover, bacterial composition of a microbiome could eventually be optimized, for example to increase production of acetate and lactate. These two acids are important modulators of health outcomes in the gut. For example acetate has been shown to prevent pathogen colonization (Fukuda et al., 2011), and lactate in the adult gut microbiome is used by butyrate-producing bacteria (Moens et al., 2016), a health-promoting SCFA (Louis and Flint, 2017).
Finally, this model could be useful to study interactions using a more complex set of species of gut microbiome species. In general these results could be important to predict the composition of microbial communities where metabolic interactions are relevant. Considering the flexibility of incorporating product equations and growth inhibitions to the model, this model could be used to find microbial consortia with desired metabolic properties such as maximized acid production.
Author Contributions
FP and DM performed experiments. FP developed the model. DG provided materials and reagents. JP-C and DG conceived the experiments, wrote, and edited the manuscript.
Funding
This work was funded by Fondecyt de Iniciacion 11130518, Fondecyt de Postdoctorado 3160525 and FONDEF IDEA ID16i10045.
Conflict of Interest Statement
Theauthors 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
We thank Dr. Eduardo Agosin for the support with equipment and methods, Paulina Torres for supporting in bioreactor cultures. We thank Dr. Martin Gotteland and Dr. Patricia García Cañete for providing materials for this study. We also thank Dr. David Mills and Lucy Joseph for access to V&E Culture Collection.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2017.02507/full#supplementary-material
References
Bäckhed, F., Roswall, J., Peng, Y., Feng, Q., Jia, H., Kovatcheva-Datchary, P., et al. (2015). Dynamics and stabilization of the human gut microbiome during the first year of life. Cell Host Microbe 17, 690–703. doi: 10.1016/j.chom.2015.04.004
Barrangou, R., Altermann, E., Hutkins, R., Cano, R., and Klaenhammer, T. R. (2003). Functional and comparative genomic analyses of an operon involved in fructooligosaccharide utilization by Lactobacillus acidophilus. Proc. Natl. Acad. Sci. U.S.A. 100, 8957–8962. doi: 10.1073/pnas.1332765100
Bartosch, S., Fite, A., Macfarlane, G. T., and McMurdo, M. E. T. (2004). Characterization of bacterial communities in feces from healthy elderly volunteers and hospitalized elderly patients by using real-time PCR and effects of antibiotic treatment on the fecal microbiota. Appl. Environ. Microbiol. 70, 3575–3581. doi: 10.1128/AEM.70.6.3575-3581.2004
Cameron, E. A., Kwiatkowski, K. J., Lee, B.-H., Hamaker, B. R., Koropatkin, N. M., and Martens, E. C. (2014). Multifunctional nutrient-binding proteins adapt human symbiotic bacteria for glycan competition in the gut by separately promoting enhanced sensing and catalysis. mBio 5:e01441-14. doi: 10.1128/mBio.01441-14
Cinquin, C., Le Blay, G., Fliss, I., and Lacroix, C. (2004). Immobilization of infant fecal microbiota and utilization in an in vitro colonic fermentation model. Microb. Ecol. 48, 128–138. doi: 10.1007/s00248-003-2022-7
Cox, L. M., Yamanishi, S., Sohn, J., Alekseyenko, A. V., Leung, J. M., Cho, I., et al. (2017). Altering the intestinal microbiota during a critical developmental window has lasting metabolic consequences. Cell 158, 705–721. doi: 10.1016/j.cell.2014.05.052
Egan, M., Motherway, M. O. C., Ventura, M., and van Sinderen, D. (2014a). Metabolism of sialic acid by Bifidobacterium breve UCC2003. Appl. Environ. Microbiol. 80, 4414–4426. doi: 10.1128/AEM.01114-14
Egan, M., O’Connell Motherway, M., Kilcoyne, M., Kane, M., Joshi, L., Ventura, M., et al. (2014b). Cross-feeding by Bifidobacterium breve UCC2003 during co-cultivation with Bifidobacterium bifidum PRL2010 in a mucin-based medium. BMC Microbiol. 14:282. doi: 10.1186/s12866-014-0282-7
Faust, K., and Raes, J. (2012). Microbial interactions: from networks to models. Nat. Rev. Microbiol. 10, 538–550. doi: 10.1038/nrmicro2832
Fu, W., and Mathews, A. P. (1999). Lactic acid production from lactose by Lactobacillus plantarum: kinetic model and effects of pH, substrate, and oxygen. Biochem. Eng. J. 3, 163–170. doi: 10.1016/S1369-703X(99)00014-5
Fukuda, S., Toh, H., Hase, K., Oshima, K., Nakanishi, Y., Yoshimura, K., et al. (2011). Bifidobacteria can protect from enteropathogenic infection through production of acetate. Nature 469, 543–547. doi: 10.1038/nature09646
Garrido, D., Ruiz-moyano, S., Jimenez-espinoza, R., and Eom, H. (2013). Utilization of galactooligosaccharides by Bifidobacterium longum subsp. infantis isolates. Food Microbiol. 33, 262–270. doi: 10.1016/j.fm.2012.10.003
Goodrich, J. K., Waters, J. L., Poole, A. C., Sutter, J. L., Koren, O., Blekhman, R., et al. (2014). Human genetics shape the gut microbiome. Cell 159, 789–799. doi: 10.1016/j.cell.2014.09.053
Greenblum, S., Turnbaugh, P. J., and Borenstein, E. (2012). Metagenomic systems biology of the human gut microbiome reveals topological shifts associated with obesity and inflammatory bowel disease. Proc. Natl. Acad. Sci. U.S.A. 109, 594–599. doi: 10.1073/pnas.1116053109
Hastie, T., and Tibshirani, R. (1990). Generalized additive models. Stat. Sci. 10, 354–363. doi: 10.2307/2246134
Hecht, A. L., Casterline, B. W., Earley, Z. M., Goo, Y. A., Goodlett, D. R., Bubeck Wardenburg, J., et al. (2016). Strain competition restricts colonization of an enteric pathogen and prevents colitis. EMBO Rep. 94:e201642282. doi: 10.15252/embr.201642282
Kanehisa, M. (2006). From genomics to chemical genomics: new developments in KEGG. Nucleic Acids Res. 34, D354–D357. doi: 10.1093/nar/gkj102
Kassinen, A., Malinen, E., Krogius, L., Palva, A., and Rinttila, T. (2004). Development of an extensive set of 16S rDNA-targeted primers for quantification of pathogenic and indigenous bacteria in faecal samples by real-time PCR. J. Appl. Microbiol. 97, 1166–1177. doi: 10.1111/j.1365-2672.2004.02409.x
Kettle, H., Louis, P., Holtrop, G., Duncan, S. H., and Flint, H. J. (2015). Modelling the emergent dynamics and major metabolites of the human colonic microbiota. Environ. Microbiol. 17, 1615–1630. doi: 10.1111/1462-2920.12599
Louis, P., and Flint, H. J. (2017). Formation of propionate and butyrate by the human colonic microbiota. Environ. Microbiol. 19, 29–41. doi: 10.1111/1462-2920.13589
Lozupone, C. A., Stombaugh, J. I., Gordon, J. I., Jansson, J. K., and Knight, R. (2012). Diversity, stability and resilience of the human gut microbiota. Nature 489, 220–230. doi: 10.1038/nature11550
Lulit, G. W., and Strohl, W. R. (1990). Comparison of growth, acetate production, and acetate inhibition of Escherichia coli strains in batch and fed-batch fermentations. Appl. Environ. Microbiol. 56, 1004–1011.
Magnúsdóttir, S., Heinken, A., Kutt, L., Ravcheev, D. A., Bauer, E., Noronha, A., et al. (2016). Generation of genome-scale metabolic reconstructions for 773 members of the human gut microbiota. Nat. Biotechnol. 35, 81–89. doi: 10.1038/nbt.3703
Marchesi, J. R., Adams, D. H., Fava, F., Hermes, G. D. A., Hirschfield, G. M., Hold, G., et al. (2016). The gut microbiota and host health: a new clinical frontier. Gut 65, 330–339. doi: 10.1136/gutjnl-2015-309990
Medina, D., Pinto, F., Ovalle, A., Thomson, P., and Garrido, D. (2017). Prebiotics mediate microbial interactions in a consortium of the infant gut microbiome. Int. J. Mol. Sci. 18:E2095. doi: 10.3390/ijms18102095
Mendoza, S. N., Cañón, P. M., Contreras, A., Ribbeck, M., and Agosin, E. (2017). Genome-scale reconstruction of the metabolic network in Oenococcus oeni to assess wine malolactic fermentation. Front. Microbiol. 8:534. doi: 10.3389/FMICB.2017.00534
Moens, F., Weckx, S., and De Vuyst, L. (2016). Bifidobacterial inulin-type fructan degradation capacity determines cross-feeding interactions between bifidobacteria and Faecalibacterium prausnitzii. Int. J. Food Microbiol. 231, 76–85. doi: 10.1016/j.ijfoodmicro.2016.05.015
Neis, E. P. J. G., Dejong, C. H. C., and Rensen, S. S. (2015). The role of microbial amino acid metabolism in host metabolism. Nutrients 7, 2930–2946. doi: 10.3390/nu7042930
Qin, J., Li, R., Raes, J., Arumugam, M., Burgdorf, S., Manichanh, C., et al. (2013). A human gut microbial gene catalog established by metagenomic sequencing. Nature 464, 59–65. doi: 10.1038/nature08821.A
Rajilić-Stojanović, M., and de Vos, W. M. (2014). The first 1000 cultured species of the human gastrointestinal microbiota. FEMS Microbiol. Rev. 38, 996–1047. doi: 10.1111/1574-6976.12075
Rakoff-Nahoum, S., Foster, K. R., and Comstock, L. E. (2016). The evolution of cooperation within the gut microbiota. Nature 533, 255–259. doi: 10.1038/nature17626
Ravcheev, D. A., Godzik, A., Osterman, A. L., and Rodionov, D. A. (2013). Polysaccharides utilization in human gut bacterium Bacteroides thetaiotaomicron: comparative genomics reconstruction of metabolic and regulatory networks. BMC Genomics 14:873. doi: 10.1186/1471-2164-14-873
Roberfroid, M., Gibson, G. R., Hoyles, L., McCartney, A. L., Rastall, R., Rowland, I., et al. (2010). Prebiotic effects: metabolic and health benefits. Br. J. Nutr. 104, S1–S63. doi: 10.1017/S0007114510003363
Rogowski, A., Briggs, J. A., Mortimer, J. C., Tryfona, T., Terrapon, N., Lowe, E. C., et al. (2015). Glycan complexity dictates microbial resource allocation in the large intestine. Nat. Commun. 6:7481. doi: 10.1038/ncomms8481
Roume, H., Heintz-Buschart, A., Muller, E. E. L., May, P., Satagopam, V. P., Laczny, C. C., et al. (2015). Comparative integrated omics: identification of key functionalities in microbial community-wide metabolic networks. NPJ Biofilms Microbiomes 1:15007. doi: 10.1038/npjbiofilms.2015.7
Sacher, J., Saa, P., Cárcamo, M., López, J., Gelmi, C. A., and Pérez-Correa, R. (2011). Improved calibration of a solid substrate fermentation model. Electron. J. Biotechnol. 14, 7. doi: 10.2225/vol14-issue5-fulltext-7
Sevelsted, A., Stokholm, J., Bonnelykke, K., and Bisgaard, H. (2015). Cesarean section and chronic immune disorders. Pediatrics 135, e92–e98. doi: 10.1542/peds.2014-0596
Shashkova, T., Popenko, A., Tyakht, A., Peskov, K., Kosinsky, Y., Bogolubsky, L., et al. (2016). Agent based modeling of human gut microbiome interactions and perturbations. PLOS ONE 11:e0148386. doi: 10.1371/journal.pone.0148386
Tailford, L. E., Crost, E. H., Kavanaugh, D., and Juge, N. (2015). Mucin glycan foraging in the human gut microbiome. Front. Genet. 6:81. doi: 10.3389/fgene.2015.00081
Tamburini, S., Shen, N., Wu, H. C., and Clemente, J. C. (2016). The microbiome in early life: implications for health outcomes. Nat. Med. 22, 713–722. doi: 10.1038/nm.4142
Thiele, I., Heinken, A., and Fleming, R. M. T. (2013). A systems biology approach to studying the role of microbes in human health. Curr. Opin. Biotechnol. 24, 4–12. doi: 10.1016/j.copbio.2012.10.001
Thiele, I., and Palsson, B. (2010). A protocol for generating a high-quality genome-scale metabolic reconstruction. Nat. Protoc. 5, 93–121. doi: 10.1038/nprot.2009.203.A
Thiele, I., Vlassis, N., and Fleming, R. M. T. (2014). FASTGAPFILL: efficient gap filling in metabolic networks. Bioinformatics 30, 2529–2531. doi: 10.1093/bioinformatics/btu321
Thomson, P., Medina, D. A., and Garrido, D. (2017). Human milk oligosaccharides and infant gut bifidobacteria: molecular strategies for their utilization. Food Microbiol. (in press). doi: 10.1016/j.fm.2017.09.001
Trosvik, P., Rudi, K., Strætkvern, K. O., Jakobsen, K. S., Næs, T., and Stenseth, N. C. (2010a). Web of ecological interactions in an experimental gut microbiota. Environ. Microbiol. 12, 2677–2687. doi: 10.1111/j.1462-2920.2010.02236.x
Trosvik, P., Stenseth, N. C., and Rudi, K. (2010b). Convergent temporal dynamics of the human infant gut microbiota. ISME J. 4, 151–158. doi: 10.1038/ismej.2009.96
Tuncil, Y. E., Xiao, Y., Porter, N. T., Reuhs, B. L., Martens, E. C., and Hamaker, B. R. (2017). Reciprocal prioritization to dietary glycans by gut bacteria in a competitive environment promotes stable coexistence. mBio 8:e01068-17. doi: 10.1128/mBio.01068-17
Tuomivaara, S. T., Yaoi, K., O’Neill, M. A., and York, W. S. (2015). Generation and structural validation of a library of diverse xyloglucan-derived oligosaccharides, including an update on xyloglucan nomenclature. Carbohydr. Res. 402, 56–66. doi: 10.1016/j.carres.2014.06.031
Vogt, S. L., Peña-Díaz, J., and Finlay, B. B. (2015). Chemical communication in the gut: effects of microbiota-generated metabolites on gastrointestinal bacterial pathogens. Anaerobe 34, 106–115. doi: 10.1016/j.anaerobe.2015.05.002
Vuoristo, K. S., Mars, A. E., Sangra, J. V., Springer, J., Eggink, G., and Sanders, J. P. M. (2015). Metabolic engineering of the mixed-acid fermentation pathway of Escherichia coli for anaerobic production of glutamate and itaconate. AMB Express 5:61. doi: 10.1186/s13568-015-0147-y
Walker, A. W., Ince, J., Duncan, S. H., Webster, L. M., Holtrop, G., Ze, X., et al. (2011). Dominant and diet-responsive groups of bacteria within the human colonic microbiota. ISME J. 5, 220–230. doi: 10.1038/ismej.2010.118
Ward, T. (2014). The Information Theoretically Efficient Model (ITEM): a model for computerized analysis of large datasets. arXiv:1409.6075
Wood, S. N. (2008). Fast stable direct fitting and smoothness selection for generalized additive models. J. R. Stat. Soc. Ser. B Stat. Methodol. 70, 495–518. doi: 10.1111/j.1467-9868.2007.00646.x
Keywords: metabolic interaction, gut microbiome diet, prebiotics, mathematical modeling, fructooligosaccharides (FOS)
Citation: Pinto F, Medina DA, Pérez-Correa JR and Garrido D (2017) Modeling Metabolic Interactions in a Consortium of the Infant Gut Microbiome. Front. Microbiol. 8:2507. doi: 10.3389/fmicb.2017.02507
Received: 07 September 2017; Accepted: 01 December 2017;
Published: 14 December 2017.
Edited by:
Michele Guindani, University of California, Irvine, United StatesReviewed by:
Angela M. Zivkovic, University of California, Davis, United StatesEmanuele Bosi, University of Florence, Italy
Copyright © 2017 Pinto, Medina, Pérez-Correa and Garrido. 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) or licensor 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: Daniel Garrido, ZGdhcnJpZG9jQGluZy5wdWMuY2w=