Abstract
We applied the so-called chemical kinetics approach to complex bacterial growth patterns that were dependent on the liquid-surface-area-to-volume ratio (SA/V) of the bacterial cultures. The kinetic modeling was based on current experimental knowledge in terms of autocatalytic bacterial growth, its inhibition by the metabolite CO2, and the relief of inhibition through the physical escape of the inhibitor. The model quantitatively reproduces kinetic data of SA/V-dependent bacterial growth and can discriminate between differences in the growth dynamics of enteropathogenic E. coli, E. coli JM83, and Salmonella typhimurium on one hand and Vibrio cholerae on the other hand. Furthermore, the data fitting procedures allowed predictions about the velocities of the involved key processes and the potential behavior in an open-flow bacterial chemostat, revealing an oscillatory approach to the stationary states.
1. Introduction
Kinetic curves deliver the dynamic fingerprint of time-evolving systems. This applies to the study of complex chemical systems [1], has been used to portray ecological population dynamics [2], and is perhaps true for bacterial systems despite their tremendous microscopic complexity. In this last sense, one may ask how much information about the underlying mechanism of bacterial growth dynamics could be extracted by evaluating bacterial growth curves. This might resemble the attempt to decipher the mechanism of a chemical system by studying its kinetic behavior.
In their attempt to predict bacterial growth kinetics in food substrates, Doona and coworkers [3–5] introduced the so-called chemical kinetics approach (CKA) to rationalize bacterial growth curves. Basically, the approach encompassed the treatment of bacteria, nutrients, and metabolites as chemical species involved in a number of coupled processes that obeyed the rules of chemical kinetics. The time-evolution of the various components was expressed by a set of coupled ordinary differential equations that were numerically resolved. Since bacterial growth is essentially nonlinear, the incorporation of positive or negative feedback allowed an integral, “all-in-one” description of evolving bacterial growth including the growth stages usually denominated as “lag,” “exponential,” “stationary,” and “death” phases [6–11] with just one model and without changing the parameters or their values. Hence the CKA relates to the exploration of the global and continuous bacterial growth dynamics and perhaps allows predictions as well as the identification and discrimination of complex mechanisms that could be at the origin of atypical growth behavior.
From the above viewpoint, the CKA is settled in between empirical laws and their extensions usually based on logistic curve fitting procedures [7–11] and the kinetic modeling at the biochemical level [12–14]. While the first approach allows only very limited insight into the mechanism behind bacterial growth dynamics, the second approach usually features a molecularly inclined view on selected bacterial processes but not on the overall growth dynamics. Thus the CKA is a mechanism-based simulation procedure that uses a coarse-grained dynamical viewpoint.
Peleg and Corradini [15] recently not only discussed the benefits of the CKA but also pointed out that the use of fixed order kinetics, as well as fixed stoichiometries in the CKA, could be a potential disadvantage in dealing with macroscopic and complex elements such as bacteria. This is readily understood by considering, for instance, the interaction between one bacterium and one nutrient molecule leading to some change in the bacteria, that is, the chemical equivalent of the conversion of the bacteria into a product. This would be the case of second-order kinetics. However, it is more likely that one bacterium interacts with many nutrient molecules to convert them into a “product.” This consideration can lead to a substantial increase in the formal kinetic order of such processes and because of the huge discrepancy in the particular orders, that is, very low order in the bacteria and very high order in the nutrient, they must not necessarily remain constant during the whole time-evolution. Hence a further evaluation of the CKA is pertinent so that it may be applicable to more complex growth patterns than treated so far.
As one example of less usual bacterial growth patterns, among others [16, 17], we previously found that the rate of later stage bacterial growth of enteropathogenic E. coli (EPEC) was correlated to the liquid-surface-area-to-volume ratio (SA/V) of the own liquid bacterial cultures [18], that is, to the scaling properties of the bacterial culturing system in dispersion (Figure 1).
The SA/V effect was related to the generation of a volatile inhibitor metabolite, namely, CO2 [19], which was produced by the bacteria during and after their conversion into an anaerobic regime [20]. Such conversion occurs because of the progressive depletion of the dissolved oxygen in the medium. The rate by which CO2 escapes from the liquid cultures depends on the SA/V of the medium and its physical escape gives rise to a relief of inhibition.
Hence we reasoned that the relief of inhibition accounted for the SA/V-dependent growth rates. Moreover, equivalent to the curve shape of typical zero-order kinetics, the last growth stage of EPEC was characterized by a quasilinear increase of the bacterial concentration in time. The suggested scenario at the origin of such growth behavior was supported by complementary experiments as well as by a simplified kinetic model, which was built in the spirit of the CKA. Such kinetic CKA-inspired model reproduced qualitatively the growth curves in terms of an interplay between inhibition (generation of the inhibitor) and the relief of inhibition (escape of the inhibitor).
In this paper, we show that our former observations of SA/V-dependent EPEC growth as well as the present recordings of E. coli JM83 and Salmonella typhimurium (also Salmonella enterica subsp. enterica ser. Typhimurium) growth can be quantitatively reproduced by the same model based on the CKA. This was in contrast to the Vibrio cholerae growth patterns that seemingly follow a distinct growth dynamics. Hence our modeling indicates not only that the SA/V-effect, as well as the particular kinetic curve shapes, is caused by a common mechanism for the different bacteria but also that changes may occur when the global growth dynamics differ.
To our knowledge, we demonstrate for the first time that the CKA combined with data fitting can discriminate between the growth dynamics of different bacteria simply by analyzing the shape of the kinetic growth curves. Moreover, the ranges of the extracted kinetic parameters allowed us to estimate the velocities of the key steps of bacterial autocatalytic growth, growth inhibition, and its relief by the escape of the inhibitor. The applied procedure also permitted predictions of the possible growth behavior in a bacterial chemostat, shedding more light on the dynamic properties of SA/V-dependent bacterial growth.
2. Material and Methods
Model calculations were performed with the simulation-adjustment software package Sa3 written by D. Lavabre, Laboratoire des IMRCP, Université Paul Sabatier, Toulouse, France, that is freely available at http://cinet.chim.pagesperso-orange.fr/index.html#. The general algorithm for the numerical integration of the differential equations was based on a semi-implicit fourth-order Runge-Kutta method with stepwise control for stiff ordinary differential equations. The minimization algorithms for the adjustment of the kinetic parameter values were of either the Powell type or random walk simulated annealing. Fitted rate parameters were automatically and iteratively returned to the numerical integration until a minimum in the residual error was reached.
The bacterial strains EPEC, E. coli JM83, Salmonella typhimurium, and Vibrio cholerae were grown under the same standard procedures with fixed initial bacterial concentrations in nutrient-rich liquid media (batch cultures) as previously described [18, 19]. The optical density (OD) was measured in a UV-VIS spectrophotometer cuvette at 600 nm and at 37°C ± 0.3°C under gentle magnetic stirring. The scaling properties were established by varying the height of the liquid column in the spectrophotometric cuvette or by obstructing the liquid/gas surface area (SA/V = 0 cm−1) with a premolded paraffin block that tightly fitted the inside of the cuvette. The total bacterial concentration was assumed directly proportional to the OD.
3. Results and Discussion
3.1. Model Description
Our chemical kinetics model (Scheme 1) [18] was extended to mimic flow conditions and the continuous input of the CO2 inhibitor. As a novelty, the numerical execution includes the possibility to adjust the kinetic parameters with the aim of fitting experimentally observed data, which is essential to perform meaningful predictions as well as comparisons between the growth kinetics of different bacteria.
 |
