Theory (cid:1) of (cid:1) the (cid:1) reactant-stationary (cid:1) kinetics (cid:1) for (cid:1) a (cid:1) coupled enzyme (cid:1) assay

(cid:1) A (cid:1) theoretical (cid:1) analysis (cid:1) is (cid:1) performed (cid:1) on (cid:1) the (cid:1) nonlinear (cid:1) ordinary (cid:1) dif-ferential (cid:1) equations (cid:1) that (cid:1) govern (cid:1) the (cid:1) dynamics (cid:1) of (cid:1) a (cid:1) coupled (cid:1) enzyme (cid:1) catalyzed (cid:1) reaction. (cid:1) The (cid:1) assay (cid:1) consists (cid:1) of (cid:1) a (cid:1) non-observable (cid:1) reaction (cid:1) and (cid:1) an (cid:1) indicator (cid:1) (observable) (cid:1) reaction, (cid:1) where (cid:1) the (cid:1) product (cid:1) of (cid:1) the (cid:1) ﬁrst (cid:1) reaction (cid:1) is (cid:1) the (cid:1) enzyme (cid:1) for (cid:1) the (cid:1) second. (cid:1) Both (cid:1) reactions (cid:1) are (cid:1) governed (cid:1) by (cid:1) the (cid:1) single (cid:1) substrate, (cid:1) single (cid:1) enzyme (cid:1) Michaelis–Menten (cid:1) reaction (cid:1) mechanism. (cid:1) It (cid:1) is (cid:1) demonstrated (cid:1)


Tabulated list of symbols (continued)
Symbol Description Expression κ 1 ratio of first-oder rate constants for non-observable reaction k −1 /k 2 κ 2 ratio of first-order rate constants for indicator reaction k −3 /k 4 σ 1 specific non-observable substrate concentration s 0 1 /K M 1 σ 2 specific indicator substrate concentration s 0 2 /K M 2 V 1 limiting rate of non-observable product E 2 k 2 e 0 1 V 2 limiting rate of indicator reaction product P k 4 s 0 1 η 1 reciprocal of specificity time for non-observable reaction V 1 /K M 1 η 2 reciprocal of specificity time for indicator reaction V 2 /K M 2

Introduction
Difficult detectable enzyme reactions are frequently coupled with easily observable reactions to be studied in enzyme kinetic experiments. An example is the enzyme E 1 reacting with the substrate S 1 to form an intermediate complex C 1 , which then, under the action of the enzyme, forms a product, the enzyme E 2 and releases the enzyme E 1 , following a Michaelis-Menten (MM) [7] reaction mechanisms. This first reaction mechanism (1) represents the non-observable enzyme catalyzed reaction: Then, the enzyme E 2 reacts with the substrate S 2 to form an intermediate complex C 2 , which forms the product P and releases the enzyme E 2 , Mechanism (2) represents the easily observable reaction, which is is generally known as the indicator reaction. In the above chemical schemes, k 1 , k −1 , k 3 , k −3 are microscopic rate constants, and k 2 , k 4 are catalytic constants. Surprisingly, there has been little mathematical investigation of coupled enzyme reactions [16,6,2] due to the nature of non-linear coupling of the differential equation system. Based on the prescribed coupling (1)-(2), we will demonstrate in the forthcoming sections that S 1 and C 1 evolve as a two-dimensional, autonomous, and singularly perturbed equation of the forṁ where ε is a small parameter for the non-observable reaction (0 < ε 1) and f, g : R 2 → R. Singular perturbation problems of this type are ubiquitous in mathematical chemistry [3,9] and biology [1] and the underlying mathematical structure of two-dimensional singularly perturbed dynamical systems that model enzyme kinetics is well understood. It is generally the case that, under certain conditions, and after a short and fast initial transient, the trajectories (solutions to (3)) in phase-space 1 , P, will remain asymptotically close to a one-dimensional slow manifold, M ε , as they approach the stable fixed point (0, 0). The presence of the slow manifold is due to an underlying disparity in the intrinsic chemical timescales. Generally, the complex C 1 will reach a quasi-steady state (QSS) very rapidly relative to the time it takes for S 1 to deplete. The presence of a slow manifold M ε is significant in that it allows for a reduction (simplification) of the model equations (3). This is because on the slow manifold, C 1 is locally a graph over S 1 , i.e., and therefore the two-dimensional problem (3) becomes one-dimensional on M ε :Ṡ We will show that the complete (coupled) system that includes both the non-observable and indicator reactions is a singularly perturbed problem of the form,Ṡ where λ is parameter for the indicator reaction, and is bounded between 0 < λ 1. Singular perturbation problems with multiple small parameters ε, λ are less prevalent in the enzyme chemical kinetics literature, although the original theorems of Tikhonov where developed in a general context and extend to multi-parameter systems in chemical systems [5]. In the forthcoming Sections2 and 3, we will derive the rate expressions for the coupled reaction, as well as expressions for the small parameters ε and λ, and show that the coupled reaction is influenced by a two-dimensional and invariant slow manifold, M * λ,ε , of the form, that attracts nearby trajectories in the four-dimensional phase space P, Consequently, due to the presence of the slow manifold M * λ,ε , the model equations of the coupled reaction (7) can be reduced to a two-dimensional system, and, in Section 4, we exploit this geometrical mathematical structure [8,10] and derive asymptotic approximations to the nonlinear ordinary differential model equations as well as establish regions of validity for the asymptotic solutions. Finally, in Section 5, we conclude with a brief discussion of the results and their relevance in possible future work.

