Mathematical models are widely used to create complex biochemical models. Model reduction in order to limit the complexity of a system is an important topic in the analysis of the model. A way to lower the complexity is to identify simple and recurrent sets of reactions and to substitute them with one or more reactions in such a way that the important properties are preserved but the analysis is easier. In this paper we consider the typical recurrent reaction scheme E + S −↽−−⇀− ES −−→ E + P which describes the mechanism that an enzyme, E, binds a substrate, S, and the resulting substrate-bound enzyme, ES, gives rise to the generation of the product, P . If the initial quantities and the reaction rates are known, the temporal behaviour of all the quantities involved in the above reactions can be described exactly by a set of differential equations. It is often the case however that, as not all necessary information is available, only approximate analysis can be carried out. The most well-known approximate approach for the enzyme mec hanism is provided by the kinetics of Michaelis-Menten. We propose , based on the concept of the flow-equivalent server which is used in Pet ri nets to model reduction, an alternative approximate kinetics for the analysis of enzymatic reactions. We evaluate the goodness of the proposed approximation with respect to both the exact analysis and the approximate kinetics of Michaelis and Menten. We show that the proposed new approximate kinetics can be used and gives satisfactory approximation not only in the standard deterministic setting but also in the case when the behaviour is modeled by a stochastic process.
1, Mathematical models are widely used to describe biological pathways because, as it is phrased in [1], they “offer great advantages for integ rating and evaluating information, generating prediction and focusing experimental directions”. In the last few years, high-throughput techniques have increased steadily, leading to the production of a huge volume of data used to derive the complex texture behind the biological/biochemical mechanisms, and creating in this way the structure needed for mathematical modelling. Indeed, many models based on the combination and the integration of various elements in order to investigate their relationships and behaviour have been devised which become more complex with the growth of available data. The complexity is reflected in the number of dynamic state variables and parameters, as well as in the form of the kinetic rate expressions.
Such complexity leads to difficulties both from the point of view of defining the model as the parametrisation becomes unfeasible and for what concerns the analysis of the model. It is often the case hence that in order to have a model which is feasilble for the analysis simplifications must be performed.
In this paper we focus our attention on the simplified, approximate treatment of a set of reactions that very often appears as building blocks of complex models. We consider the reactions
k1
E + S −←k−→−1 ES −k→2 E + P
(
System of enzymatic reactions can be described by Petri nets [2] (Figure 1
shows the Petri net corresponding to the reactions given in (
The concept of flow equivalent server has already been used in [4] where a
complex signal transduction model was considered. In that paper we have shown
that this concept can be applied not only to the small set of reactions given in
(
The paper is organised as follows. Section 2 provides the necessary background, Section 3 describes the concept of the flow equivalent server and Section 4 presents the results of the comparison between the approximation approaches. We conclude with a discussion and an outlook on future works in Section 5. 2
In 1901 Henri [5] proposed a partly reversible reaction scheme to describe the
enzymatic process. According to this scheme the enzyme E and the substrate
S form, through a reversible reaction, the enzyme-substrate complex ES. This
complex can then give rise to the product P through an irreversible reaction
during which the enzyme is freed and can bind again to other molecules of the
substrate. This scheme is summarised in (
There are two typical approaches to associate a quantitative temporal
behaviour to the reactions in (
The deterministic approach describes the temporal behaviour of a reaction
with a set of ordinary differential equations (ODE). For the reactions in (
In the stochastic approach a continuous time Markov chain (CTMC) is used to describe the process. Each state of the chain is described by a vector of integers in which the entries give the quantities of the molecules, which, accordingly, assume discrete values. These discrete values are resulting either directly from molecule count or from discretization of continuous values. Reactions are modeled by transitions between the states. For example, from state | x1, x2, x3, x4| where x1, x2, x3 and x4 are the quantities of the molecules E, S, ES and P , respectively, there is a transition to state | x1 − 1, x2 − 1, x3 + 1, x4| with rate k1x1x2 which corresponds to the binding of one molecule E with one molecule S to form one molecule of ES. It is easy to see that even for small models the corresponding CTMC can have a huge state space whose transition rate structure is non-homogeneous. Exact analytical treatment of the se chains is often unfeasible and in most cases simulation is the only method that can be used for their analysis. 2.1
Under some assumptions, the temporal, quantitative dynamics of the mechanism
described by the reactions in (
= k1[E][S] − [ES](k−1 + k2) = 0
dt
from which, considering that the total amount of enzyme is conserved, i.e. [E] +
[ES] = [E]0, the quantity of ES can be expressed as
where the term kM = k−1+k2 is called the Michaelis-Menten constant. Applying
k1
(
[E]0[S] [ES] = k−1+k2 + [S] = k1 [E]0[S] kM + [S] vMM = k2[E]0[S] Accordingly, after the “plateau” level of ES is reached, the kinetic parameters k1, k−1 and k2 together with [S] and the initial total quantity of the enzyme, [E]0, determine the overall rate of the production of P .
Applying the approximate kinetics of Michaelis and Menten, the differential equations describing the reactions become d[E] dt d[S] dt d[P ] dt
In this section we derive an alternative approximate kinetics for the analysis of
enzymatic reactions, based on the concept of the flow equivalent server. This
technique was originally proposed in the context of the steady-state solution of
queueing networks [3, 8, 9] and can be adapted to our purposes with a proper
interpretation of the assumptions on which it is based. The idea behind this
concept is to consider the reactions given in (
Figure 2 depicts the Petri net corresponding to the reactions of (
In practice, this corresponds to the construction and to the (steady state) solution of the continuous time Markov chain (CTMC) that corresponds to the sub-model in isolation. Providing the speed of the composit e transition in the tabular form highlighted by Figure 3 is convenient for cases where the domain of the function is “small”, but may be impractical in many com mon situations. Despite the computational complexity of the approach, we must notice that the equilibrium assumption of the flow equivalent method is used only to obtain an approximate characterisation of the throughput for different sets of initial conditions and does not mean that the equivalent speed can only be used for steady state analysis.
The concept of flow equivalent server described above is used traditionally in a stochastic setting. However, it can be applied in a deterministic setting as well using arguments that are summarized by the following points. The complexity of the approach in the stochastic setting becomes prohibitive when the amount of the substrate S becomes very large. On the other hand, this is the case in Given k1, k−1, k2, and #E
S1 S2 · · · Sn vF ES(S1) vF ES(S2)
· · · vF ES(Sn)
Fig. 3. Flow Equivalent Server characterisation which the stochastic (or at least the average) behaviour of the model is conveniently captured by a set of ODE, i.e., by a deterministic model. Moreover, in the case of our model, the equilibrium solution of the set of differential equations corresponding to the short-circuited model is simple enough to obtain an analytic expression for the speed of the composite transition as it is described in the following.
We assume that the initial condition is [E]0 = M1, [S]0 = M2, [ES]0 = 0,
and [P ]0 = 0. We will denote the steady state measures of the compounds by
[E], [S], [ES] and [P ]. In the short-circuited version of the reactions given in (
dt d[S](t)
dt d[ES](t) dt = 0 = −k1[E][S] + k−1[ES] + k2[ES] = 0 = −k1[E][S] + k−1[ES] + k2[ES] = 0 = +k1[E][S] − k−1[ES] − k2[ES] which are three dependent equations; – and the following equations expressing conservation of mass
[E] + [ES] = M1, [S] + [ES] = M2
In (
k2 [E] + [S] + kM − p([E] − [S])2 + 2kM ([E] + [S]) + kM2
Accordingly, the set of ordinary differential equations describing the reactions
given in (
k2 [E] + [S] + kM − p([E] − [S])2 + 2kM ([E] + [S]) + kM2 k2 [E] + [S] + kM − p([E] − [S])2 + 2kM ([E] + [S]) + kM2 which explicitly reflects the assumption of the conservation of E and the observation that substrate S is transformed into product P .
In this section, we first compare in Section 4.1 the MM and FES approximate kinetics from the point of view of the speed they assign to the production of P as function of the reaction rates (k1, k−1, k2) and the concentration of the enzyme and the substrate ([E], [S]). Subsequently, in Sections 4.2 and 4.3 we compare the quantitative behaviour of the approximations to that of the full model in the deterministic and in the stochastic setting, respectively.
It is easy to check that as the quantity of the substrate tends to infinity the
two approximate kinetics lead to the the same speed of production. In both cases
for the maximum speed of production we have
vmax = lim vMM = lim vF ES = k2[E]
[S]→∞ [S]→∞
(
A typical way of illustrating the approximate Michaelis-Me nten kinetics is to
plot the production speed against the quantity of the substrate. Figure 4 gives
such illustrations comparing the speeds given by the two approximate kinetics.
Reaction rate k2 is either 0.1, 1 or 10 and reaction rates k1 and k−1 are varied
in order to cover different situations for what concerns the ratio k1/k−1. Two
different values of [E] are considered. The limit behaviours expressed by (
In this section we compare the different kinetics in the deterministic setting. Once
the initial quantities and the reaction rates are defined, the systems of differential
equations given in (
For the first experiments we choose such parameters with which the two approximate kinetics result in different speeds of production. Based on Figure 4 this is achieved whenever the quantity of the enzyme is comparable to the quantity FES with 1,1,0.1 MM with 1,1,0.1 FES with 10,1,0.1 MM with 10,1,0.1 FES with 1,10,0.1 MM with 1,10,0.1
FES with 1,1,1 MM with 1,1,1 FES with 10,1,1 MM with 10,1,1 FES with 1,10,1 MM with 1,10,1 FES with 1,1,10 MM with 1,1,10 FES with 10,1,10 MM with 10,1,10 FES with 1,10,10 MM with 1,10,10 1
FES with 1,1,0.1 MM with 1,1,0.1 FES with 10,1,0.1 MM with 10,1,0.1 FES with 1,10,0.1 MM with 1,10,0.1
FES with 1,1,1 MM with 1,1,1 FES with 10,1,1 MM with 10,1,1 FES with 1,10,1 MM with 1,10,1 FES with 1,1,10 MM with 1,1,10 FES with 10,1,10 MM with 10,1,10 FES with 1,10,10 MM with 1,10,10 on the left side and with [E] = 10 on the right side; reaction rates are given in the legend in order k1, k−1 and k2 of the substrate. Accordingly, we set [E]0 = [S]0 = 10. For the full model [ES]0 needs to be set too, and we choose [ES]0 = 0. This choice does not help the approximations. They assume that the total enzyme concentration [E]0 + [ES]0 is immediately distributed between [E] and [ES], thus making possible an immediate (consistent) production of P . On the contrary, in the full model the production of [ES] takes time and thus the speed of the production of P must start from 0, growing to a high value only later. Figures 5 and 6 depict the quantity of the product and the speed of its production as functions of time for two different sets of reaction rates. In both figures the kinetics based on flow equivalence provides precise approximation of the production of P . The Michaelis-Menten kinetics instead fails to follow the full model, but this is not surprising as the derivation of this kinetics assumes small amount of enzymes. It can also be seen that high values of k1/k−1 (Figure 6) lead to worst approximation in case of Michaelis-Menten kinetics. On the right hand side of the figu res one can observe that for the full model the speed of producing P is 0 at the beginning and then it increases fast to the speed foreseen by the FES approximation. 0 10 20 50 0 10 20 30 40 50
A second set of experiments is illustrated in Figures 7 and 8. We choose sets of parameters with which the speed of production of the MM and FES approximations are similar. In these cases both approximations are close to the reference behaviour. Still, it can be seen that for high values of k1/k−1 (Figure 8) the approximation provided by the Michaelis-Menten kine tics is slightly less precise.
In the following we turn our attention to the cases in which both the approximations are less reliable. In Figure 9 we plotted the case k1 = 0.1, k−1 = 0.1, k2 = 0.1, [E] = 1, [S]0 = 1 and [ES]0 = 0. As mentioned earlier, with [ES]0 = 0 the initial production speed in the original model is 0 while it is immediately high in the approximate kinetics. With low values of k1 and k−1, the time taken by the system to reach the quasi-steady-state situation ass umed by the approximate kinetics is quite long. For this reason there is a longer initial period in which P is produced by the approximations at a “wrong” speed. Furthe rmore, decreasing k1 and k−1 would lead to a longer period in which the approximate kinetics are not precise (see Figure 9).
Another way of “disturbing” the approximations is to dynami cally change the quantity of the substrate in the system. In the original model, because of the intermediate step yielding ES, the speed of producing P changes only after some delay. On the contrary, the approximations react immediately. The harsher the change in the quantity of the substrate the larger is the difference between the original model and the approximations. This phenomenon is reflected in the 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 model by adding the following term to the differential equation that describes the quantity of the substrate: 10(U (t − 5) − U (t − 5.1)) − 10(U (t − 10) − U (t − 10.1)) (12) where U denotes the unit-step function. The effect of (12) is to add 1 u nit of substrate to the system in the time interval [5, 5.1] and to take away 1 unit of substrate from it in the time interval [10, 10.1]. The resulting behaviour is depicted in Figure 10. The approximations change the speed of producing P right after the change in the quantity of the substrate while the original model reacts to the changes in a gradual manner. Naturally, if the quantity of the substrate undergoes several harsh changes then the MM and the FES kinetics can result in bad approximation of the full model. 4.3
In the following we compare the different kinetics in the stochastic setting, by analysing the corresponding CTMCs. In particular, we determine by means of Petri Nets & Concurrency { 67 simulation the average and the variance of the quantity of the product as function of time. The simulations were carried out in Dizzy [11].
The reaction rates for the first set of experiments are k1 = k−1 = k2 = 1. As in the previous section, this choice allows to test a situation where the speed of the two approximations are different. For the same reason, we choose the same initial quantity for the enzyme and the substrate [E]0 = [S]0 = 1. In the stochastic setting the discretization step, denoted by δ, has to be chosen as well. This choice has a strong impact because as the granularity with which the concentrations are modeled is increased, the behaviour of the CTMC tends to the deterministic behaviour of the corresponding ODE. Figures 11 and 12 depict the average and the variance of the quantity of the product with δ = 0.01 and δ = 0.001, respectively. In both figures the approximate kinetics based on flow equivalence gives good approximation of the original average behaviour while the Michaelis-Menten approximation results in too fast pro duction of P . On the right side on the figures one can observe that also the variance is approximated better by the FES approximation. Fig. 12. The average (left) and the variance (right) of the quantity of the product as function of time with k1 = 1, k−1 = 1, k2 = 1, [E]0 = [S]0 = 1 and δ = 0.001
For the second set of experiments we set k1 = 10 and k−1 = k2 = 1 and as initial states we choose again [E]0 = [S]0 = 1. In this case too, as it was shown in Figure 4, the speeds of production of P as predicted by the MM and FES approximations are quite different. Figures 13 and 14 depict the resulting behaviour for two different values of δ. As in case of the deterministic setting, the Michaelis-Menten approximation suffers from the increa sed k1/k−1 ratio and becomes less precise than before. The FES based approach still results in good approximation for both the average and the variance of the production. 0 1 2 5 0 1 2 3 4 5 Fig. 14. The average (left) and the variance (right) of the quantity of the product as function of time with k1 = 10, k−1 = 1, k2 = 1, [E]0 = [S]0 = 1 and δ = 0.001 5
In this paper we have considered the approximate treatment of the basic enzymatic reactions E +S −↽−−⇀− ES −−→ E +P . In particular, an approximate kinetics, based on the concept of flow equivalent server, has been proposed for its analysis. This FES approximate kinetics has been compared to both the exact model and to the most common approximate treatment, namely, the Michaelis-Menten kinetics. We have shown that the FES kinetics is more robust than the one of Michaelis-Menten.
The FES approximation for the basic enzymatic reactions is computationally convenient due to the fact that it has been possible to find an analytic expression for the speed of the composite reaction in this case. While it is very unlikely for this to be true in the case of more complex kinetics, the method is very general and we will study it further within this context to see if it is possible to find other functional expressions for the speed of the composite reaction. One direction of research will be computing the flow equivalent characterization of the kinetics for a number of specific parameter sets and then of constructing the functional representations via interpolation.