The kinetic scheme consists of two coupled autocatalytic processes, one autoinhibitory step and the relief of inhibition, all together giving rise to considerable dynamic complexity and possible nonintuitive behavior. Step (I) denotes the autocatalytic growth of the aerobic bacteria (B) limited by O2. It accounts for the common “lag” and “exponential” phase. Process (II) stands for the conversion of aerobic to anaerobic bacteria (Ban). Step (III) is the autocatalytic growth of anaerobic bacteria limited by the nutrient (N) and inhibited by CO2. It represents the SA/V-dependent and quasilinear growth stage. Step (IV) stands for the escape of CO2 from the system giving rise to a relief of inhibition and causing the dependence on the SA/V.
Since all processes are coupled, the kinetic scheme gathers together the different growth stages over the whole time-evolution into one single model. The corresponding set of reaction fluxes and differential equations are given in Table 1.
3.2. Initial Conditions and Parameter Estimation
The apparent 1 : 1 stoichiometry between B and O2 in step (I) of the kinetic scheme is hardly realistic with respect to a quantitative data treatment. It is more likely that many O2 molecules interact with one bacterium to cause an event such as the indicated self-replication. For that reason, we have used a normalized initial concentration of = 4.0×10−8 M giving rise to a unity reaction order in [O2] of model step (I) (please see the Appendix for further details).
In the model, the escape rate of CO2 from the system was considered independent of the other components of the mixture. A semiquantitative reproduction of former measurements of the CO2 escape at different SA/V [17] delivers estimates for the values of the first-order rate constant () to be in the range between 10−5 and 10−4 s−1 for SA/V ratios between 0.14 and 2.5 cm−1.
3.3. Data Fitting of SA/V-Dependent E. coli Growth
The adjustment of the kinetic parameters to fit the experimentally observed growth curves was started manually with the “lag” and “exponential” growth stage and then extended to the entire kinetic curves. Subsequently, the automatic fitting procedures were carried out.
Starting with the first attempts, we noticed a very good agreement (flexibility) between the model and the experimental data leading to excellent fitting results with respect to each individual experiment, with residual errors as small as (Figure 2(a)). In these attempts, all kinetic parameters for the individual fitting of each of the kinetic curves were allowed to vary until the minimum residual error between simulated and experimental data was reached. Table 2 shows the obtained parameter variations and the parameter mean values for the 6 fitted kinetic curves.
(a)
(b)
The data given in Table 2 indicates relatively small variations of the parameter values although each experiment was fitted independently. Note that the automatic fitting procedure has “chosen” parameter values for that are in good agreement with the experimental data of the CO2 escape that we have measured for different SA/V ratios.
An exception to the small parameter variations is represented by the control parameter that was implemented to tune the inhibition strength. In four experiments this parameter remained at a value of around 108 M−1, but, in the experiments performed at SA/V = 0.8 cm−1, the order of magnitude was 6 and in the case of SA/V = 0.4 cm−1 it was 9. A less pronounced situation occurred for in which the order of magnitude is usually −1 except at SA/V = 0.57 cm−1 where its value was −2 and at SA/V = 0.8 cm−1 where its value was −3. These variations could be the consequence of a parameter coupling between and and the further relationship between and in the rate expression .
Table 2 gives insight into the relative velocities of the assumed key processes for the SA/V-dependent EPEC growth. For instance, the autocatalytic growth of Ban (Ban + N + CO2) is predicted to occur significantly slower than that of B (B + B), / ≈ 0.02. Furthermore, the process appears to be the slowest step in the model that apparently determines the overall rate.
Although the parameter variations were found to be rather small during the individual fitting procedures, a more challenging and realistic task for the CKA should be to keep the parameters , , , and at fixed values for the 6 experiments and only to allow variations in in order to reflect the differences in the SA/V ratios of the experimental series, thus assuming that all the other parameters are invariant regardless of changes in the scaling properties of the system.
The results of such attempt where exclusively has been adjusted during the automatic fitting procedure and all the other parameters were kept at fixed values, that is, the same , , , and for all 6 experiments, are shown in Figure 2(b). This data fitting was employed by selecting arbitrarily one experiment (SA/V = 0.67 cm−1) to start with, keeping at a reasonable order of magnitude in accordance with the experimental observations of the CO2 escape, performing the adjustment of , , , and , and then using the adjusted parameters for SA/.67 cm−1 for the other experiments by exclusively adjusting .
As anticipated by the previous fitting results, we obtained an excellent fitting of the SA/V = 0.67 cm−1 growth curve. But also for the other experiments between SA/V = 0.5 and SA/V = 1.0 cm−1 a very reasonable data reproduction was obtained. In contrast, the fitting attempt of the SA/V = 0.4 cm−1 curve with the most horizontal SA/V-dependent growth phase failed. This could be due to experimental data scattering as well as hidden parameters or may reflect the need to include additional model steps. Nevertheless, the fitting of five individual experiments by only changing can be considered as satisfactory.
Figure 3 shows the fitting of two experiments performed under significantly extended time and with a completely obstructed surface (SA/V = 0 cm−1) as well as an entirely free surface (SA/V = 0.71 cm−1) of the bacterial cultures by using culture volumes of 1.4 cm3 in both cases. Although the fitted parameters are not identical to those given in Figure 2, since dealing with two independent experimental series, they all lie very close to the observed parameter ranges of Table 2.
These results reveal that the model can respond properly to the two extreme culture surface conditions of “open” versus “closed” causing a marked difference in the later growth rates. Interestingly, the curve shape of the SA/V = 0 cm−1 case resembles the typical S-shaped portrays of bacterial growth kinetics [6], which was reproduced in the present case by keeping at zero. Such situation also occurs when the culture surface is not completely obstructed but the SA/V becomes sufficiently small. In fact, this is the case for one standard in vitro culturing condition, that is, the growth in a static tube filled with culture medium. Hence, based on our observations, the later stage growth usually appears stationary (zero growth rate) but only due to a small SA/V, that is, because of the very slow escape of the inhibiting CO2.
3.4. Data Fitting of SA/V-Dependent E. coli JM83, Salmonella typhimurium, and Vibrio cholerae Growth
Our current experiments showed that the SA/V growth effect is not restricted to EPEC. In fact, Vibrio cholerae, E. coli JM83, and Salmonella typhimurium, at least, also respond to changes in the scaling properties during the later stage growth phases. However, a simple eye inspection of the shape of these growth curves shown in Figure 4 indicates differences with respect to the formerly recorded EPEC growth curves. Specifically, the later growth stage in the other bacterial species is curved while that of EPEC appears linear.
(a)
(b)
(c)
We subjected this experimental data to our modeling approach. Again, for the fitting of each bacterial system, only was changed to account for the different SA/V. Since all systems were different also the extracted kinetic parameter values differ between EPEC, Vibrio cholerae, and E. coli JM83. Table 3 provides an overview of the fitting results in comparison to the previous EPEC fitting. In general, the obtained values for E. coli JM83 are closer to EPEC than those of Salmonella typhimurium.
As indicated in Figure 4, a reasonably good experimental data fitting could be obtained for E. coli JM83 and Salmonella typhimurium while an acceptable fitting for Vibrio cholerae was impossible.
The fitting results of EPEC, E. coli JM83, and Salmonella typhimurium versus Vibrio cholerae are particularly interesting in two aspects:(1)from the viewpoint of bacterial phylogeny there is a difference between EPEC, E. coli JM83, and Salmonella typhimurium on one side and Vibrio cholerae on the other side [21]. We assume that these phylogenetic differences have an impact on the mechanism responsible for the growth kinetics. This implies that our model can discriminate between EPEC “like” and “unlike” growth dynamics as given by the circumstance that—based on a fixed model structure—in one case a reasonable set of parameter values can be obtained while in the case of Vibrio cholerae it cannot be obtained;(2)our model, despite its apparent flexibility, cannot fit any number of bacterial growth curves, indicating a certain specificity that, to our opinion, only the CKA can reveal. Such specificity was already observed during the model design, for instance, when we tested different degrees of autocatalysis that ended up as quadratic in the proposed model. The potential behavior of a cubic autocatalytic process in model step (I), 2B instead of B, resulted in the impossibility of fitting any of the experimental EPEC data, even for the individual fitting attempts.
Hence the model in its present structure not only proves valid because of its capacity to fit various experimentally observed data but also appears specific by apparently rejecting the fitting of EPEC “unlike” growth dynamics.
3.5. Prediction of Chemostat Behavior
We have evaluated the behavior of our model under hypothetical flow conditions with the continuous inflow of the reactants B and O2 and the corresponding outflow of the reaction mixture simulating a continuously fed and well-mixed bacterial chemostat [22]. The parameter values and initial conditions were taken from the fitting of the experiment SA/V = 0.71 cm−1 in Figure 3 as a prototypical case. As shown in Figure 5, some interesting transient behavior was observed before the concentrations reached stationary values. At specific flow rates, all concentrations are approaching their stationary states in a damped oscillatory manner. Note that this predicted behavior lasts very long.
Sustained oscillations or kinetic bistability by systematic variations of the flowrate was not observed or identified in the given parameter space. However, the simulations indicate complex dynamic behavior under flow conditions that could be of importance for practical purposes when EPEC-like organisms are grown in a bacterial chemostat. Moreover, the oscillatory approach to stationary bacterial concentrations could have an impact on any open bacterial growth system, such as the in vivo duplication of pathogenic bacteria and their expression of virulence.
The flow term extension of the model equations was also used to simulate the continuous input of the inhibitory CO2 into the semiopen batch system. As expected, a significant decrease in the rate of the later stage growth phase occurs. In contrast, an initial oversaturation of the cultures with CO2 under batch conditions is predicted to have no effect on the simulated bacterial growth curves due to the escape of CO2 before the SA/V-dependent growth phase was reached.
4. Conclusions
In this paper, we have applied the chemical kinetics approach (CKA) to rationalize and reproduce still unrecognized bacterial growth that was dependent on the liquid-surface-area-to-volume ratio (SA/V) of the bacterial cultures. Like in its previous applications, for instance, in the fitting of Staphylococcus aureus growth kinetics [3], the CKA was confirmed as a useful tool to understand the bacterial growth dynamics from a kinetic point of view as well as to suggest predictions of complex behavior under open-flow conditions.
As a novel feature, we found that the CKA combined with data fitting can serve to discriminate between different growth dynamics as well as to recognize similarities in growth mechanisms. This was demonstrated by the grouping of EPEC, E coli JM83, and Salmonella typhimurium on the one hand and Vibrio cholerae on the other hand. Furthermore, we have shown that the CKA can also apply to more complex growth dynamics with the possible extraction of rate parameter values or their relations that may be used to assess the relevance of the key processes that give rise to the curve shapes of the bacterial time-evolution.
It is important to underscore that the applied fitting procedures were entirely different from classical curve fitting because in the present case we performed an adjustment of the kinetic parameters of a fixed set of differential equations and did not use arbitrary equations as in usual curve fitting approach. The fitting proceeded by an automatic error minimization between the simulated and the experimental data. Usually, if a model is “wrong” no satisfactory parameter adjustment can be achieved. This case could apply to the Vibrio cholerae system. Hence a successful fitting, that is, a small residual error between experimental and simulated data like in the cases of EPEC, E coli JM83, and Salmonella typhimurium, indicates that the chosen model was a reasonable proposal.
Upcoming studies will be devoted to further experimental insight into the SA/V effect such as the examination of corresponding bacterial growth under open-flow conditions and a finer-grained model outline based on additional biochemical as well as kinetic information.
Appendix
According to experimental observations [18], the initial concentration of dissolved O2 in the culture medium of a typical growth experiment drops to nearly zero when the “exponential growth phase” ends and the “quasilinear,” SA/V-dependent, growth stage starts. In terms of a 1 : 1 stoichiometry between B and O2, the bacterial concentration at this point of O2 depletion should be equal to the initial concentration of the limiting O2. We estimate the initial O2 media concentration roughly as 10−4 M based on Henry’s law [23]. Given that the OD value for the transitions between the “exponential” and the “quasilinear” stage for all experiments is OD ≈ 0.5, where [O]2 ≈ 0 M, and taking the available information given in Table 4, the ratio [O2]0/[B] = 1.0×10−4/6.25×10−13 = 1.6×108 is associated with the apparent kinetic order in O2 of model step (I), that is, with the number of O2 molecules that are required to drive one bacterium to a change.
To avoid numerical difficulties by using a very small initial [B] and a very high kinetic order in O2, the initial concentration of B was normalized by using the concentration equivalence [.] between the bacteria and consumed O2 at OD = 0.5 extrapolated to [] at OD = 0. Hence for the simulations the initial bacterial concentration was normalized to 4.0 × 10−8 M giving rise to a unity reaction order in [O2] of model step (I).
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
Joaquín Sánchez gratefully acknowledges support by CONACYT, Mexico, Project CB-2010-01/156257. The authors also thank Dr. Jean-Claude Micheau for stimulating discussions.