Derivation of the governing equations for the coupled enzyme reactions
Applying the law of mass action to coupled reaction mechanism (1)-(2) yields seven rate equationsė where lowercase letters represent concentrations of the corresponding uppercase species. Typically, in vitro laboratory enzyme assays present the following initial conditions By examining the system of rate equations (11), the coupled reaction mechanism obeys three conservation laws, the enzyme E 1 , substrates S 1 and S 2 conservation laws, Mathematically speaking, the solution trajectory to (11) must lie on the intersection of the hyperplanes defined in (13), which means the original sevendimensional problem can be reduced to a four-dimensional problem. Using (13a) and (13b) to decouple the enzyme concentrations, the redundancies in the system (11) are eliminated to yielḋ where e 1 (t), e 2 (t) and p(t) are readily calculated once s 1 (t), c 1 (t), s 2 (t) and c 2 (t) are known.

Rate expressions for the non-observable enzyme catalyzed reaction
The rate equations (14a)-(14b) are uncoupled from (14c)-(14d). Interestingly, these rate equations have the same structure to those of the single substrate, single enzyme reaction following the MM mechanism. Therefore, it is possible to derive rate equations to model the coupled enzyme catalyzed reaction, and estimate its kinetic parameters using the general theory of the reactantstationary assumption [4] (RSA). The rate equations for the non-observable reaction are identical to those of the single substrate, single enzyme reaction following the MM mechanisms.
3.1 Review of the single substrate, single enzyme MM reaction Revisiting the analysis for the single substrate, single enzyme reaction, we find that there there can be a rapid buildup of c 1 during an initial fast transient of the non-observable reaction step. After this rapid buildup (where the rate of depletion of c 1 approximately equals its rate of formation) c 1 is assumed to be in a QSS,ċ The timescale t c1 is the time associated with the initial transient buildup of where K M1 = (k −1 + k 2 ) /k 1 is the Michaelis constant for the non-observable reaction. The quasi-steady-state assumption (QSSA, 15), in combination with (14a)-(14b), leads to the derivation of the well-known rate expressions (17a (a)-(b)), where V 1 = k 2 e 0 1 is the limiting rate of the non-observable reaction. Hence, the rate equations (14a)-(14b) are reduced to an differential-algebraic equation systems with one single differential equation for s 1 .
Since (17) is only valid after the initial transient, t c1 , it is necessary to define a boundary condition for s 1 at t = t c1 . To find this condition, it can be assumed that there is a negligible decrease in s 1 during the initial transient. This is known as the RSA, and is expressed as The RSA provides an initial condition for (14a) under the variable transformationt → t − t c1 . The mathematical expression (17b) is the MM equation, and the system (17a)-(17b) governs the dynamics of the substrate s 1 and complex c 1 of the observable reaction under the QSS and RSA. The closed-form of (17b) with the initial condition (18), is known as the Schnell-Mendoza equation [12] and is written explicitly in terms of the Lambert-W function From the perspective of asymptotic theory, Schnell and Mendoza [12] have provided a piecewise solution for the single enzyme reaction in terms of a fast transient solution for s 1 , which we will denote as s f 1 , as well as a QSS solution for s 1 , which we denote as s ss 1 : From the earlier work of Segel [13], we have a fast solution "c f 1 " for the complex c 1 , as well as a QSS solution, "c ss 1 ", Collectively, equations (20) and (21) constitute an asymptotic solution that, with fulfillment of some appropriate qualifiers (that we will discuss shortly), serves as an accurate approximation to the true solution of (14).
In addition to the time scale t c1 , which quantifies the length of the initial fast transient (build-up of c 1 ), the time it takes for the majority of the substrate s 1 to diminish is given by t s1 . Although there are several methods for estimating the significant timescales of chemical reactions [15], we employ the heuristic method proposed by Segel [13], and approximate the depletion time to be effectively the total depletion of s 1 (s 0 1 ) divided by the maximum rate of the substrate of depletion after t c1 : Generally speaking, t s1 is a reasonable measure of how long it takes for the non-observable reaction to reach equilibrium. In fact, given (22), it is clear that t s1 is a strict lower bound on the total duration of the reaction.
3.2 Geometrical picture of the enzyme catalyzed reaction, and conditions for the validity of asymptotic solutions While the asymptotic solutions are useful in that they can be used to make certain predictions about the behavior of the system, asymptotic theory fails to yield a more visual/geometric understanding of the problem. Likewise, the asymptotic solutions (20)-(21) do not inherently elucidate their own validity. To paint a complete picture of the mathematical structure behind the single-enzyme, single-substrate reaction, we turn to dynamical systems theory, and analyze this problem from phase-space, P. From the phase-space perspective, after the initial buildup of c 1 , the phase-space trajectory of the single-substrate, single-enzyme reaction (14a-b) hugs a slow manifold, 'M ε ,' and is asymptotic to M ε in the approach to equilibrium. The time it takes for the trajectory to reach the slow manifold is approximately t c1 , while the time it takes for the trajectory to equilibrium is approximated by t s1 . The condition for the validity of the asymptotic solution resides in how well the c 1 -nullcline approximates the slow manifold, M ε and also how straight the phase-space trajectory is in its approach to the slow manifold during the initial fast transient. The former of these conditions is known as the QSSA, and the latter is of course the geometrical interpretation of the RSA. We note that, chemically speaking, if the trajectory is close the slow manifold M ε , then the complex C 1 is assumed to be in a QSS for which the difference and the rate of C 2 depletion is approximately equal to the rate C 2 formation. Mathematically, it was originally proposed that the QSSA was valid if t c1 t s1 , In other words, it was assumed that the c 1 -nullcline should be considered a good approximation to the slow manifold M ε if the timescale accounting for the build-up of c 1 was small compared to the timescale accounting for the duration of the reaction.
As for the validity of the RSA, Segel [14] proposed that one could assume little change in s 1 (an almost straight phase-space trajectory towards the slow manifold) if the depletion of s 1 over the timescale t c1 is minimal: Since |ṡ 1 | ≤ s 0 1 e 0 1 , the strict inequality given in (24) translates to, Through scaling analysis, Segel [13] went on to show that the RSA, ε 1, determines single-handedly the validity of the asymptotic solutions (20) and (21). Introducing the dimensionless variablesŝ 1 = s 1 /s 0 1 andĉ 1 = c 1 /c 1 , Segel and Slemrod [14] demonstrated that, with respect to the dimensionless timescale τ = t/t c1 , equations (14a)-(14a) scale as, Thus, it is apparent from the dimensionless equations (26)-(27b), that if ε 1, then not only will the RSA hold, but the QSSA (which assumes that the c 1 -nullcline is a good approximation to M ε ) also holds. In fact, the RSA, ε 1, is more restrictive than separation of timescales. After some algebraic calculations, the separation of timescales, t c1 /t s1 1, can be written as: where K S1 = k −1 /k 1 , and K 1 = k 2 /k 1 . For the RSA to be valid, the condition must be satisfied, which is more stringent than condition (28), and hence dictates the conditions under which equation (17b) or (19) can be applied. For this reason, it is nowadays considered that MM expressions are valid under the RSA (see Figures 1a and 1b), rather than the QSSA [11].

