Bioinformatics Advance Access originally published online on August 25, 2007
Bioinformatics 2007 23(20):2754-2759; doi:10.1093/bioinformatics/btm400
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Metal reduction kinetics in Shewanella
1BACTER Institute and 2Department of Mathematics and Department of Biochemistry, University of Wisconsin – Madison,Wisconsin 53706, USA
*To whom correspondence should be addressed.
 |
ABSTRACT |
|---|
 TOP ABSTRACT 1 INTRODUCTION2 METHODS3 RESULTS4 DISCUSSIONACKNOWLEDGEMENTSREFERENCES |
|---|
Motivation: Metal reduction kinetics have been studied in cultures of dissimilatory metal reducing bacteria which include the Shewanella oneidensis strain MR-1. Estimation of system parameters from time-series data faces obstructions in the implementation depending on the choice of the mathematical model that captures the observed dynamics. The modeling of metal reduction is often based on Michaelis–Menten equations. These models are often developed using initial in vitro reaction rates and seldom match with in vivo reduction profiles.
Results: For metal reduction studies, we propose a model that is based on the power law representation that is effectively applied to the kinetics of metal reduction. The method yields reasonable parameter estimates and is illustrated with the analysis of time-series data that describes the dynamics of metal reduction in S.oneidensis strain MR-1. In addition, mixed metal studies involving the reduction of Uranyl (U(VI)) to the relatively insoluble tetravalent form (U(IV)) by S.alga strain (BR-Y) were studied in the presence of environmentally relevant iron hydrous oxides. For mixed metals, parameter estimation and curve fitting are accomplished with a generalized least squares formulation that handles systems of ordinary differential equations and is implemented in Matlab. It consists of an optimization algorithm (Levenberg–Marquardt, LSQCURVEFIT) and a numerical ODE solver. Simulation with the estimated parameters indicates that the model captures the experimental data quite well. The model uses the estimated parameters to predict the reduction rates of metals and mixed metals at varying concentrations.
Contact: rl8q{at}cms.mail.virginia.edu
Supplementary information: Supplementary data are available at Bioinformatics online.
|
1 INTRODUCTION |
|---|
|
TOPABSTRACT1 INTRODUCTION2 METHODS3 RESULTS4 DISCUSSIONACKNOWLEDGEMENTSREFERENCES |
|---|
The biotransformation of wastes containing uranium, radionuclides and heavy metals is cost-effective as compared to traditional methodologies (Lovley and Phillips, 1992; Macaskie et al., 1997). Dissimilatory metal reducing bacteria (DMRB) can reduce various metals and radionuclides, including sediment-abundant Fe(III), Mn(III/IV) and aqueous species of U(VI), Cr(VI), Co(III) and Tc(VII) (Gorby and Lovley, 1992; Gorby et al., 1998; Lloyd and Macaskie, 1996; Lovley, 1993; Lovley et al., 1991; Nealson and Saffarini, 1994; Roden and Zachara, 1996; Wildung et al., 2000). The kinetics of metal reduction has previously been formulated as Monod models based on time-dependent data from bacteria (Gaudy and Gaudy, 1980; Monod, 1949; Rittmann and McCarty, 2001). These time-courses contain an important amount of information about the structure and regulation of the underlying biological system. The task of extracting the information involves the use of a suitable modeling framework that captures suitably the dynamics of the data. In addition to linear models, there exists an infinite variety of non-linear structures that are potential candidates for optimal data representation. Biology does not have a rich repertoire of first-principle laws and the best-suited functions for describing higher-level biological phenomena are simply unknown. As a result, biology often resorts to approximations that are locally anchored in mathematical theory and whose global appropriateness is supported by experience. As the range of approximations is virtually unlimited, modeling involves criteria such as data fit, interpretability and mathematical tractability rather than deep theory. The descriptive models are often non-linear and formulated as a system of differential equations where the optimization of parameters is more complex than linear programming where a fast and unique solution is often guaranteed. Non-linear methods such as regression, genetic algorithms, simulated annealing or dynamic programming often lead to challenging convergence issues (Voit et al., 2006).
Power-law representations are a convenient mathematical construct to model the behavior of biological systems. Here, we discuss a method to efficiently model metal reduction pathways that can be linearized. Many parts of biochemical systems are indeed linear (Stanbury et al., 1983) and thus by using power law representations for the underlying model, linear segments of pathways can be estimated with efficient linear regression techniques (Lall and Voit, 2005). In case of mixed metals, the power law model is solved through a non-linear regression. In our specific case, the number of parameters is not very high and the non-linear regression (using lsqcurvefit in Matlab) found a good solution.
|
2 METHODS |
|---|
|
TOPABSTRACT1 INTRODUCTION2 METHODS3 RESULTS4 DISCUSSIONACKNOWLEDGEMENTSREFERENCES |
|---|
Models of biochemical systems are often based on non-linear ordinary differential equations (ODE'). Given availability of time-series data for system variables of interest, it should be possible to use a non-linear search algorithm, which solves the equations at each iteration and ultimately yields optimal parameter values. However, such an approach is not feasible for large problems and alternative techniques must be explored. As a result, the decoupling procedure is applicable. This procedure uses the slopes as estimates of the true differentials on the left-hand side of the ODE model. As a result, the non-linear ODE is reformulated from one involving n non-linear ODE' to a larger system of n x N algebraic equations where each set may be treated as an independent regression task in N equations (Lall and Voit, 2005).
The concepts of Biochemical System theory (BST) have been used as they have been discussed in the literature several times (Savageau, 1976; Torres and Voit, 2002; Voit, 1991, 2000). In this modeling environment, all processes are represented as products of power law functions which are mathematically derived from Taylor's theorem of numerical analysis that are applied to variables in logarithmic space. Hence, each step Vi which involves at most n dependent (state) variables and m independent (control) variables takes the format
|
| (1) |
i is the rate constant and the exponent fij is the kinetic order that quantifies the effect of variable Xj on Vi.
BST offers alternative representations among which the most widely used are the generalized mass action (GMA) and the S-system representation. In GMA models, each reaction is represented by a product of power law functions which describes all effluxes and influxes and the general GMA structure is
|
| (2) |
The S-system represents the collection of influxes into a given pool with a single power law term and the collection of effluxes from the pool with a second power law term. Hence, the generic S-system structure is
|
| (3) |
As our method focuses on linear chains, the ß1-term can be estimated from the generic format of a linear model as
|
| (4) |
On substitution of numerical values for all variables and slopes at N time steps
|
| (5) |
S.oneidensis is a gram negative facultative anaerobic bacterium found in variety of environments such as freshwater lakes, marine sediments, subsurface formation and variable depths in redox stratified aquatic systems. It is an important model organism due to its metal reducing and bioremediation capabilities. It possesses diverse respiratory capacities: in addition to aerobic respiration it can anaerobically respire inorganic and organic substrates such as fumarate, nitrate, nitrite, thiosulfate, Fe(III), Mn(III), Cr(VI) and U(VI).
For the purpose of this model, time-series data from cells grown in the presence of aqueous species of Fe(III), Co(III), Cr(VI), U(VI) and Tc(VII) as electron acceptors were used. The metal complexes were diluted to 500 µM in 20 mM sodium bicarbonate buffer (pH 6.8) containing 10 mM of either sodium lactate or sodium acetate as the electron donor. The 10 ml reaction mixtures were added to 26.3 ml glass pressure tubes which were made anoxic by bubbling for 10 min with N2:CO2 (80:20). Washed cells which were 
1 ml were added into the anoxic reaction mixture (2 x 108 cells/ml, final concentration) through the rubber stoppers using needles and syringes which were made anoxic by sparging with N2:CO2. Samples were taken at selected time points using degassed needles and syringes and assayed for metal reduction. The details of these experimental procedures can be found in the literature (Liu et al., 2002). For mixed metals, uranium reduction studies with mixed iron hydroxides were conducted in batch cultures under non-growth conditions in 60 ml serum vials which contained 50 ml of modified bicarbonate medium. Serum vials were purged with N2:CO2 at a ratio of 80:20 and after sterilization uranyl acetate and lactate were added from sterile anaerobic stock solutions to give final concentrations of 400 µM and 10 mM, respectively. Solution U(VI) concentration and soluble Fe(II) concentration were determined spectrophotometrically at 575 and 562 nm. Further details on these experimental procedures are in the literature (Wielinga et al., 2000). The metal reduction pathway, which includes the electron transport system (ETS) has been illustrated in Figure 1. Table 1 indicates the reactions of metal reduction by S.oneidensis (Liu et al., 2002).
|
|
The S-system equations for the metal reduction are set up symbolically in accordance with guidelines that have been documented (Voit et al., 2000). The formation of aqueous complex Fe(II) citrate is driven by the concentration of the Fe(III) citrate and by lactate as the electron donor. As there is only a production term for the Fe(II) citrate, the dynamics is given as:
|
| (6) |
X1 = concentration of Fe(II) citrate in micromolar; X7 = concentration of Fe(III) citrate in micromolar and X8 = concentration of Lactate in micromolar.
The finite difference method computes the slopes (Si(tk)) of the time course and these slopes become the estimators of the true differentials.
|
| (7) |
Taking logarithms of both sides yields
|
| (8) |
The N equations of this type represent the N time points considered in the regression. On using the N equations simultaneously, the linear regression determines optimal values for the kinetic parameters 
1, g11 and g12.
The regression in this step was limited to the variables X7 and X8 which was determined by our specific knowledge of the metal reduction pathway. In case nothing about the contributing variables and the regulation of this step is known, the power-law representation will include all variables that can affect this step, for instance, as inhibitors or cofactors. While the regression will still be linear, it will contain more variables. In principle, this would not cause problems but it can lead to redundancies and possible overfitting. As a result, it is advantageous to limit each regression to those variables that are known to affect the process under investigation.
The dynamics of Fe(II)NTA is similar to Fe(II)citrate and is given by:
|
| (9) |
X2 = concentration of Fe(II)NTA in micromolar; X9 = concentration of Fe(III)NTA in micromolar and X10 = concentration of Lactate in micromolar.
Logarithm transformation yields the next set of algebraic regression equation at N time points.
|
| (10) |
In a manner similar as above, linear regression yields estimates for 2, g21 and g22.
The reduction of species of Co(III)EDTA–, (UO2)
, 
and Tc(O4)– was driven by their respective concentrations, the concentration of the lactate and in the case of 
, the concentration of hydrogen in combination with the lactate. It should be noted that the concentration of the lactate as a function of time is different in each of the metal reduction reactions as the experiments for the degradation of the different metals were conducted independently of each other.
As a result, the dynamics of degradation of the aqueous species of Co(III)EDTA–, (UO2)
, 
and Tc(O4)– are given as:
|
| (11) |
X3 = concentration of Co(III)EDTA– in micromolar; X11 = concentration of Lactate in micromolar.
|
| (12) |
X4 = concentration of UO2
in micromolar; X12 = concentration of lactate in micromolar.
|
| (13) |
X5 = concentration of 
in micromolar; X13 = concentration of lactate in micromolar and X14 = concentration of hydrogen in micromolar.
|
| (14) |
X6 = concentration of 
in micromolar; X15 = concentration of lactate in micromolar.
Logarithmic transformation yields sets of algebraic regression equations at N time points for each of Equations (9) through (14) and a subsequent linear regression results in estimates of parameters ß3, h31, h32, ß4, h41, h42, ß5, h51, h52, h53, ß6, h61 and h62.
Table 2 indicates the reactions for mixed metals where the reduction of uranyl (U(VI)) to the relatively insoluble tetravalent form (U(IV)) by S.alga (BR-Y) occurs in the presence of aqueous Fe(III)complexes. The dominant sink for Fe(II) appears to be magnetite (Fe11Fe2111O4) which is a solid that has been noted in the microbial reduction of ferric ion (Benner et al., 2002).
|
The S-system equations for the mixed metals are given by:
|
| (15) |
|
| (16) |
X16 = concentration of C4H6O6U in micromolar; X17 = concentration of Fe2+ in micromolar; X18 = concentration of CH3CHOHCOO– in micromolar; X19 = concentration of Fe(OH)3 in micromolar.
The parameters ß16, h161, h162, 17, g171, g172, ß17, h171 and h172 are estimated with a non-linear least squares formulation that consists of an optimization algorithm (Levenberg–Marquardt) and a numerical ODE solver (Runga–Kutta, ODE45) and is implemented in Matlab.
|
3 RESULTS |
|---|
|
TOPABSTRACT1 INTRODUCTION2 METHODS3 RESULTS4 DISCUSSIONACKNOWLEDGEMENTSREFERENCES |
|---|
Linear regression, in all cases, results in a unique solution. The data was subjected to random noise that was added uniformly from a 10% range about each experimental point. Due to the small number of measurements, the experimental variability was represented better in our approach to test the reliability of results. The artificial addition of noise is a quasi-sensitivity analysis that is beneficial because, in a given set of time-series data the (linear) regression result is unique despite noise within the particular data. By adding random noise, we obtain parameter estimates for slightly modified time courses which reduces the dependence on one set of experiments.
Experimental data were processed by adding 10% random noise which was done to create models that are relatively unaffected by experimental noise. Regressions with the processed data were performed using Matlab and the results are shown in Figure 2. The figure shows plots of the kinetic orders versus the logarithm of the rate constant for metal reduction. The kinetic orders g12, g22 are estimated between 0 and 6, h32, h42, h52 and h53 estimated between 0 and 5 and g11, g21, h31, h41, h51, h61 and h62 estimated between 0 and 2. A comparison of these results shows the degree to which each kinetic order is possibly affected by experimental noise. For instance, in Figure 2, g11 is very stably estimated for a range of artificial noise scenarios, while g12 is not. For the parameters in the higher range of 0–6, a statistical analysis revealed that, in most cases, a 95% confidence interval resulted in parameters falling between 0 and 3. Histogram plots were generated for the parameters using the statistical package R to determine whether they were normally distributed or not. For normal distribution, the 95% confidence interval for the parameters was obtained by using the mean and SD (µ ± 2
) whereas for skewed distribution the quantiles were used. Combinations of parameter values were selected from ranges of kinetic orders and rate constants obtained and models were generated using these values with an ODE solver in Matlab. From the various sets of parameter values it was determined which of these generated models fitted the experimental data. These parameter values have been listed (see Supplementary Material). The results indicated that the S-system model captures the dynamics of the metal reduction pathway quite well (Fig. 3). The model uses the estimated parameters to predict the reduction rates of metals at varying concentrations, and it is observed from these fits that the reduction rates had the following rate trend: Fe(III)citrate > Fe(III)NTA > Co(III)EDTA– > 
> 
> 
(Fig. 4). It is to be noted that there is a high degree of variability in the reduction rates for the different metals. As an example, at a metal concentration of 500 µm, the reaction rate of cobalt is more than 8 times greater than uranium which is a significant increase as compared to the reaction rate of Fe(III)NTA which is 5% greater than cobalt which is a minor increase. These plots can be used directly or extrapolated to determine reduction rates at the concentration levels at which the metal contaminants are present in soil or groundwater which would be useful in the design of bioremediation systems for the clean-up of the toxic waste. For mixed metals, the S-system model captured the reduction of uranyl in the presence of varying ratios of goethite and ferrihydrite (Fig. 5). It is seen from these fits that the rate and extent of uranyl reduction was inversely related to the fraction of ferrihydrite present, and the amount of Fe2+ released to the solution was in proportion to the ferrihydrite present. The optimized parameter values for these fits were determined using a non-linear least squares regression method (Levenberg—Marquardt algorithm). As a result, the S-system model was able to capture the observed dynamics of the mixed metal reduction pathway quite well, where the aqueous Fe(III)complex has played an inhibitory role in the reduction of Uranium. The model uses the estimated parameters to predict the reduction rates of mixed metals at varying ratios of goethite to ferrihydrite. It is seen from these fits that the reduction rates decreases with increasing ratios of ferrihydrite which inhibits the degradation of uranium (Fig. 6). Color versions of Figures 2, 4 and 6 are included in the Supplementary Material.
|
|
|
|
|
As is the case with single metals, the plot is useful for determining reduction rates at concentration levels at which uranium in the presence of aqueous Fe(III)complex are present in aqueous or solid media.
|
4 DISCUSSION |
|---|
|
TOPABSTRACT1 INTRODUCTION2 METHODS3 RESULTS4 DISCUSSIONACKNOWLEDGEMENTSREFERENCES |
|---|
Kinetic models that have evaluated the feasibility of microbial metal reduction have followed Monod-based kinetics, which were used in the design and control of modern microbiological wastewater system (Rittmann and McCarty, 2001; Tchobanoglous and Burton, 1991). The present work uses the S-system representation of the Biochemicals Systems Theory (Lall and Voit, 2005). S-system models can be designed correctly from the topology and regulation of the network. The proposed fitting method is useful for models within BST, as the procedure leads to linear regressions in logarithmic coordinates. However, if the procedure were applied to rate functions such as Michaelis–Menten rate laws, the estimation is somewhat simplified but not linear. Hence, the question arises whether BST models are sufficiently accurate to capture the dynamics of natural systems. It has been suggested by various applications in the past that BST models are of similar accuracy as other models (Alvarez-Vasquez et al., 2004; Curto et al., 1998). In this work, it has been demonstrated that, with accurate measurements that have been assayed for metal reduction, an S-system model provides reasonable parameter estimates that is sufficiently accurate to capture the metal reduction dynamics.
Due to the relatively small number of data points, we added random noise to the data points, which avoids overdependence on a specific set of measurements. We pursued a strategy based on the replacement of differentials with estimated slopes. Simultaneous fitting with slopes was employed, which led to reasonable fits of the time-dependent data. In the case of mixed metals, we used the estimated parameter values for non-linear searches that produced a global well-fitted model satisfying all topological constraints within the pathway. It is to be noted that the experimental data used in this work has been previously modeled using the Monod rate laws (Liu et al., 2002). The Monod-based models have the disadvantage, however, that they become statistically overparameterized resulting in large uncertainties in their parameters. Uncertainties in the parameters arise from the high correlation between them, making them extremely sensitive to experimental error. It was also noted that in many cases the SDs for the Monod parameters were larger than the parameter values themselves. It is vital that one carefully considers whether to evoke complex mechanistic rate functions versus the simpler representation in the S-system model as the former may involve dozens of parameters for a single enzyme-catalyzed step (Savageau, 1976; Schulz, 1994). The model provides insights into the variability of reduction rates for the metals at different concentrations in the soil or groundwater. Such information is useful, as it can be used in conjunction with a reactive transport model to predict groundwater plume lengths and estimate whether a plume has or will migrate off-site. It also enables the prediction of concentrations in the plume at various lengths from the source of the contamination. Subsequently, it can be determined whether the rate of remediation is sufficient to contain the spread of contaminants and what action needs to be taken to protect receiving streams and drinking water sources. For future successful analysis, it is important to develop improved methods that support and automate inverse tasks based on time series of dynamic biological models. This will require a comparative assessment of smoothers and different types of search algorithms. Methods to generate good initial guesses to initiate the non-linear searches can also help identify the right basins of attraction.
|
ACKNOWLEDGEMENTS |
|---|
|
TOPABSTRACT1 INTRODUCTION2 METHODS3 RESULTS4 DISCUSSIONACKNOWLEDGEMENTSREFERENCES |
|---|
This work was supported by the BACTER Institute through a grant from the Department of Energy as part of the Genomics:GTL program (DE-FG02-04ER25627). The authors would like to thank Tim Donohue and Laura Vanderploeg at the University of Wisconsin, Jim Fredrickson and Chongxuan Liu at the Pacific Northwest National lab and Kelvin Gregory at Carnegie Mellon University for their support during the project.
Conflict of Interest: none declared.
|
FOOTNOTES |
|---|
Associate Editor: Alfonso Valencia
Received on March 19, 2007; revised on July 31, 2007; accepted on August 1, 2007
|
REFERENCES |
|---|
|
TOPABSTRACT1 INTRODUCTION2 METHODS3 RESULTS4 DISCUSSIONACKNOWLEDGEMENTSREFERENCES |
|---|
Alvarez-Vasquez F. Integration of kinetic information on yeast sphingolipid metabolism in dynamical pathway models. J. Theor. Biol (2004) 226:265–291.[CrossRef][Web of Science][Medline]
Benner SG, et al. Reductive dissolution and biomineralization of iron hydroxide under dynamic flow conditions. Environ. Sci. Technol (2002) 36:1705–1711.[Medline]
Curto R. Mathematical models of purine metabolism in man. Math. Biosci (1998) 151:1–49.[CrossRef][Web of Science][Medline]
Gaudy AFJ, Gaudy ET. Microbiology for Environmental Scientists and Engineers (1980) New York: McGraw-Hill, Inc.
Gorby YA, Lovley DR. Enzymatic uranium precipitation. Environ. Sci. Technol (1992) 26:205–207.
Gorby YA, et al. Microbial reduction of Co(III)EDTA in the presence and absence of manganese(IV) oxide. Environ. Sci. Technol (1998) 32:244–247.
Lall R, Voit EO. Parameter estimation in modulated, unbranched reaction chains within biochemical systems. Comput. Biol. Chem (2005) 29:309–318.[Web of Science][Medline]
Liu C. Reduction kinetics of Fe(III), Co(III), U(VI), Cr(VI) and Tc(VII) in cultures of dissimilatory metal-reducing bacteria. Biotechnol. Bioeng (2002) 80:637–649.[CrossRef][Web of Science][Medline]
Lloyd JR, Macaskie LE. A novel PhosphorImager-based technique for monitoring the microbial reduction of technetium. Appl. Environ. Microbiol (1996) 62:578–582.
Lovley DR, Phillips EJP. Bioremediation of uranium contamination with enzymatic Uranium Reduction. Environ. Sci. Technol (1992) 26:2228–2234.
Lovley DR. Dissimilatory metal reduction. Annu. Rev. Microbiol (1993) 47:263–290.[CrossRef][Web of Science][Medline]
Lovley DR, et al. Microbial uranium reduction. Nature (1991) 350:413–416.[CrossRef]
Macaskie LE, et al. Bioremediation of uranium-bearing wastewater: biochemical and chemical factors influencing bioprocess application. Biotechnol. Bioeng (1997) 53:100–109.[Medline]
Monod J. The growth of bacterial cultures. Annu. Rev. Microbiol (1949) 3:371–393.[CrossRef][Web of Science]
Nealson K, Saffarini D. Iron and manganese in anaerobic respiration: environmental significance, physiology, and regulation. Annu. Rev. Microbiol (1994) 48:311–343.[CrossRef][Web of Science][Medline]
Rittmann B, McCarty P. Environmental Biotechnology: Principles and Applications (2001) Columbus, OH: McGraw-Hill Higher Education.
Roden EE, Zachara JM. Microbial reduction of crystalline iron(III) oxides: influence of oxide surface area and potential for cell growth. Environ. Sci. Technol (1996) 30:1618–1628.
Savageau MA. Biochemical Systems Analysis: a Study of Function and Design in Molecular Biology (1976) Advanced Book Program, Reading, MA, USA: Addison-Wesley Publication Company. xvii,379.
Schulz RA. Enzyme Kinetics: from Disease to Multi-enzyme Systems (1994) Cambridge, UK: Cambridge University Press.
Stanbury JB, et al. The Metabolic Basis of Inherited Disease (1983) Auckland, New York: Mcgraw-Hill Book Company.
Tchobanoglous G, Burton FL. Wastewater Engineering (1991) 5th. New York: McGraw Hill.
Torres NV, Voit EO. Pathway Analysis and Optimization in Metabolic Engineering (2002) New York: Cambridge University Press. xiv,305.
Voit EO. Canonical Nonlinear Modeling: S-system Approach to Understanding Complexity (1991) New York: Van Nostrand Reinhold. xii,365.
Voit EO. Computational Analysis of Biochemical Systems: a Practical Guide for Biochemists and Molecular Biologists (2000) New York: Cambridge University Press.
Voit EO, et al. Regulation of glycolysis in Lactococcus lactis: an unfinished systems biological case study. IEE Proc. Syst. Biol (2006) 153:286–298.[CrossRef][Medline]
Wielinga B, et al. Inhibition of bacterially promoted uranium reduction: ferric (hydr)oxides as competitive electron acceptors. Environ. Sci. Technol (2000) 34:2190–2195.
Wildung RE, et al. Effect of electron donor and solution chemistry on products of dissimilatory reduction of technetium by Shewanella putrefaciens. Appl. Environ. Microbiol (2000) 66:2451–2460.
CiteULike 
Connotea 
Del.icio.us What's this?
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||