Using singular perturbation theory to determine kinetic parameters in a non-standard coupled enzyme assay

We investigate how to characterize the kinetic parameters of an aminotransaminase using a non-standard coupled (or auxiliary) enzyme assay, where the peculiarity arises for two reasons. First, one of the products of the auxiliary enzyme is a substrate for the primary enzyme and, second, we explicitly account for the reversibility of the auxiliary enzyme reaction. Using singular perturbation theory, we characterize the two distinguished asymptotic limits in terms of the strength of the reverse reaction, which allows us to determine how to deduce the kinetic parameters of the primary enzyme for a characterized auxiliary enzyme. This establishes a parameter-estimation algorithm that is applicable more generally to similar reaction networks. We demonstrate the applicability of our theory by performing enzyme assays to characterize a novel putative aminotransaminase enzyme, CnAptA (UniProtKB Q0KEZ8) from Cupriavidus necator H16, for two different omega-amino acid substrates.


Introduction
Enzyme assays are an important tool for characterizing enzymes. In the classic assay, a single reaction converts a substrate into a product, using an enzyme as a catalyst, and the product is measured over time to estimate the initial reaction velocity (Bisswanger 2014). A mathematical analysis allows one to use this information to deduce the kinetic properties of the enzyme in question (Murray 2002). When characterizing an enzyme for which the reaction product is difficult to observe, an auxiliary enzyme can be introduced to convert the product of the primary enzyme reaction into a chemical that can more easily be measured (e.g. NADH) (Storer and Cornish-Bowden 1974;Rudolph et al. 1979). This is called a coupled enzyme assay. However, the problem is more complicated if the reaction catalysed by the auxiliary enzyme has a product which is also a substrate for the primary enzyme.
In this paper, we are interested in understanding how to characterize a primary enzyme in a non-standard coupled enzyme assay, where one of the products of the auxiliary enzyme is a substrate for the primary enzyme and where we explicitly account for the reversibility of the auxiliary enzyme reaction. This is motivated by the characterization of the class of enzymes known as aminotransaminases, which catalyse the transfer of keto and amine groups between organic compounds (Slabu et al. 2017), and are an important class of biocatalysts in synthetic biology for reasons we discuss below. In Fig. 1 we show a schematic of the reaction network we consider in this paper. While we refer to the specific organic compounds pyruvate and alpha-L-alanine in Fig. 1 and through the paper in order to be consistent with the experiments we perform, we emphasize that our analysis is more general. In particular, the general organic-group-transfer property of aminotransaminases means that 'pyruvate' can be generalized to 'keto acid', and 'alpha-L-alanine' can be generalized to 'amino acid' for an aminotransaminase primary enzyme. Since ammonia is produced when an amino acid is dehydrogenated, the presence of ammonia is a consequence of using an aminotransaminase as the primary enzyme. This means that the reaction network we show in Fig. 1 is particularly applicable to coupled enzyme assays with the primary enzyme being an aminotransaminase.
The technique to deal with a standard coupled enzyme assay with one primary enzyme and one auxiliary enzyme was first outlined within a mathematical framework in McClure (1969). We direct the interested reader to Rudolph et al. (1979) for a history of the theoretical work on standard coupled enzyme assays, as well as a practical guide to performing these types of assays. In this paper, we investigate when the standard approach identified in McClure (1969) is appropriate for the non-standard assay we consider, and how to modify the approach when it is not. We then implement our theory by performing assays on a previously-uncharacterized enzyme for two different substrates.
Calculating enzyme properties through enzyme assays is an inherently dynamic process; the steady state is reached only when one of the substrates vanishes. In this dynamic process, it is helpful to measure the indicator chemical during a regime where the reaction velocity is approximately constant, corresponding to a linear increase of the indicator chemical in time (we henceforth refer to this as a 'linear growth' regime). Often this means estimating the reaction velocity through initial rate experiments, Pyruvate (s 1 ) Substrate (s 2 ) Aldehyde (s 3 ) Alpha-L-alanine (s 4 ) NAD + (s 5 ) NADH (s 6 ) Ammonia (s 7 )ṽ 1 v 2ṽ−2 Fig. 1 A diagram of the reaction network we investigate in this paper. Each box denotes a chemical compound involved in the network, and each arrow denotes a reaction to and from the specified group of products and substrates, respectively. Each reaction is catalysed by an enzyme; the primary enzyme catalyses the reaction on the right, marked with velocityṽ 1 , and the auxiliary (secondary) enzyme catalyses both directions of the reversible reaction on the left, marked with forward velocityṽ 2 and backwards velocityṽ −2 . Note that our model does not rely on pyruvate being converted into alpha-L-alanine, it being valid for the primary enzyme being a general aminotransaminase. Hence, in this diagram 'pyruvate' can be generalized to 'keto acid', and 'alpha-L-alanine' can be generalized to 'amino acid' with the goal of avoiding complications due to reversible reactions. Linear growth is helpful as a sanity check that things are proceeding as they should (typically using the 'eyeball norm'). Moreover, experimental noise can increase the error when attempting to implement an accurate fit in nonlinear systems -restricting oneself to a linear growth regime significantly reduces this issue.
For these reasons, we use singular perturbation theory (Bender and Orszag 2013;Kevorkian and Cole 2013) to analyse the mathematical systems we derive, and to understand how to determine the kinetic parameter values of the primary enzyme in such a system. This provides a significant reduction of the complexity of the system, and thus allows us to minimize issues associated with experimental noise, as we aim to determine functional forms for the measurable reaction velocities. Moreover, as exhibited in recent examples (Dalwadi et al. 2018a, b;Kumar and Josić 2011;, singular perturbation theory is particularly well suited to understanding systems of chemical reactions exhibiting a significant separation of timescales, so we will use a similar approach here. In our analysis, we assume that the timescales of enzyme complex formation are much shorter than the other reaction timescales in our problem, and that the effective reaction velocities are governed by Michaelis-Menten-type laws. The suitability of the Michaelis-Menten equations in many scenarios, including but not limited to standard coupled enzyme assays, has been investigated in recent work by Schnell and colleagues, for example Stroberg and Schnell (2016); ; , using singular perturbation theory.
To demonstrate how our theoretical results can be used to characterize enzymes, we perform several coupled enzyme assays to characterize a putative omegaaminotransferase (also known as an aminotransaminase or an omega-amino-acid amino transferase) enzyme, CnAptA, for two omega-amino acid substrates. Omegaaminotransferases belong to class III aminotransferases; these catalyse the transamination of omega-amino acids such as beta-alanine or gamma-aminobutyric acid by transferring the amino group to a keto acid, using pyridoxal 5 -phosphate (PLP) as a cofactor (Mehta et al. 1993). Aminotransferases are increasingly important biocatalysts for the synthesis of industrially relevant chiral compounds and in the pharmaceutical industry for the production of optically pure amines and amino alcohols necessary for the synthesis of important drugs (Sayer et al. 2013). Moreover, aminotransferases are particularly important in synthetic biology because they are the key to creating biosustainable production routes to many important platform chemicals. That is, carboxylic acids make up eight of the top twelve platform chemicals selected by the US Department of Energy that can be derived from biomass (Werpy and Petersen 2004). As current industrial methods to produce these acids often involve fossil fuels, biological production routes provide environmentally sustainable alternatives.
The outline of this paper is as follows. In Sect. 2, we briefly recap the mathematical analysis of the standard enzyme assay, which allows us to put the main results of this paper into an appropriate context. In Sect. 3, we derive a mathematical model for the non-standard coupled enzyme assay of interest in this paper, which has the reaction network shown in Fig. 1. We then explore how the non-standard assay differs from the standard one and determine how the measured reaction velocity relates to the kinetic parameters of the primary enzyme in the non-standard case. We do this using an asymptotic analysis to systematically reduce the complexity of the system, and we present the main distinguished limit (i.e. the dominant balance of most practical relevance) of the system in Sect. 3.1. We present the remaining distinguished limit in Sect. 3.2, where we find that the distinguished nature of this limit is apparent only in a non-measurable chemical during the linear regime or during the depletion regime, and that the measurable system behaviour in the linear regime is captured by a sub-limit of the strong reverse reaction regime treated in Sect. 3.1. We summarise our theoretical results in Sect. 4, outlining the procedure that should be implemented to use experimental data to infer the kinetic parameters of the primary enzyme. In the same section, we compare our proposed method with a naive nonlinear fit, verifying that our method is preferable. In Sect. 5, we demonstrate that our results can be used effectively in characterizing enzymes by performing the coupled enzyme assay modelled in this paper and applying our theoretical results. Finally, in Sect. 6, we discuss our results in the context of their significance for experimental approaches involving non-standard coupled enzyme assays.

Standard coupled enzyme assay
We start by briefly recapping the mathematical model for the standard coupled enzyme assay, presented in broad dimensional terms in McClure (1969). In this section, we present the relevant chemicals involved in abstract terms before relating them to the chemicals labelled in Fig. 1 for the non-standard coupled enzyme assay we consider in the next section. Our goal in this section is to recap how to determine the kinetic parameters of the primary enzyme in the system shown in Fig. 2, where we assume knowledge of the kinetic parameters of the auxiliary enzyme. Experimentally, we are able to measure the concentration ofs 6 , and so our overarching goal is to understand how to use this measurement to infer the kinetic parameters that characterize the primary enzyme. In general, the kinetic parameters used to characterize an enzyme are intrinsic to the ordinary differential equations (ODEs) assumed to govern this system, using Michaelis-Menten-type laws to quantify each reaction velocity (Menten and Michaelis 1913). The network in Fig. 2 implies the following seven dimensional ODEs We use the convention that dimensional quantities are marked with a tilde, whereas dimensionless quantities are not. The initial conditions arẽ which corresponds to a system initially containing the minimal number of substrates necessary for both reactions to take place:s 1 ,s 2 , ands 5 , then instantaneously adding the primary and auxiliary enzyme such that the entire system is well-mixed. We assume that the affinity a given substrate has for an enzyme is independent of the other substrates. Therefore, the reaction velocities follow Michaelis-Menten kinetics generalized to multiple substrates (Alberty 1953) This problem was recently investigated through a singular perturbation analysis in  to determine when the Michaelis-Menten reaction-type formulation of the problem is valid by considering the intermediate enzyme complexes formed during each reaction. We will assume throughout this paper that we are in this regime, essentially assuming the reactant-stationary assumption holds, i.e. that the primary substrates are approximately constant as the primary intermediate complexes form.
In (3), we use the convention that a subscript refers to a relationship with a particular chemical, and a bracketed superscript refers to a particular reaction. For example, the parameterK (2) 4 refers to the Michaelis constant fors 4 in the reaction velocityṽ 2 . In general,α andβ can be viewed as control parameters that can be varied between experiments; all remaining parameters should remain constant between different experiments (at fixed temperatures and pHs), such assumptions being appropriate to the experimental approach we describe in Sect. 5. Hence, we are interested in inferring the three parameters with a bracketed superscript of 1, which arek (1) ,K (1) 1 , andK (1) 2 , from measurements ofs 6 .
It is apparent from the system (1)-(3), as well as from Fig. 2, that the concentrations fors 3 ,s 6 , ands 7 decouple from the rest of the system. We include them here to keep the same notation as for Sect. 3. Moreover, we are able to reduce the system to two ODEs by noting the five following linearly independent conserved quantities: s 1 +s 3 =α,s 2 +s 3 =β,s 5 +s 6 =γ ,s 2 +s 4 +s 6 =β,s 6 =s 7 . (4) We delay the system reduction to discuss first the following standard assumptions made when carrying out a coupled enzyme assay, as this will aid in deciding with which ODEs to work. First, we assume that the maximum rate of reaction mediated by the auxiliary enzyme is much faster than for the primary enzyme, so thatk (2) k (1) . As the maximum rate of each reaction is proportional to the concentration of the corresponding enzyme, this can be ensured by using a high enough ratio of auxiliary to primary enzyme. Second, we will initially saturate the system with large quantities ofs 5 , so thatγ K (2) 5 , thus reducing the complexity ofṽ 2 in (3). Third, we will initially saturate the system with eithers 1 ors 2 (so thatα K (1) 1 orβ K (1) 2 ) if we are interested in determiningK (1) 2 orK (1) 1 , respectively. This reduces the complexity ofṽ 1 in (3). This final assumption is not required to facilitate a solution of the reduced system, though it does reduce the number of fitting parameters. In Sect. 3, we only make the first two assumptions in our analysis. Not making the third assumption allows us to consider the separate cases of saturatings 1 ands 2 at the same time.
We will present the method for determiningK (1) 1 , essentially the case with a large initial saturation ofs 2 . The approach to determiningK (1) 2 is equivalent for the standard coupled enzyme assay discussed in this section, but this will not be the case for the non-standard assay we consider in Sect. 3. If we initially saturate the system with large quantities ofs 2 , it is convenient to use the following two ODEs to describe the system: This is the standard system to consider when using coupled enzyme assays, investigated in the literature in McClure (1969), as discussed in Sect. 1.

Fig. 3
The measurable chemicals 6 as the numerical solution to (1), the asymptotic solution (8) (prior to the depletion dynamics), and the long-time asymptotic solution (9) (also prior to the depletion dynamics). We use the parameter values in Table 2, choosingα = 2.5 mM and usingβ =γ = 50 mM to focus on the effect of a limiting initial concentration ofs 1 , rather thans 2 ors 5 To implement our systematic asymptotic analysis, we form dimensionless variables as in Table 1, and obtain the following system where typical values and definitions of the dimensional and dimensionless parameters are given in Tables 2 and 3, respectively. The most important definition here is = k (1) /k (2) which, as discussed above, is small in coupled enzyme assays. We also define a =K (2) 4 /α and K (1)

=K
(1) 1 /α. In the limit of → 0, with t = O(1), (6) has leading-order solution which corresponds to the dimensional solutioñ for our measurable chemical. The long-time limit of (8) is Therefore, the kinetic parametersk (1) andK (1) 1 can be inferred by measuring the long-term (by which we meanK (2) 4 /k (2) t α/k (1) , for reasons noted below) production rate ofs 6 for different values ofα. A similar procedure, but with saturating s 1 instead ofs 2 , allows the inference ofk (1) (again) andK (1) 2 . The solution (9) tells us that the long-term production rate will beαk (1) /(K (1) 1 +α), as can be seen in Fig. 3, and that we requiret K (2) 4 /k (2) befores 6 production will appear to be at a constant rate (K (2) 4 /k (2) ≈ 200 s in Fig. 3). However, this constant rate will not occur indefinitelygiven that ds 1 /dt = O(a ), the production rate ofs 6 will start to deviate from constant when t = O(1/a ), equivalent tot = O(α/k (1) ) = O(2.5 × 10 4 s) in Fig. 3. In the language of matched asymptotic expansions, this corresponds to an outer scaling, a timescale over which (what we term) the depletion dynamics occur. For brevity, we do not consider the depletion dynamics further for the standard enzyme assay, but we will for the non-standard enzyme assay.

Non-standard coupled enzyme assay
We now consider the main problem we are concerned with in this paper -the determination of the kinetic parametersk (1) ,K (1)

, andK
(1) 2 of the primary enzyme in the system shown in Fig. 1, where we assume knowledge of the kinetic parameters of the auxiliary enzyme. Experimentally, we are able to measure the concentration of NADH (i.e.s 6 ) in the system, and so our overarching goal is to understand how to use this measurement to infer the kinetic parameters that characterize the primary enzyme for a range of different substrates. The differences between the standard coupled enzyme assay in Sect. 2 and the problem we consider in this section are that the latter involves a reversible auxiliary reaction and thats 1 is a product/substrate of the forward/backward auxiliary reaction, rather than the auxiliary reaction being unidirectional ands 1 being decoupled from the auxiliary reaction in Sect. 2. The reaction network we consider is particularly relevant to the class of enzymes known as aminotransferases, 1 as discussed in the Introduction.
The governing equations for the network shown in Fig. 1 consist of the following seven ODEs where we define each variable in Table 1. As in the previous section, we use the convention that dimensional quantities are marked with a tilde, whereas dimensionless s 4 Alpha-L-alanines 4 = (k (1)K (2) 4 /k (2) )s 4 s 5 NAD +s 5 =γ s 5 quantities are not. Here, the three reaction velocities arẽ each following Michaelis-Menten kinetics generalized to multiple substrates (Alberty 1953), under the assumption that the affinity a given substrate has for an enzyme is independent of the other substrates.
In conjunction with the governing equations (10), we use initial conditions whereα,β,γ > 0. This corresponds to a system initially containing pyruvate (s 1 ), an omega-amino acid substrate (s 2 ), and NAD + (s 5 ), with no other substrates present in the reaction network, then instantaneously adding the primary and auxiliary enzyme such that the entire system is well-mixed. Hence, there are thirteen dimensional parameters. We reiterate that our goal is to infer the three parametersk (1) ,K

, andK
(1) 2 that correspond to the primary enzyme, assuming that we know the remaining ten parameters, which correspond to the auxiliary enzyme and the initial chemical concentrations. We provide typical parameter values in Table 2.
As would be expected on physical grounds, the inclusion of a reverse reaction for the auxiliary enzyme causes less NADH to be made over a given time (Fig. 4a). Moreover, the reverse reaction can result in a linear phase that is very short (Fig. 4b), making the standard approach of calculating a steady production rate difficult to implement in practice. To this end, we will characterize the different possible types of NADH Value determined by choice ofk (2) and assay environment.
For the auxiliary enzyme, we use values taken from Yoshida and Freese (1965) in our numerical simulations. We give an idea of the range of these parameter values by also providing data from four other organisms: Enterobacter aerogenes (Chowdhury et al. 1998), Mycobacterium tuberculosis (Hutter and Singh 1999), Phormidium lapideum (Sawa et al. 1994), and Rhodopseudomonas capsulata (Hutter and Singh 1999) production through a systematic asymptotic analysis-see Bender and Orszag 2013; Kevorkian and Cole 2013 for general descriptions of the approach. Before we nondimensionalize the system, we note from (10) and (12) the five following linearly independent conserved quantities: s 1 +s 4 =α,s 2 +s 3 =β,s 5 +s 6 =γ ,s 2 +s 4 +s 6 =β,s 6 =s 7 . (13) Thus, we may immediately reduce the number of ODEs from seven to two. However, the asymptotic analysis we will carry out is more intuitive if we keep most of the ODEs; the only ones we disregard are those for ds 3 /dt and ds 7 /dt in (10). For the former, we do this since the aldehyde produced from the primary enzyme (s 3 ) decouples from the rest of the system, so we are able ignore this ODE henceforth, noting thats 3 can be determined from (13). For the latter, the concentration of ammonia is always equivalent to the concentration of NADH (s 7 ≡s 6 ), so we will replaces 7 bys 6 henceforth.
We form the following dimensionless variables For the auxiliary enzyme, we use values taken from Yoshida and Freese (1965) in our numerical simulations. We give an idea of the range of these parameter values by also providing data from four other organisms: Enterobacter aerogenes (Chowdhury et al. 1998), Mycobacterium tuberculosis (Hutter and Singh 1999), Phormidium lapideum (Sawa et al. 1994), and Rhodopseudomonas capsulata (Hutter and Singh 1999) The measurable chemicals 6 as the numerical solution to (10) in a a linear plot and b a log-log plot. We use the parameter values in Table 2, withk (−2) = 0.5 mM s −1 ,α = 2.5 mM, andβ = 2.5 mM which we also summarize in Table 2. The timescale here (i.e.K (2) 4 /k (2) ) is that of initial alpha-L-alanine (s 4 ) and NADH (s 6 ) formation, and is 200 s for the parameters in Table 2. For the substrate dimensionless variables we scale pyruvate, beta-alanine, and NAD + with their initial conditions and alpha-L-alanine and NADH with the order of magnitude of alpha-L-alanine present during the regime of linear NADH growth, which our analysis will show to be approximately 50 µM. From (14), we obtain the dimensionless governing equations where the three dimensionless reaction velocities are defined as The initial conditions of the dimensionless system are We provide the definitions of the dimensionless parameters in Table 3, from which we note that , δ, and K are all much smaller than 1. Additionally, we note that the value of k (−2) , representing the ratio of the backward to forward rates controlled by the auxiliary enzyme, can vary across many orders of magnitude. We reiterate that our goal is to understand how to relate the measurement of NADH (s 6 ) in the linear production regime to the values of the kinetic parameters governing the primary enzyme. Our approach is to comprehensively investigate how the system behaves for different magnitudes of k (−2) , with a focus on understanding when a linear regime for NADH production is appropriate and what the measurement of NADH in this regime will tell us about the primary enzyme. To this end, we perform an asymptotic analysis that exploits the small parameters , δ, and K (noting that the three limits all commute), with a focus on analysing how the asymptotic size of Table 3 The dimensionless parameters in our system, obtained using the data in Table 2 =k (1) The 'approximately equals to' signs arise from uncertainty in the kinetic parameters of the primary enzyme k (−2) affects the system behaviour. We emphasize that the parameter can always be made to be small by suitable choices of the relative concentrations of the primary and auxiliary enzyme, which scale withk (1) andk (2) , respectively, and likewise for δ by choosingγ , the initial concentration ofs 5 , to be significantly larger thanK (2) 5 , the Michaelis constant fors 5 in the forward auxiliary reaction. Since the smallness of these parameters reduces the complexity of the system, it is experimentally favourable to choose the relative enzyme concentrations and initial conditions to make and δ small. We note that the smallness of K is slightly different, however. While the typical parameter values we give in Table 2 unanimously result in small values of K , there seems no fundamental reason why it should be small. To ensure the generality of our analysis, we therefore also present results for K = O(1) in "Appendix C".
It turns out that there are two distinguished asymptotic limits in this system (so that all other regimes with , δ, and K each small are sub-cases of these two), whereby 1, respectively. We investigate the former in the next section, discussing its relevance to our inference problem, and we investigate the latter in Sect. 3.2, where we show that it can be considered a sub-limit of the former case in terms of experimental measurements. Moreover, we show that the latter regime exhibits the same constant 'long-time' production rate of NADH as the standard coupled enzyme assay case we summarized in Sect. 2. In Table 4, we provide a summary of the results we derive in the remainder of this section.

Strong reverse reaction: k (−2) = O(1/K)
When k (−2) = O(1/K ), the auxiliary-enzyme-controlled reaction in the forward direction is much stronger than for the reverse direction. This scaling results in a distinguished limit where the nonlinearity of the reverse reaction is important when t = O(1), which is where the linear growth regime occurs for the standard coupled enzyme assay. We now exploit the limits , δ, K 1 with k (−2) K = O(1). We note that while small values of would result in large values of K (−2) 6 , we will treat K (−2) 6 = O(1) (as we do with the remaining dimensionless parameters in the system) Table 4 A summary of the solutions we derive in this paper for each timescale in the two distinguished limits of the problem, as well as the dimensionless parameters that control the concentration of s 6 in each asymptotic region We give the equation number for the definition of each parameter unless it is given in Table 3 as this choice keeps more terms at leading order (i.e. constitutes a distinguished limit) and is the case we expect from experimental parameters.

Linear-phase regime
We first consider the linear-phase regime, which occurs over the timescale t = O(1). By 'linear-phase', we mean that s 6 will exhibit linear growth for some time period which we can approximate. In the language of matched asymptotic expansions, this corresponds to linear growth in the intermediate matching region between this timescale and the longer depletion timescale (which we will show corresponds to In the limits set out in Sect. 3.1, the leading-order version of (15) for t = O(1) is with initial conditions (17). We may immediately deduce that, over this timescale, where we define We can therefore reduce the system (17)-(18) to the following single ODE where we define  (23) and (24), respectively. We use parameter values Although we cannot solve (21) explicitly, we may deduce from (21) that the earlytime behaviour is and the long-time behaviour is These asymptotic solutions are good approximations of the numerical solution to (21) in their regions of validity (Fig. 5). It will also be helpful to note that Hence when measuring s 6 experimentally, we will not observe linear growth in time immediately (rather, (23) implies that we will observe a t 2 relationship), but this linear relationship will start to develop after a lag period. To help understand further the solution of the ODE (21), we first note that we can reduce its dependence from three parameters to two by making the scaling s 6 = K Y . This results in the system A comparison of numerical and asymptotic solutions of the reduced ODE (26) for aV = 10 −2 , γ = 10 0 , bV = 10 −2 , γ = 10 2 cV = 10 4 , γ = 10 2 , dV = 10 2 , γ = 10 0 . The onset of the linear regime is set as 2.5 times the right-hand side of (30), and is denoted by a cross/asterisk when the first/second term on the right-hand side is maximal where we defineV Solutions to (26) generally exhibit a transient region before settling to linear growth for a wide range of parameter values (Fig. 6). The long-time linear growth exhibited by the ODE (26) is and it is of interest to understand when this linear growth becomes apparent. One way to do this is to examine the menagerie of distinguished asymptotic limits of the ODE (26), which are interesting in their own right. Moreover, the analytic solutions afforded by an asymptotic analysis can be used to fit time-series data from experimental results more easily than numerical solutions of a system of ODEs. We investigate the subsidiary distinguished asymptotic limits in terms of the two parametersV and γ in "Appendix A", and we show a comparison of numerical and asymptotic solutions to (26) in Fig. 6. We illustrate the four distinguished limits that we obtain from the analysis of "Appendix A" in Fig. 7. A result from this analysis with particular importance for accurately interpreting the data is that there are two distinct linear regimes in the distinguished limit I.  Fig. 7 The distinguished asymptotic limits of the reduced ODE (26) for asymptotically large or small values ofV and γ , denoted by dark solid lines and labelled using Roman numerals contained within a circle. Note that the distinguished limit I corresponds toV 1 and γ = O(1), the distinguished limit II corresponds toV γ = O(1) withV 1 and γ 1, the distinguished limit III corresponds toV /γ 2 = O(1) with V 1 and γ 1, and the distinguished limit IV corresponds toV 1 and γ = O(1). The sublimits of the distinguished asymptotic limits are noted using letters contained within a square While the duration of this transient region can be calculated numerically for given parameter values or approximated analytically from the asymptotic results we derive in "Appendix A", it is helpful to obtain an approximate analytic estimate for when we expect to see a linear relationship for the general case. To this end, we first note that where f 1 and f 2 satisfy the following ODEs which yield the following bounds We note that the long-time limit of the lower bound in (29) is the long-time limit (27) of the ODE (26) since the nonlinear term in the ODE becomes negligible when Y 1. Mathematically, the lower bound of (29) appears linear when (1 + γ )t 1. Combining these two constraints, we expect the NADH concentration to grow linearly in time when As 1/(1 + γ ) is bounded above by 1, it is the second of these constraints that can lead to very large times before the linear regime is encountered. This problem is exacerbated when γ 1 andV = V /K (−2) 6 1, as can be seen by the significantly large times until linearity for Fig. 6b. Moreover, we note that the two constraints in (30) approximately balance whenV = O(γ 2 ) for large γ (essentially Regime III in Fig. 7), and whereV = O(1) for small γ . This represents a boundary between the two different constraints, where the large-V case (on the right of Fig. 7) corresponds to t 1/(1 + γ ). When we denote the onset of the linear regime in Fig. 6, we use an asterisk when (1 + γ )/V > 1/(1 + γ ) and a cross when (1 + γ )/V < 1/(1 + γ ).
In each case, we take for definiteness the onset to occur at 2.5 times the value of the right-hand side of (30).

Physical interpretation
To understand the physical implications of our analysis from the linear-phase regime, it is helpful to reframe our results in terms of dimensional quantities. Hence, the long-time dimensional linear behaviour of NADH, given in (24) for the dimensionless problem, iss at leading order, where Here, ω can be thought of as a measure of the relative strengths of the forward and backward auxiliary reactions: ω → 0 + when the backward reaction is much stronger, but ω → 1 − when the forward reaction is much stronger. Additionally,Ṽ can be thought of as the 'natural' observed strength of the reaction, as this is the reaction strength when the backward auxiliary reaction is unimportant, as we show in Sect. 3.2. Moreover, ω can be calculated from the experimental setup, so the standard practice of inferring the unknown parameters for a single enzyme can be used, as long as the observed reaction velocity is modified by a factor of 1/ω. Thus, we have determined the reaction velocity in the linear-phase regime in terms of the system parameters. Our remaining task is to understand when we expect to observe this regime. The dimensional version of the lower constraint for the linear regime, given in (30) for the dimensionless problem, ist max ωK Hence, this gives us the approximate time we expect to wait before observing a linear relationship between time and NADH production. The constraint (33) provides a warning against particularly low initial concentrations of the primary enzyme, pyruvate, or the substrate, as this will cause a very large lag time. One can use (33) to specify a required lag time, and determine constraints for the initial concentrations discussed above. However, while we have identified and characterized the linear regime of NADH production, this regime is unfortunately transient; it can only last until substrate depletion starts to affect the system at leading order, which occurs when t becomes of O(1/ ). We investigate this regime in the next subsection.

Depletion regime
To investigate the depletion regime, we introduce T = t and analyse the limit → 0 with T = O(1). Moreover, from (24) and (25), we see that we must also make the scalings S 4 = s 4 / , and S 6 = s 6 / , resulting in the system where the three reaction velocities are now defined as The matching conditions from the earlier timescale yield the following 'initial conditions' Naively taking the limits of , δ, K → 0 in (34) to obtain a leading-order system would yield a duplication of information. To avoid this, we must form appropriate linear combinations of the governing equations in order to obtain enough information for a leading-order system. Using this approach, exploiting conserved quantities where possible, we obtain the following differential-algebraic system at leading-order s 5 + cS 6 = 1.
Thus, we have reduced the problem to that of one ODE with four additional algebraic relationships, and the initial conditions (36). To further simplify this system, it is convenient to use (37b-d) to write then to transform (37a) into a separable ODE for s 1 , obtaining Solving (39) allows the remaining variables to be deduced algebraically from (37b-e), (38). We are able to integrate (39) to obtain the following implicit solution where The term with coefficient bK (1) Thus, we have characterized the depletion regime. One can use (37c-d), and (40) to obtain the following functional forms for the chemical concentrations in the depletion regime The asymptotic solution for S 6 in (41) agrees well with the numerical results for a range of values of k (−2) in both the linear-phase regime of the previous section and the depletion regime of this section (Fig. 8).
While the transient-linear regime of the previous section is much easier to fit and has the benefit of being easier to identify from time data, this linear behaviour may only be exhibited for a short duration. Given experimental data on the NADH concentration, we can therefore estimate the three unknown parameters k (1) , K

, and K
(1) 2 via a nonlinear parameter fitting of (41). In practice, it is likely that either s 1 or s 2 will be saturating, corresponding to either K (1) 1 or K (1) 2 being small, thus reducing the number of parameters that need to be fitted at any one time. It is worth noting that, unlike in the standard assay considered in Sect. 2, the form of S 6 is not symmetric in the parameters K (1)

and K
(1) 2 . This is to be expected, since the symmetry in the reaction network is broken in the non-standard assay.
To identify when we may be able to use the transient-linear regime to fit data, we now estimate an upper bound for this regime. While it is possible to determine this by obtaining an early-time solution to the system (37) to calculate a correction to the long-time solution from the earlier time regime (24), the resultant expression is cumbersome and not particularly insightful. Instead, we use the following more crude approximation for an upper limit to the linear regime: obtained from the requirement to stray from linearity in this depletion regime. In terms of dimensional quantities, (42) Over the timescale t = O(1/a ), the levels of S 6 increase and the levels of s 2 and s 5 decrease. From (41), we see that the maximum amount of NADH produced is bounded above by an O(1) quantity. Hence, we see from (34)-(35) that there will be additional depletion dynamics if s 5 = O(δ). This condition is not guaranteed to occur, and will depend on the values of the kinetic parameters, as well as on the initial conditions of the system. We do not consider this depletion possibility further, as it has no bearing on our main goal of characterizing the primary enzyme. We consider the remaining distinguished limit in the system, that of a weak reverse reaction, in the next section.

Linear-phase regime
The distinguished limit where the reverse reaction is weak occurs when k (−2) = O( K ), and it will be helpful to introduce ρ := k (−2) /( K ) = O(1) for notational purposes. This scaling results in a distinguished limit in which the nonlinearities of the primary reaction and the reverse auxiliary reaction are important when t = O(1/ ). This causes the limit to be distinguished for s 4 over a timescale of t = O(1/ ).
In this case, the leading-order version of (15) for t with initial conditions (17). The solution to (44) is where V is defined in (20). Hence, s 6 initially scales with t 2 before eventually scaling with t. Thus, when t = O(1) for a weak reverse reaction, the observed NADH production is essentially equivalent to that of Sect. 2, even though pyruvate is both a substrate of the primary reaction and a product of the auxiliary reaction. This is because there is not enough alpha-L-alanine being created to produce a large enough amount of pyruvate to affect the NADH production, and the reverse reaction is too weak to siphon away the NADH that is produced by a significant amount. In dimensional terms, the observed long-time concentration of NADH in this regime is essentially the sub-limit of (31) as ω → 1, withṼ and ω defined in (32). Additionally, we are able to obtain a lower bound for this linear regime from (45), namely that t 1 for the linear regime to hold. In dimensional terms, this corresponds tõ

Depletion regime
As t increases further, the substrate depletion and the reversible reaction from NADH to NAD + start to affect the problem at leading order. To investigate this, we make the scalings t = T / and S 6 = s 6 / . Then, the leading-order version of (15) is While (48) may appear complicated at first glance, it reduces readily to allow an implicit analytic solution for the system. To see this, we first use (48c) to reduce (48a) and deduce that at leading order over this timescale. Using (49a) in (48b) yields a separable ODE for s 2 , with implicit solution which can also be written in terms of the Lambert W function. This tells us that the leading-order substrate concentration will start to decrease before the pyruvate (Fig. 9), even if the substrate is saturating i.e. K (1) 2 → 0. Finally, forming linear combinations of (48b)-(48e) to yield conserved quantities allows us to deduce that From our analytic solutions, we are able to infer that the weak reverse reaction only affects the amount of s 4 in the system over this timescale; the other chemicals are independent of the reverse reaction despite it appearing at leading order. This change of scaling affects the asymptotic analysis of the depletion dynamics between the strong and weak reverse reaction cases; care is needed in scaling appropriately if taking k (−2) → 0 in (37)-(41). We are able to deduce the long-time dynamics in terms of the implicit solution for s 2 given in (49b). As the experimental calculation for the reaction velocity is taken when the growth of NADH is linear, it is helpful to note when and how this breaks down in this long-time limit. To do this, we expand (49b) for small time to note that Moving back to dimensional quantities, we find that Fig. 9 A comparison of the substrate depletion between the strong and weak reverse reactions, obtained from numerical solutions to the full system (15). We use the parameter values in Table 3, withα =β = 2.5 mM, k (−2) = K = 5 × 10 −4 for the weak reverse reaction and k (−2) = 1/K = 20 for the strong reverse reaction (52) The simplest way to make this constraint less limiting is to have as small a value of =k (1) /k (2) as possible. Another way is to increaseβ, though we note that this also has the effect of making it difficult to determineK (1) 2 . Therefore, to minimize the correction term (and hence increase the duration of time at which the NADH growth is linear, thus reducing experimental error), we should decreaseα and increaseβ. The former of these is perhaps a little surprising, though we note that we must also ensure that pyruvate depletion does not affect anything at leading order, soα K (2) 4 ≈ 0.05 mM must also hold.
There will be further asymptotic regions and depletion dynamics at play here at longer timescales, such as if s 5 becomes of O(δ). However, these additional regions are of limited relevance to our analysis here.

Guide to using our theoretical results
We now summarize the main results from Sects. 3.1 and 3.2, and give a step-by-step guide of how to use these to characterize enzymes. We have found that a linear-growth regime for NADH is possible for a strong reverse reaction, but this regime can be quite short and is sandwiched between two nonlinear-growth regimes. The NADH concentration in the linear regime is given in (31) in terms of dimensional quantities. Although the linear-growth regime in the weak reverse reaction case is also sandwiched between two nonlinear-growth regimes, it lasts for a much longer duration in this case and so results from the linear regime are easier to use in practice.
Due to our analytic asymptotic results, we are able to deduce upper and lower time constraints for the linear regimes, which we give in (33) and (43) in terms of dimensional quantities for the strong reverse reaction case and in (52) for the weak reverse reaction case. We are also able to obtain an implicit closed-form solution for the NADH concentration away from this linear regime, which would allow the unknown quantities to be inferred through a nonlinear parameter fitting, if the linear regime was too short for accurate measurements to be obtained. We now provide a step-by-step guide to using our results with time data for the concentration of the indicator chemical in each case, using different values of the initial substrate concentrationsα andβ. We also present a summary of these guidelines as a flow chart in Fig. 10, and note that a summary of the solutions we derive in this paper in provided in Table 4.
Step 1 Characterize the auxiliary enzyme.

Using the weak reverse reaction results
In this case, the NADH concentration should exhibit linear growth for a significant period of time, with the relevant bounds given in (52). We are able to reduce the lower bound in (52) by adding more of the auxiliary enzyme, thereby reducing the lag period. Additionally, although the upper bound in (52) contains the parameters we are trying to infer, we can see that adding less of the primary enzyme will increase it, extending the length of time until the depletion dynamics kick in.
From the experimental data in the linear growth regime, the procedure is to calculate the average slope of NADH concentration. This slope is equal toṼ , defined in (32). The kinetic parameters of the primary enzyme can be inferred by fitting the average slope obtained from the experimental data to the expression we give forṼ in (32), using nonlinear data-fitting to minimize the sum of the difference between the data and functional form.

Using the strong reverse reaction results
In this case, the NADH concentration should exhibit linear growth for some period of time, with the relevant bounds given in (33) and (43). While the upper bound of (43) can be extended by adding less of the primary enzyme, it is not necessarily true in this case that adding more of the auxiliary enzyme will decrease the lower bound of (33). If the linear regime is able to be observed for an appropriate length of time, then the appropriate approach is to calculate the average slope of NADH concentration from the experimental data in the linear growth regime. This slope is equal toṼ ω, defined in (31). Here, ω is a function of the kinetic parameters of the auxiliary enzyme and the initial chemical concentrations, whereasṼ depends on the kinetic parameters of the primary enzyme and the initial chemical concentrations. Hence, the kinetic parameters of the primary enzyme can be inferred by fitting the average slope obtained from the experimental data to the expression we give forṼ ω in (32), using nonlinear data-fitting to minimize the distance between the data and functional form.
If the linear regime is not able to be observed for an appropriate length of time, there are two options, depending on whether better data is available fort = O(K (2) 4 /k (2) ) or fort = O(α/k (1) ). In the former case, one should use the results from Sect. 3.1, where the data-fitting has been reduced to the single ODE (21). In terms of dimensional quantities, the ODE (21) can be re-written as +s 6 ,s 6 (0) = 0.
While there are four different parameter groupings in (53), the kinetic parameters of the primary enzyme are all contained withinṼ , as described in (32). Therefore, the data fitting to this ODE is likely to be significantly simpler than data-fitting to the initial system of ODEs (1). If better data is instead available fort = O(α/k (1) ), then the results of Sect. 3.1.3 should be used. Here, the functional form for the NADH concentration is given by (38), (40), and (41). In this case, the data-fitting should involve minimizing the sum of the difference between the time data and the functional form.

Comparison with a naive nonlinear solver
We now compare the asymptotic method we describe in this paper with a naive nonlinear parameter fitting. To do this, we first generate 'true' time courses for NADH from the system (10)-(11) for four different values ofα, then add Gaussian noise with mean 0 µM and standard deviation σ µM, and we use σ as a control parameter. An example of these time courses with σ = 1 is given in Fig. 11. The asymptotic parameter fitting is carried out by fitting the asymptotic solution (31) to the true time courses using the MATLAB function 'lsqcurvefit' for the parametersk (1) andK (1) 1 . The naive nonlinear parameter fitting is carried out using the same MATLAB function 'lsqcurvefit' on the ODE system (10)-(11), fitting the parameter valuesk (1) ,K (1)

, andK
(1) 2 . However, as we generate the true time courses using a large value ofβ in order to focus onK (1) 1 rather thanK (1) 2 , we ignore the value ofK (1) 2 generated from this. In general, the naive nonlinear fit is worse than the asymptotic fit in determining k (1) and better than the asymptotic fit in determiningK (1) 1 (Fig. 12). However, the asymptotic fit appears to be more consistent in its predictions than the naive nonlinear fit, as the naive nonlinear fit can predict wildly incorrect values ofk (1) even for small noise. Moreover, we note that the computational time required for the naive nonlinear fit is around 1000 times larger than for the asymptotic method.

Application to experimental data
The above concludes our mathematical analysis of the problem. To showcase how these results can be used in practice, we performed enzyme assays to characterize an (a) (b) Fig. 12 A comparison between the asymptotic method presented in this paper, and a naive nonlinear fit of the system (10)-(11). We show the root mean squared error (RMSE) for both methods in estimating the parameters ak (1) and bK (1) 1 . The dashed line shows the value of the parameter being estimated, so a RMSE well below this line signifies good accuracy whereas a RMSE near or above the line signifies poor accuracy. We generate the 'true' data from the system (10)-(11), then add Gaussian noise with mean 0 µM and standard deviation indicated on the x-axis. Each marker represents the RMSE calculated from 100 total realisations. We use the parameter values from Table 2, as well ask (1) = 0.02 µM s −1 ,k (2) = 10 µM s −1 , k (−2) = 10 µM s −1K (1) 1 = 1.5 mM,K (1) 2 = 1 mM,β = 50 mM, andγ = 50 mM. We also useα = 0.5, 1, 2.5, 5 mM to generate four different 'true' time course curves for each different standard deviation of noise enzyme for two different substrates using a previously uncharacterized putative omegaaminotransferase, cloned and purified in a heterologous Escherichia coli host. The aminotransferase is CnAptA (UniProtKB Q0KEZ8), named here for the first time in the literature and identified in the genome of Cupriavidus necator H16, an industrially relevant, facultative chemolithoautotrophic chassis microorganism of the Synthetic Biology Research Centre (SBRC-Nottingham). We describe the experimental protocol used to purify this enzyme in "Appendix B".

Enzyme assay protocol
For these assays, the corresponding reaction network is shown in Fig. 1. The primary enzyme is CnAptA, which converts pyruvate and an additional omega-amino acid substrate into a corresponding aldehyde and alpha-L-alanine. Our assay solution also included beta-nicotinamide adenine dinucleotide sodium salt (NAD + ) as an electron acceptor for the auxiliary enzyme reaction, and we monitored the formation of NADH over time, produced from the auxiliary reaction. We used a commercially available alanine dehydrogenase (A-DH) from Bacillus cereus (alanine dehydrogenase, recombinant; CAS-No 9029-06-5, Sigma-Aldrich Company Ltd.) as the auxiliary enzyme, and we performed the experiment at pH 10, the pH condition recommended by the supplier of this auxiliary enzyme, this condition resulting in an undetectable reverse reaction. Our analysis suggests that the primary enzyme will be easier to characterize, since the bounds on the validity of the linear regime are less severe. Hence, we can either use the sub-limit as ω → 1 of the results in the distinguished limit (31), or we can use the sub-limit directly deduced in Sect. 3.2. The approach for this is Performing these assays at different acidities, which may be more relevant for certain bacterial cell environments, would result in a strong reverse reaction, which could be dealt with using the analytic representation for the distinguished limit given in (31). The approach for these cases is summarised in Sects. 4 and 4.2.
We performed the enzyme assays in 96 well microtiter flat bottom transparent plates (Costar ® Sigma-Aldrich). The reaction mixture contained different concentrations of the substrate, 2.5 mM NAD + , 5 mM pyruvate, 5 µl A-DH [0.28 mg/ml], 5 µM of pyridoxal 5 -phosphate (PLP), and 5 µl of CnAptA [0.0245 mg/ml]. We then used an appropriate amount of sodium carbonate buffer to obtain a final volume of 200 µl. The two different substrates we used were 3-aminobutanoic acid (Sigma-Aldrich 97%, CAS No 541-46-6) and 5-aminovalerate (ACROS Organics,. To detect the NADH present in the solution, we measured the light absorbance at 340 nm using a TECAN Infinite ® M1000 PRO plate reader, at 30 • C, the absorbance being linearly proportional to the concentration of NADH in the solution. We monitored the light absorbance at 340 nm at 60 s intervals for 1 h. After subtracting off the base level of absorbance at the start of each experiment and the absorbance increase monitored in the control experiment with no substrate, we converted these clean absorbance data into an NADH concentration, shown for each substrate in Figs. 13a and 14a. We expect our TECAN data to be accurate to within a few µM, and this error margin is why there sometimes appears to be 'negative' concentrations for early time. As our analysis predicts, there is a clear linear regime after a lag period. We obtained a constant reaction velocity from this linear regime, using the average slope between 1800 and 3600 s. From this constant reaction velocity, we used a nonlinear data fit to the functional form of the solution (31)-(32) with ω → 1 (as discussed in Sects. 4 and 4.1) by minimizing the sum of the difference between the data and functional form squared, using the inbuilt function fminsearch in MATLAB, a simplex search method for finding minima (Lagarias et al. 1998). From Figs. 13b and 14b, we see that our nonlinear fits generally provide accurate representations of the data. For the enzyme CnAptA, our nonlinear fit predicts values ofk (1) = 0.026 µM s −1 andK (1) 2 = 1.3 mM for the substrate 3aminobutanoic acid and values ofk (1) = 0.021 µM s −1 andK (1) 2 = 0.93 mM for the substrate 5-aminovalerate, assuming that the levels of pyruvate we use are saturating. In terms of the catalytic efficiencies, we see thatk (1) /K (1) 2 = 2.0 × 10 −5 s −1 for 3aminobutanoic acid andk (1) /K (1) 2 = 2.3 × 10 −5 s −1 for 5-aminovalerate. Therefore, we note that while 5-aminovalerate has a higher affinity (a lowerK (1) 2 ) to and catalytic efficiency for CnAptA than 3-aminobutanoic acid, the maximum reaction velocity is higher for 3-aminobutanoic acid than for 5-aminovalerate. However, these differences are fairly small -there is not a large difference in enzyme behaviour between the two substrates.

Discussion
We have used singular perturbation theory to determine how to infer the kinetic parameters of the primary enzyme from the measured NADH concentration in a coupled enzyme assay, with reaction network shown in Fig. 1. In particular, this will allow putative aminotransferases to be quickly characterized using coupled enzyme assays. Our analysis also allows us to obtain estimates for when we expect to see a linear growth in NADH concentration, the typical regime measured in enzyme assays since it allows a quick sanity check of experimental results and is much easier to fit to data. We have outlined how to use our theoretical results to infer the kinetic parameters of the primary enzyme from experimental data in Sect. 4, where we also compared our method with a a naive nonlinear fit. Moreover, we have used our results to characterize a novel putative aminotransferase for two different substrates in Sect. 5.
Our analysis shows that there are two distinguished asymptotic limits in the system: the important one for NADH (the measurable chemical) being the strong reverse reaction case, considered in Sect. 3.1; the other is the weak reverse reaction case, considered in Sect. 3.2. We demonstrate that all the information that can be inferred from the NADH concentration in the weak reverse reaction case is a sub-limit of the strong reaction case, even though there are differences in the concentrations of alpha-L-alanine, an intermediate chemical, between these two cases. In our analysis, we exploit the smallness of three dimensionless parameters: , δ, and K , with a focus on analysing how the asymptotic size of k (−2) affects the system behaviour and treating all other parameters in the system as O(1) to retain their generality. As the parameters and δ can always be made to be small through experimental choice, and we also present results for K = O(1) in "Appendix C" (which shows a slight change in the depletion regimes), our work represents a comprehensive asymptotic analysis of the possible behaviours of the non-standard coupled enzyme assay we consider in this paper.
Additionally, we note that our results suggest that one cannot deduce whether the strong reaction case pertains simply by looking for a change in reaction velocity when the amount of auxiliary enzyme is varied. This can be seen by (31), where we see that the only two variables that involve the experimentally controlled auxiliary enzyme concentration,k (2) andk (−2) , appear as a ratio, thus removing any dependence on this enzyme concentration. Thus, it is not always immediately clear from observations whether one is in the weak or strong reverse reaction regime, unless the auxiliary reaction has previously been completely characterized. Due to the weak reverse reaction being a sub-limit of the strong for the experimentally measurable NADH concentration, as discussed above, we advise working with the results derived in the strong reverse reaction regime in Sect. 3.1.
In the derivation of our model, we assume that the kinetic parameters of the putative enzyme are unchanging over time. This may not be valid for enzymes that are unstable or suffer a loss of activity over time, for example through incubation in a nonoptimal buffer at certain temperatures. Accounting for this change in enzyme activity would require additional information about the type and effect of the degradation occurring.
Although it would be possible to perform a brute force fitting of the kinetic parameters of the primary enzyme using the governing equations (10), such an approach would be time-consuming and would not provide the general analytic expressions we derive through our asymptotic analysis. Moreover, we are able to use our analytic expressions to estimate when we will be in the linear growth regime, and for how long we expect this regime to last, in terms of the experimental parameters. This provides a useful sanity check when it comes to interpreting the data.
Historically, the measured reaction velocity is often converted into kinetic parameters for the enzyme using a Lineweaver-Burke plot (Lineweaver and Burk 1934). However, this procedure is generally not the best way to calculate these values as it amplifies noise in the measurements (Fell 1997). Nonlinear data fitting is orders of magnitude better than it was at the time of Lineweaver and Burk, so it is much more accurate nowadays to fit the data to analytic forms such as those we have deduced using singular perturbation theory. For the case we consider in this paper, measurements with varying initial substrate concentrations give us a two-dimensional array of data from which to fit three variables. Additionally, choosing one of these substrate concentrations to be much larger than the estimated Michaelis constantK (1)

orK
(1) 2 reduces the problem to a one-dimensional array of data from which to fit two variables.
While we have considered a specific type of non-standard coupled enzyme assay here, the asymptotic techniques we use are far from restricted to this system. We hope that the techniques used in this paper will be applied to different reaction networks, to help characterize enzymes that can only be assayed in non-standard reaction network topologies.

A Distinguished asymptotic limits of equation (26)
In this Appendix, we investigate the distinguished asymptotic limits of (26), thus allowing us to characterize the possible behaviours of its solution.

A.1 Case I
The first distinguished limit we discuss (Region I in Fig. 7) occurs whenV 1 and γ = O(1) for the ODE (26). Here, there are two different timescales in the problem: t = O(1) and t = O(1/V ), for which Y = O(V ) and Y = O(1), respectively; the limit is distinguished only over the latter. Therefore, our brief analysis of the earlier timescale consists simply of the asymptotic solution For the later timescale, the leading-order governing equation is and hence we are in the full quasi-steady case. The equation (55) has solution which matches appropriately to the earlier-time solution. We note that the large-time limit of (56) is Y ∼V t/(γ + 1). This is different to the large-time limit of (54) (when t = O(1)), which is Y ∼V t. Therefore, there are two different linear regimes here, and one should take great care in this distinguished limit to ensure that one is aware of in which linear regime the measurements are being taken: failure to do so could result in a significant error of interpretation. When γ → ∞ in this distinguished limit, we move into sub-limit A (Fig. 7). Here, the earliest timescale is the same as (54), but the later timescale splits into two, the earlier with t = O(1/(V γ )) and Y = O(1/γ ), the later with t = O(γ /V ) and Y = O(1), as could be inferred from the large-γ limit of (56). On the timescale and on the timescale When γ → 0 in this distinguished limit, we move into sub-limit D (Fig. 7). Here, the full leading-order solution in the later timescale is simply the large-time limit of the earlier timescale. Hence, the leading-order solution is simply (54), uniformly valid for all time.

A.2 Case II
The next distinguished limit (Region II in Fig. 7) occurs whenV 1 and γ 1, with V γ = O(1), in the ODE (26). Here, there are two different timescales in the problem: t = O(1) and t = O(γ 2 ), for which Y = O(1/γ ) and Y = O(1), respectively; the limit is distinguished only over the earlier of these. Taking Y =Ȳ /γ , the leading-order equation for the earlier timescale is The Riccati equation (59) has solution where Ai(z) and Bi(z) are the standard Airy functions. The solution (60a) exhibits the following large-t behaviourȲ ∼ V γ t as t → ∞.
The later timescale in this problem occurs when t = O(γ 2 ) = O(1/V 2 ) and Y = O(1). From (26), this timescale yields the following leading-order equation with solution This exhibits the large-t behaviour Hence, we may deduce that in this regime the variable Y initially behaves as Y ∼ V t 2 /2, then as Y ∼ V t/γ , then finally as Y ∼V t/γ , and therefore exhibits at least two points of inflection, as well as three distinct power laws. The solutions of the distinguished limit in Case II tend to those of the sub-limit A (Fig. 7) asV γ → 0. The early timescale, when t = O(1) with solution (60a), splits into two, with solutions (54) and (57). The late timescale solution (63) is the same as in sub-limit A, as given in (58).

A.4 Case IV
The final distinguished limit (Region IV in Fig. 7) occurs whenV 1 and γ = O(1) in the ODE (26). Here, the main balance occurs when t = O(1) and Y = O(V ). In this case, the leading-order version of the ODE (26) is which has solution In the limit of γ → ∞, the solution (70) reduces to that in sub-limit C (Fig. 7), given by (68), with the important timescale t = O(1/γ ). In the limit of γ → 0, the solution (70) reduces to that of sub-limit D (Fig. 7), given by (54). This completes our comprehensive asymptotic analysis of the distinguished limits of the ODE (26).

C Case in which K = O(1)
In this Appendix, we provide the governing system in the depletion regime for the case with a strong reverse reaction when K = O(1), as the t = O(1) case proceeds in the same manner as when K is small. Here, substrate depletion starts to affect the system at leading order when t becomes of O(1/ ). To investigate this, we introduce T = t, and we must also make the scalings S 4 = s 4 / , and S 6 = s 6 / . These scalings result in the system where the three reaction velocities are now defined as v 1 = s 1
Naively taking the limit as → 0 in (34) to obtain a leading-order system would yield a duplication of information. To avoid this, we must form appropriate linear combinations of the governing equations in order to obtain enough information for a leading-order system. Using conserved quantities where possible, we obtain the following differential-algebraic system at leading-order s 5 + cS 6 = 1.
It is still possible to follow the analysis of Sect. 3.1 and rearrange this system into a single ODE, though this is too cumbersome to be particularly enlightening.