Scaling analysis of the indicator reaction
Returning our attention to the coupled reaction, we use the same approach developed by Segel [13], and calculate a lower bound, 't s2 ', on the time it takes for the indicator reaction to occur, where k M2 = (k −3 +k 4 )/k 3 is the Michaelis constant for the indicator reaction. Introducing the dimensionless variablesĉ 1 andŝ 1 , as well as the additional dimensionless variablesĉ 2 ,ŝ 2 , into equations (14c)-(14d), it is easy show that, with respect to the dimensionless time t * = t/t s2 , equations (14c-14d) become, with parameters λ, σ 2 , κ 2 and α defined as: Equations (32a)-(32b) imply that if λ 1, then c 2 will be a fast variable with respect to s 2 on the t * timescale, and we can proceed to approximate the solution to the indicator reaction by means of geometric singular perturbation theory by calculating the timescale, t c2 , which roughly accounts for the duration of the initial fast transient of the indicator reaction. First, we impose a RSA on s 2 : Equation (14d) is now linear under the RSA of s 2 : where

The exact solution to the linearized equation is
and therefore the initial fast transient is defined naturally by Dimensional analysis of the indicator reaction supports the equation (37) derived via linearization. Rescaling (14c) and (14d) with respect tot = t/t c2 and s 2 ,ĉ 2 ,ŝ 1 ,ĉ 1 yields, and hence the indicator reaction is a fast-slow system with respective timescalest and t * if λ 1. Furthermore, the identity of λ as a natural small parameter of the indicator reaction can also be derived using the heuristic approach of Segel [13], We remark, that equation (36) does not suggest that there is a rapid accumulation of c 2 over [0, t c2 ], as the accumulation of c 2 over this interval is entirely dependent on h(t). Thus, t c2 is not chemically analogous to t c1 of the nonobservable reaction, in that it is not a measure of how long it takes c 2 to reach the QSS phase.

Asymptotic analysis of the coupled enzyme system
Now that we have obtained expressions for relative timescales t c1 , t s1 , t c2 , t s2 and small parameters ε, λ, we are in a position to derive an asymptotic solution that accurately approximates the indicator reaction. Introducing again the dimensionless variables,ĉ 1 ,ŝ 1 ,ŝ 2 andĉ 2 , equations 14 transform as, where δ is the ratio of the substrate depletion timescales δ = t s1 /t s2 , and the constants ω 1 , ω 2 are bounded below by 1 and denote: From the point of view of scaling, we can derive solutions for two particular cases, which are define by the scale of δ: (i) Case 1: indicator reaction is as fast or faster than non-observable reaction (δ ∼ 1 or δ 1), and Case 2: indicator reaction is slower than non-observable reaction (δ 1).

Case 1: δ ∼ 1 or δ 1
To begin our analysis, we will assume that the indicator reaction happens as fast or faster than the non-observable reaction, and therefore we analyze the model equations (14) with respect to the timescale T . If δ is order unity or larger, then c 2 is a fast variable on the timescale T (as well as on t * ). We begin the analysis by setting λ = 0, for which there exists a three-dimensional critical manifold, M 0 , Several remarks are in order. First, M 0 is hyperbolic and thus structurally stable to perturbations in λ. To see this, we study the adjoined equation, and treat s 1 , c 1 and s 2 as constants. The exact solution of (45) is and therefore c 2 = ψ(s 1 , c 2 , s 2 ) is a stable root of (45). As a result of this hyperbolicity, the perturbed system (i.e., equation (40) with 0 < λ 1) will contain invariant manifolds, M λ , which are smooth perturbations to M 0 : Second, experimental initial conditions lie on the critical manifold M 0 and, as a result, we can assume that when 0 < λ 1 that (with "d(·, ·)" denoting distance) and hence Next, we set ε = 0, from which we define a second critical but one-dimensional manifold M where φ(s 1 ) is given by In fact, the equation c 1 = φ(s 1 ) defines the one-dimensional critical manifold M (1) 0 corresponding to M ε of Figures (1a) and (1b). It is a straightforward calculation to show that M (1) 0 is hyperbolic, thus when 0 < ε 1 there will exist one-dimensional slow manifolds that are embedded within the two-dimensional invariant manifold 2 W s (M ε ): Collectively, the manifolds M λ and M ε define a two-dimensional invariant manifold, M * λ,ε : We remark further that M * λ,ε is embedded in M λ , and equivalently defined as the intersection (see Figure 2 for a three-dimensional schematic drawing), where Γ is also an invariant three-dimensional manifold in P that contains Near the slow manifold M * λ,ε , the dynamics of the complete system are approximately:ṡ Thus, from the point of view of geometric singular perturbation theory, the dynamics of the coupled system are determined by the intersection of invariant, slow manifolds. Furthermore, with the geometric framework in place, we can now supply the timescales t c1 and t c2 with a phase space interpretation: − The timescale t c1 is roughly the time it takes the phase space trajectory to reach Γ . − The timescale t c2 is approximately the time it takes the phase space trajectory to reach M λ . − With experimental initial conditions, the phase space trajectory starts within O(λ) of M λ , and the initial fast transient is characterized by t c1 .
As a direct result from the geometric portrait of the problem in phase space, we see that the timescales t c1 and t c2 are approximate measures on how long it will take for C 1 and C 2 to reach a QSS. Moreover, since experimental conditions imply that trajectories start close to M λ , it can be assumed without a significant loss in accuracy, that the complex C 2 effectively starts in a QSS.
In the case that C 2 starts in a QSS, i.e., with experimental initial conditions, the coupled reaction will consist of only one fast transient that occurs over t c1 . This fast transient occurs while C 1 rapidly approaches a QSS, and, after t c1 , it can be assumed that both C 1 and C 2 are in a QSS. In contrast, if initial conditions are not experimental (i.e., if c 2 (0) > 0), then fast transients in both the indicator and non-observable reactions will occur as both C 1 and C 2 simultaneously approach a QSS. Moving forward, we now proceed to derive approximations to the mass action equations under the assumption that initial conditions are experimental. The term "s 0 1 − s 1 − c 1 " in the numerator of (49) is the total amount of free enzyme (we will denote this as e f 2 ), and it is straightforward to show thaṫ If t < t c1 , then equation (58) is approximately, and therefore, based on the linear theory of (59) over the τ timescale. Inserting (60) into (49) yields Thus, from (61), if separation of timescales t c1 and t s1 is sufficiently large, i.e., then the total production of the complex c 2 is asymptotically negligible during the initial fast transient of the non-observable reaction. This condition is sufficient, but not necessary. With experimental initial conditions, the initial production rate of c 2 is 0 (i.e.,ċ 2 (0) = 0), it can be assumed that the production of c 2 over t c1 is negligible if t c1 is extremely short in duration. For simplicity, we will take c 2 = 0 for t ≤ t c1 . In addition, from the scaling we see that that if δ 2 = t c1 /t s2 is sufficiently small, then the depletion of s 2 on the τ timescale is negligible. As with the condition for the production of c 2 on the τ timescale (62), this condition (64) is sufficient but not necessary, sinceṡ 2 (0) = 0 and therefore the depletion of s 2 on this regime will be small if t c1 is short. As a result, we will take s 2 = s 0 2 for t < t c1 . Geometrically speaking, on the t * and T scales, when t ≥ t c1 , the nonobservable reaction is close to the manifold Γ and the distance from the trajectory to the manifold M * λ,ε becomes arbitrarily small as t → ∞. Thus, after the initial fast approach towards M * λ,ε (characterized by t c1 ), the production and depletion of the intermediate complex c 2 is approximately This expression can be simplified further by noting that the quantity s 0 1 −φ(s 1 ) is bounded from above and below, and therefore assuming ε 1, we have, where Δs 1 = s 0 1 − s 1 . Applying (66) to (14c) and inserting (67) into (14c) yields the leading order approximation toṡ 2 , It is worth noting, that the approximation (68) has a form similar to that of the MM equation. After t c1 , the change in the substrate s 1 is expressible in terms of the Lambert-W function, and equation (69) can then be inserted into (68): The resulting equation, (70), which serves as an approximation to the dynamics of the indicator reaction after the trajectory has reached the slow manifold M * λ, , can be integrated directly using separation of variables, where the term θ(t) is given by In our approximation, since Δs 1 = 0 for t ≤ t c1 , it follows that c 2 = 0 and s 2 = s 0 2 when t = t c1 and hence the asymptotic approximations match. Thus, as a result of the scaling and asymptotic analysis, we obtain the following approximation for the indicator reaction: The power of the asymptotic solution (73) resides in the fact that it is valid over a large range of δ, and is accurate even if the indicator reaction is as fast or faster than the non-observable reaction (see Figures 3a and 3b for a comparison between the numerical and asymptotic solutions).

Case II: δ 1
We now consider the case when δ 1, in which case the non-observable reaction is completed much faster than the indicator reaction. Based on the scaling given in (40c)-(40d), it is apparent that if δ −1 λ ∼ 1, and δ 1, then s 2 is a slow variable on the T timescale. In the case that s 2 is slow over the timescale t s1 , we can approximate s 2 as The inequality given (82) translates to which is a natural condition given the scaled equations (63a)-(63b), and prescribes a condition for the validity of (74) as well as the asymptotic approximation (80)-(80b).

Discussion
Scaling analysis of the mass action equations that model the kinetics of a coupled enzyme assay has revealed two small parameters, ε and λ. These parameters, when sufficiently small, i.e., establish a disparity amongst four timescales: τ, T,t and t * . This disparity gives rise to fast/slow temporal dynamics of the coupled enzyme assay and establishes the presence of slow manifolds in the phase space, P, of the model equations. From the theory of geometric singular perturbations, a geometric framework, in the form of a intersecting slow manifolds, has been established from which coupled enzyme kinetics can be analyzed. From this structure, three manifolds that govern the long-time dynamics of the reaction have been identified: (1) a three-dimensional slow manifold M λ , (2) a three-dimensional invariant manifold Γ , and (3) a two-dimensional slow manifold M * λ,ε . The latter of these manifolds is embedded in the former two, and is defined as their intersection, with M * λ,ε being the dominant slow manifold as t → ∞. Moreover, due to the specific coupling of the governing equations for the coupled reaction, a fourth slow manifold, M ε , exists and lies in a two-dimensional manifold W s (M ε ). The Cartesian product of M ε with R 2 ≥0 defines the three-dimensional manifold Γ : From a chemical point of view, the manifold M λ defines a surface in phase space on which the complex C 2 is in a QSS in which the production of C 2 is roughly equal to the disassociation of C 2 . Likewise, on the invariant manifold Γ , the complex C 1 is in a QSS. These manifolds (Γ, M λ ) intersect along a third invariant manifold, M * λ,ε , on which both complexes C 1 and C 2 are in QSS.
With the manifolds Γ and M λ are associated two fast timescales, t c1 and t c2 , that account for the initial rapid approach of the phase space trajectory towards M λ and Γ respectively. With experimental initial conditions, the phase space trajectory starts within O(λ) from M * λ . The major significance of this observation is that, when experimental initial conditions are present, the model reduction of (14) to the differential-algebraic equation, is valid when t = 0 provided λ 1. Thus, when experimental initial conditions are prescribed, there is only one fast transient that occurs over t c1 , and is a rough measure of the length of time it takes the phase space trajectory to reach M * λ,ε . Two additional timescales, t s1 and t s2 , have been derived and account for respective amount of time it takes for the substrates s 1 and s 2 to completely deplete: Asymptotic solutions to the rate equations (14), as well as conditions for their validity have been derived that are accurate when δ (the ratio of the substrate depletion timescales for non-observable and indicator reactions, δ = t s1 /t s2 ) is δ 1, δ 1, and δ ∼ 1. The asymptotic solutions for the substrate concentrations s 1 and s 2 are given in terms of Lambert-W functions, and, in the limiting case when the indicator reaction is completed much later than the non-observable reaction (for which δ 1), the asymptotic approximation for s 2 is homologous to the Schnell-Mendoza equation for single-enzyme, singlesubstrate reaction following the MM-type mechanisms: The asymptotic solutions and their associated conditions for validity, timescales, small parameters, and geometric structure of the problem through the lens of phase space serve as the major contributions of this paper. These solutions are useful not only in the insight they provide pertaining to the dynamics of the coupled reaction, but can also possibly be used to analyze, at least from a theoretical platform, the inverse problem [17], in which substrate curves are generated with the hope of estimating the kinetic parameters K M1 , V 1 , K M2 and V 2 . We also expect that the theoretical framework introduced in this paper is applicable to other type of coupled enzyme catalyzed reactions.