Critical Parameters for Singular Perturbation Reductions of Chemical Reaction Networks

We are concerned with polynomial ordinary differential systems that arise from modelling chemical reaction networks. For such systems, which may be of high dimension and may depend on many parameters, it is frequently of interest to obtain a reduction of dimension in certain parameter ranges. Singular perturbation theory, as initiated by Tikhonov and Fenichel, provides a path towards such reductions. In the present paper, we discuss parameter values that lead to singular perturbation reductions (so-called Tikhonov–Fenichel parameter values, or TFPVs). An algorithmic approach is known, but it is feasible for small dimensions only. Here, we characterize conditions for classes of reaction networks for which TFPVs arise by turning off reactions (by setting rate parameters to zero) or by removing certain species (which relates to the classical quasi-steady state approach to model reduction). In particular, we obtain deﬁnitive results for the class of complex-balanced reaction networks (of deﬁciency zero) and ﬁrst-order reaction networks.


Introduction
The modelling of chemical reaction networks frequently leads to high-dimensional parameter dependent systems of ordinary differential equations (ODEs).Even in the presence of a well-established structure theory for large classes of reaction networks, reducing the dimension of such systems is desirable for several reasons: From a quantitative perspective in the laboratory, parameter identification is frequently unfeasible for the full system but might be possible for a reduced equation.The Michaelis-Menten system and generalizations can be seen as examples of this; see e.g.Segel and Slemrod [42], Keener and Sneyd [32].From a qualitative vantage point, one strategy to prove special features such as the existence of periodic solutions, or multistationarity, is to prove such features for a reduced system and show that they persist for the full system in some parameter range.For a recent example of this strategy, see [17].Thus it is of general interest to identify parameter domains where a systematic reduction is possible.
Typically (although not exclusively) the reduction procedures are based on singular perturbation theory as developed by Tikhonov [43] and Fenichel [21].In the present paper, we will discuss singular perturbation reductions and critical parameters that permit reductions of this kind.The focus will be on characterizing such critical parameters that correspond naturally to structural features of the chemical reaction network.
A frequently used approach to finding appropriate parameters for singular perturbation scenarios goes back to a classical paper by Heineken et al. [28].The method relies on an adroit scaling of suitable variables (based on an intuitive understanding of the processes in the reaction network) and ideally leads to a system with slow and fast variables to which Tikhonov's and Fenichel's theorems are applicable.From another perspective, a singular perturbation approach for systems with prescribed slow and fast reactions was discussed by Schauer and Heinrich [40].More recently, a complete characterization of the parameter values (called Tikhonov-Fenichel parameter values, briefly TFPVs) which give rise to singular perturbations, and of their critical manifolds, was obtained in A. Goeke's dissertation [23] and the ensuing papers [25][26][27] by Goeke et al.Moreover for polynomial or rational systems, an algorithmic path exists toward determining these parameter values.The theory was applied to a number of reaction networks, including standard reaction networks from biochemistry [32], and for these all possible singular perturbation reductions could be determined.In addition, it turned out that the algorithmically determined TFPVs for these systems readily admit an interpretation in terms of chemical species concentrations and reaction rates: Frequently these TFPVs correspond to a "switching off" of certain reactions, or a removal of certain chemical species.This is the vantage point for the present paper.Since there is a natural limit to any algorithmic approach for systems with large numbers of variables or parameters, generalizing such structural insights is of interest.
From a mathematical as well as from a chemical perspective it seems desirable to understand whether (and how) special properties of reaction networks imply the existence of particular classes of TFPVs.The purpose of the present paper is to contribute toward this understanding.We will focus on reaction networks with mass-action kinetics, hence on polynomial differential equations.Our goal is to employ the structure of chemical reaction networks to obtain heuristics for finding TFPV candidates (respectively, candidates for scaling) in a first step and then, in a second step, proceed to verify the TFPV property for some reasonably large and relevant classes of reaction networks.We make substantial use of the structure theory going back to Horn and Jackson [29], Feinberg [12] and others.In terms of chemical reaction networks, we are concerned with slow and fast reactions, on the one hand.On the other hand, we investigate the provenance of quasi-steady state phenomena for chemical species, and their naturally associated "slow-fast" systems.Our main results apply in particular to weakly reversible reaction networks of deficiency zero.
Specifically, in Section 3 we first consider TFPVs that arise from turning off reactions, and identify graphical means for their identification (Theorems 3.11 and 3.12).We also provide an explanation of why TFPVs in many cases belong to proper coordinate subspaces (Proposition 3.3).In particular, we obtain a complete characterization for weakly reversible systems of deficiency zero.Continuing, in Section 4, we characterise sets of species (so-called LTC species sets) that "shut down" the reaction network when the corresponding variables are zero (hence, the species are present in zero concentration).Such species sets naturally lead to slow-fast systems (in a weak sense), and we further investigate their relation to linear first integrals and give conditions for when an LTC species set is the support of a linear first integral (Proposition 4.7).We proceed to discuss conditions for TFPV for systems on stoichiometric compatibility classes (Proposition 4.11 and its corollaries).Finally, we briefly consider combining the approaches to turn off certain reactions and to remove certain species.
The paper is organized as follows.Section 2 contains preliminaries on reaction networks, TFPVs and (in a weak sense, formally) slow-fast dynamical systems.Section 3, in the context of reaction networks, discusses TFPVs defined by rate parameters.Section 4 builds on Section 2 and connects results of Section 3 to the classical scaling approach and slow-fast systems.The results are illustrated by examples.In particular, we make use of many standard textbook reaction networks, such as the Michaelis-Menten reaction network (albeit its reversible version) that serves as benchmarking example.

Preliminaries
We let R, R ≥0 , R >0 denote the sets of real, non-negative real and positive real numbers, respectively.Also, we let N 0 denote the set of non-negative integers.Given m ∈ N 0 , a coordinate subspace of R m is defined by The support supp(x) of x ∈ R m is the set of all indices i with x i = 0, i = 1, . . ., n.For y = (y 1 , . . ., y n ) ∈ N n 0 and x = (x 1 , . . ., x n ) ∈ R n ≥0 , we define x y = n i=1 x y i i .If M = (m 1 . . .m k ), m i ∈ N n 0 , i = 1, . . ., k, is an (n × k)-matrix, then we define x M as the vector (x m 1 , . . ., x m k ) ∈ R k ≥0 .

Reaction networks
We consider spatially homogeneous chemical reaction networks with constant thermodynamical parameters and kinetics of mass-action type.The mathematical theory of these reaction networks was initiated and developed in seminal work by Horn and Jackson [29], and Feinberg [12].We will refer to Feinberg's recent monograph [13] as a basic source.First we introduce the notion of a reaction network and fix some terminology.
Definition 2.1.A mass-action reaction network over a set of species consists of non-negative integer linear combinations in X , and κ labels edges by positive real numbers.Isolated nodes, but not self-edges, are allowed.We refer to the nodes as complexes, to the edges as reactions, and to the labels as rate parameters.Every species is assumed to be in some complex with a positive coefficient.Throughout we let d be the cardinality of Y and m the cardinality of R.
A reaction network G is a subnetwork of another reaction network G with species set X , if G is a subdigraph of G.
We enumerate the set of complexes in some way, and thus write Y j = n i=1 y ij X i with y ij ∈ N 0 .The y ij 's are referred to as stochiometric coefficients.A labelled reaction between the complexes Y j , Y ℓ is written as Here Y j is called a reactant complex and Y ℓ a product complex.Note the reversal of the subindex of κ in the labels.A numbering of the elements of R by 1, . . ., m, provides an ordering of R and we identify the collection of κ ℓj with a vector κ ∈ R m >0 , ordered in the same way as R, such that is the i-th reaction.We will use this convention without further reference.The zero complex 0 is allowed by definition.Reactions with reactant 0 are called inflow reactions, and account for production or influx of species.
As a reaction network is given as a directed graph, terminology and properties from graph theory apply.Moreover, special terminology has been developed, parallel to terminology in graph theory.We will refer to a reaction network where all connected components of the digraph are strongly connected as weakly reversible, and otherwise apply standard terminology.
The evolution of the species concentrations in time is modelled by means of a system of ODEs, assuming mass-action kinetics.Denote by x(t) = (x 1 (t), . . ., x n (t)) the vector of concentrations of the species X 1 , . . ., X n at time t.Define the complex matrix by consisting of the stoichiometric coefficients of the complexes, and let y 1 , . . ., y d denote its columns.We let B be the reactant matrix with i-th column y j if Y j is the reactant of the i-th reaction, and N ∈ R n×m the matrix, referred to as the stoichiometric matrix, with i-th column given by y ℓ − y j if Y j −−→ Y ℓ is the i-th reaction.With this notation, the system of ODEs becomes: where reference to t is omitted and κ ∈ R m >0 .The sets R n >0 and R n ≥0 are positively invariant for (1) [45].Furthermore, there is a useful decomposition of the right-hand side of (1) in terms of the Laplacian of the reaction network.The Laplacian matrix A(κ) = (a ij ) 1≤,i,j≤d ∈ R d×d is given by where (See also Feinberg [13,Subsection 16.1] for further background.) System (1) often admits stoichiometric first integrals.These are non-zero linear forms Note that α 1 , . . ., α n might be chosen as integers, since N has integer entries.
Definition 2.2.The image of the stoichiometric matrix N is the stoichiometric subspace, and the intersection of every coset of this subspace with the non-negative orthant is a stoichiometric compatibility class (SCC).
The dimension (respectively, codimension) of the mass-action reaction network is by definition the dimension (respectively, codimension) of the stoichiometric subspace.In principle, system (1) might admit further linear first integrals.However, the following result says that it does not happen for realistic networks.

Lemma 2.3 ([14]
).If every connected component of a reaction network has exactly one terminal strongly connected component, then every linear first integral of (1) is stoichiometric.
Example 2.4.Consider the mass-action reaction network with species X 1 , X 2 and reactions The corresponding ODE system is given by ẋ1 The reaction network has no linear first integrals for generic κ, but when κ 1 = κ 2 , the vector (1, 1) defines one.This reaction network has one connected component, but two terminal strongly connected components, namely {0} and {2X 1 }.
−−→ X 3 , the dimension is 1, but for all κ there are two linearly independent linear first integrals: a stoichiometric linear first integral φ 1 = x 1 + x 2 + x 3 , and a non-stoichiometric, While it is possible at the outset to reduce the dimension of system (1) via these first integrals, for the purpose of the present paper it seems appropriate to keep the representation (1) until at a later stage.

Tikhonov-Fenichel parameter values (TFPVs)
Throughout, when refering to singular perturbation reduction, we mean this in the sense of Tikhonov [43] and Fenichel [21].In order to identify parameters that give rise to singular perturbation reductions the following approach was taken in Goeke's dissertation [23] and the subsequent papers [25,26].
Consider a parameter-dependent ODE system, with h(x, π) polynomial in x and π.We let D 1 h(x, π) and D 2 h(x, π) denote the partial derivatives with respect to x and π, respectively.Given π ∈ Π, we denote by V(h(•, π)) the zero set of x → h(x, π), and let n − s * be the generic dimension (with respect to π) of the vector subspace generated by the entries of h(x, π) for x ∈ Ω.In addition, we require that the generic rank of In the setting of mass-action reaction networks, this subspace is equal to the stoichiometric subspace under the hypotheses of Lemma 2.3.In this case, s * is the codimension of the reaction network, according to Definition 2.2.
The existence of singular perturbation reductions coincides with the existence of Tikhonov-Fenichel parameter values (TFPVs).Definition 2.5.A TFPV for dimension s (s * < s < n) of system (4) is a parameter π ∈ Π, such that the following hold:

π).
(iii) There exists x 0 ∈ Z such that all non-zero eigenvalues of D 1 h(x 0 , π) have negative real part.We let Π s ⊆ Π denote the set of TFPVs for dimension s > s * .
Note that provided (i) holds, then (ii) and (iii) are together equivalent to (ii') There exists x 0 ∈ Z such that D 1 h(x 0 , π) has exactly n − s non-zero eigenvalues (counted with multiplicity), which additionally have negative real part.
The conditions imply that the critical manifold Z is locally exponentially attracting.We have the following characterization [23,26].
Proposition 2.6.Given a parameter π ∈ Π and any smooth curve ε → ϕ(ε) in the parameter space Π with ϕ(0) = π, the system admits a singular perturbation reduction in the sense of Tikhonov and Fenichel if and only if π is a TFPV.
Thus one may think of a TFPV as a ("degenerate") parameter set from which singularly perturbed systems emanate.
If Π and Ω are semi-algebraic sets, then Π s is a semi-algebraic set as well [26].In any case, the Zariski closure of Π s exists and is denoted by (5) W s := Π s Zar .
An alternative characterization of TFPVs and the basis for an algorithmic approach to TFPVs is the following [23,26].Consider the characteristic polynomial of D 1 h(x, π).Then, given s * < s < n, a parameter value π is a TFPV with locally exponentially attracting critical manifold Z (depending on π) of dimension s, if and only if the following hold for some x 0 ∈ Z: (2) all roots of χ(τ, x 0 , π)/τ s have negative real part.
(vi) The system ẋ = h(x, π) admits s independent local analytic first integrals at x 0 .Therefore, a starting point for computing TFPVs is as follows: With h(x 0 , π) = 0 and σ n (x 0 , π) = • • • = σ n−s+1 (x 0 , π) = 0, one sees that (x 0 , π) is a solution to n + s > n equations for x ∈ R n , given π.In turn, this allows to obtain conditions on π for general polynomial systems via elimination theory.
The validity of the hypotheses for Tikhonov's and Fenichel's theorems depend on the ambient space, and thus may change when passing to an invariant subspace.As a consequence, the notion of TFPV may also depend on the ambient space.For reaction networks this observation is relevant when passing to SCCs.

Slow-fast systems and scalings
In the present paper, we call a smooth system of the form on an open subset of R s × R r × R, with a parameter ε in a neighborhood of 0, a slow-fast system.
A classical approach to a rigorous foundation of quasi-steady state phenomena in chemical reaction networks goes back to Heineken et al. [28]: In order to obtain a slow-fast system from (4) some variables of the system that satisfy a compatibility condition are scaled by a small parameter.We outline a simplified version of this technique: Given a smooth curve ε → π * + ερ + . . . in the parameter space (with π * not necessarily a TFPV), we obtain a system (7) h(x, π * + ερ with "small" parameter ε.Note that As for the compatibility condition, we follow [34]. is called an LTC index set for (7), and the set of corresponding variables {x i 1 , . . ., x ir } an LTC variable set, if (8) h (0) (x) = 0, whenever (The acronym stands for "locally Tikhonov consistent" [34].)If the ODE system models a reaction network, then the corresponding species set {X i 1 , . . ., X ir } is called a set of LTC species.If the concentrations of all the species in an LCT set are all zero, then no reaction can take place.Note that (8) cannot be fulfilled if there are inflow reactions in the reaction network, as h (0) (x) contains a non-zero constant monomial.
For an LCT index set {i 1 , . . ., i r }, define and collect the remaining variables in u 2 .Partitioning and rewriting h * (x, ε) =: g(u 1 , u 2 , ε), one obtains a system ), with g 1 , g 2 being polynomials, provided that g 1 , g 2 are so, arriving at the slow-fast system as in ( 6), (10) u * In the singular perturbation reduction following Heineken et al. [28], one applies Tikhonov's theorem to (10), upon verifying the necessary conditions.In the literature, a frequently used shortcut is to directly solve g 1 (u 1 , u 2 , ε) = 0 (with small ε) for u 1 and substitute the result into the second equation of (9).We will refer to this procedure as classical QSS reduction.Note that without further analysis, e.g.verifying the hypotheses for Tikhonov's theorem, this is a purely formal procedure.
Obviously, any superset of an LTC species set is also an LTC species set, but minimal LTC species sets are of primary interest: For non-minimal LTC species sets the partial derivative D 1 g 1 necessarily has zero columns, so a local resolution of the implicit equation g 1 = 0 cannot exist.
Since solutions of (10) are bounded on compact subsets of their maximal existence interval, one finds u 1 = O(ε) on these compact subintervals.But it is not guaranteed that system (9) admits a local (n − r)-dimensional invariant manifold close to u 1 = 0 for small positive ε, hence there remains the question whether a singular perturbation reduction exists.Thus, LTC variable sets provide candidates for Tikhonov-Fenichel reductions, but these need further investigation.Moreover, even in the singular perturbation setting, there may not be a connection to TFPVs.We will get back to this later.
In some cases, direct application of Tikhonov-Fenichel does not work, but singular perturbation reduction with a critical variety of higher dimension is possible.For instance, Schneider and Wilhelm [41] considered a scenario, where the fast part of ( 10) admits non-trivial first integrals.In such a setting, the partial derivative D 1 g 1 cannot have full rank, but if the rank is full on every level set of the first integrals, and the non-zero eigenvalues have negative real parts, then reduction works.(Conversely, the local existence of such first integrals is also necessary [25, Prop.2].) For reaction networks it is of interest to understand whether first integrals of system (9) (and thus of ( 7)) carry over, upon scaling, to the fast system at ε = 0, and to a possible reduction.We first note an obvious fact.
which is a true but uninteresting fact due to u2 = 0.However, we do have: Proposition 2.9.Let φ(u 1 , u 2 , ε) denote a smooth first integral of system (9).Then, the following hold: , and φ is a first integral of (10), and φ(u * 1 , u 2 , 0) is (constant or) a first integral of the fast system.
Since stoichiometric first integrals are of particular importance for reaction networks, we note a special case: Let φ = m 1 u 1 + m 2 u 2 be a first integral of ( 9) with row vectors m 1 , m 2 .If m 2 = 0, then φ is also a first integral of the fast system.(Note that this is a rather restrictive condition.)Example 2.10.Consider the reversible Michaelis-Menten reaction network [32], ( 11) where X 1 chemically is a substrate, X 2 an enzyme, X 3 an intermediate complex, and X 4 a product, formed by conversion of the substrate X 1 .Using (1), we obtain the ODE system The right-hand side vanishes for instance when Tikhonov's theorem is not directly applicable to this slow-fast system, since the rank condition is not satisfied.But a step-by-step approach yields a reduction to dimension one: With the first integral φ 1 = x 2 + x 3 one obtains a three-dimensional system for x 1 , x * 2 and x 4 , for which a singular perturbation reduction to dimension two exists.Then with Proposition 2.9 the first integral φ 2 = x 1 + x 3 + x 4 allows a further reduction to dimension one.
3 TFPVs for reaction networks

General considerations
While the notion of TFPV applies to all parameter dependent polynomial (and more general) vector fields, special properties of reaction networks impose restrictions.We give an elementary illustration of this fact.This relation defines a cone in the parameter space.On the other hand, by (2), every linear 2 × 2 system describing a first order reaction network with two species (hence, the reaction network has the complexes X 1 , X 2 and possibly 0), takes the form with non-negative κ ij (κ ij is zero if the corresponding reaction does not exist).The determinant condition on the Jacobian of the system simplifies to due to non-negativity.In addition, the existence of stationary points requires conditions on κ 13 and κ 23 .Evaluating the TFPV conditions, one sees that they all admit an interpretation in the reaction network framework: Certain reactions are being "switched off".Furthermore, the conditions for TFPVs to exist yield very simple irreducible components of W s , namely coordinate subspaces.
For a number of standard reaction networks in biochemistry (in particular those described in the first chapter of Keener and Sneyd [32]), all TFPVs were determined algorithmically in Goeke's dissertation [23] and in the subsequent papers [26,27].It turned out that all of these admit an interpretation as a degenerate scenario in reaction network terms, via "switched off" reactions or missing species (and in some cases a combination of these).Based on these observations, and employing the theory of reaction networks, we will investigate conditions on reaction networks that guarantee the existence of singular perturbation scenarios.
For the reaction networks discussed in [26,27], one finds that every irreducible component of W s (see (5)) is just a coordinate subspace.Indeed, it is not easy to find (realistic) systems where some component of W s is not a coordinate subspace.This may be the case when non-stoichiometric first integrals exist for only some κ: Let s * be the codimension of the reaction network (the number of independent stoichiometric first integrals).Assume the set Π of κ's that give rise to extra linear first integrals is a proper algebraic variety and hence has measure zero (as for the reaction network in (3) in Example 2.4).Any point in Π is a candidate for a TFPV in dimension s > s * , if furthermore the critical manifold intersects the non-negative orthant and is attracting.Going back to network (3), the set Π consists of TFPVs and is characterized by the condition κ 1 = κ 2 , as one easily verifies that there exists a linearly attracting critical manifold.
An artificial way to construct further examples where the set of TFPVs is not included in a coordinate subspace, is to consider any parametrized polynomial system for which the dimension of the set of stationary points is larger than s * for some choice of parameters, and furthermore, all negative monomials of the i-th polynomial are multiples of x i .The latter is enough to constructively interpret the system as arising from a mass-action reaction network [11], though the networks obtained in this way are typically not realistic.The following example is generated in this way.
Example 3.2.Consider the following mass-action reaction network Generically, the variety of stationary points consists of the point (0, 0) and has dimension s * = 0.However, when κ 1 κ 4 = κ 2 κ 3 , then the variety has dimension one and consists of the line Additionally, a direct computation shows that the critical manifold is attracting for (x 1 , x 2 ) ∈ R 2 ≥0 .Hence, κ 1 κ 4 = κ 2 κ 3 defines a set of TFPVs for dimension one.In this case, there are no linear first integrals.
We now turn to system (1), and first establish conditions to ensure that every TFPV lies in some proper coordinate subspace, thus every irreducible component of W s is contained in some coordinate subspace.If we require the critical manifold to intersect the positive orthant, the existence of TFPVs κ ∈ R m >0 is easily precluded for important classes of reaction networks.In preparation for Proposition 3.3 we introduce some objects and some notation.Let N ′ ∈ R s * ×m consist of s * linearly independent rows of the stoichiometric matrix N .Moreover let E ∈ R m×q be a matrix whose columns are the extreme rays of the polyhedral cone ker(N ) ∩ R m ≥0 , and for λ ∈ R m denote by diag(Eλ) the matrix with the entries of Eλ in the diagonal, and zeros off-diagonal.Consider the matrix ( 12) (with B the reactant matrix, see Subsection 2.1).Finally, let Λ be the set of λ ∈ R q ≥0 such that Eλ ∈ R m >0 .(The particular choice of N ′ will be irrelevant.) Let G be a mass-action reaction network of codimension s * .With the notation introduced above, assume G belongs to one of the following cases: (a) The set of positive stationary points (b) The reaction network is injective [18], hence the coefficient σ n−s * (x, κ) of τ s * of the characteristic polynomial of the Jacobian of system (1) is a polynomial in x and κ with only non-negative coefficients.
(c) For all λ ∈ Λ, at least one of the minors of N ′ diag(Eλ)B ⊤ is non-zero.
Then, there are no TFPVs κ ∈ R m >0 for which some irreducible component Z of the critical variety intersects the positive orthant.
Proof.(a) The parametrization gives that the dimension of V κ is s * for all κ ∈ R m >0 .(b) For κ to be a TFPV, we require σ n−s * (x 0 , κ) = 0 for some x 0 ∈ Z, which occurs only if κ or x 0 belong to a coordinate subspace.(c) The condition implies that the Jacobian of (1) has rank n − s * at any positive stationary point [37].
Remark 3.4.We make a few observations regarding the relevance of the criteria in Proposition 3.3.
• Condition (c) holds for surprisingly many networks and is computationally easy to verify.
When Λ = R q >0 , which occurs often, then condition (c) holds if there is one minor with all non-zero coefficients of the same sign.
• Many realistic reaction networks admit parametrizations in the sense of Proposition 3.3(a): Among these are reaction networks admitting toric steady states [38], complex-balancing equilibria [8,12,29] (see also Subsection 3.3), and there are many reaction networks for which parametrizations can be found using linear elimination of some variables in terms of the rest [19,39] (see [6] for a short account on how to find parametrizations).
• Injective reaction networks admit at most one equilibrium in each SCC [9,18], and several criteria, in addition to the one stated in Proposition 3.3(b), have been established.These criteria involve graphical conditions [1,2] and sign vectors [36].
To include TFPVs with critical manifold intersecting the positive orthant, it is appropriate (and necessary in the cases covered in Proposition 3.3) to deviate from the convention in Definition 2.1 and allow the rate parameters to be zero, thus change the parameter range of κ to R m ≥0 .Passing from generic κ ∈ R m >0 to a special κ ∈ R m ≥0 may be seen as considering a subnetwork of the original reaction network.To indicate this, we make the following definition.
We denote by G(κ) the subnetwork obtained from G by removing the reactions with indices in {1, . . ., m} \ supp(κ), that is, the i-th reaction is removed if κ i = 0, for i = 1, . . ., m. Isolated nodes are not removed from G(κ), and hence G and G(κ) have the same set of complexes and species.
>0 is a TFPV for dimension s and thus U ⊆ W s .Since U is Zariski dense in C, its Zariski closure is C and it follows that C ⊆ W s .Proposition 3.6(b) does not imply that all rate parameters in the coordinate subspace are TFPVs given that one is a TFPV, but only that this is the case in an open set relative to the coordinate subspace.The next example illustrates this.
If we now consider κ = (κ 1 , 0, 0, κ 4 , 0, 0), the codimension of G( κ) is also s = 1, and the positive part of the stationary variety consists of the attracting line x 2 = κ 1 κ 4 .Hence, κ is a TFPV for dimension one.The minimal coordinate subspace containing κ is not an irreducible component of W 1 , as it is a proper Zariski closed set of the coordinate subspace C.
In what follows we consider TFPVs for two classes of reaction networks, namely first order reaction networks and complex-balanced reaction networks.Due to special properties of the Laplacian matrix, TFPVs for complex-balanced reaction networks can be identified.Our results build on the understanding of the kernel of A(κ) in (2).Therefore, we first review key results about Laplacian matrices and especially their kernel, using a graphical approach.

Some properties of Laplacian matrices
In this subsection we recall and review some properties of Laplacian and compartmental matrices.For the following known facts refer e.g. to Jacquez and Simon [30,Subsection 4.1].
The Laplacian matrix of a directed graph (and thus the Laplacian A(κ) of a reaction network) is a compartmental matrix.We recall some notions.
• A quadratic matrix with real entries is called a compartmental matrix if all its off-diagonal entries are ≥ 0 and all its column sums are ≤ 0.
Lemma 3.8.Let L(σ, τ ) be a compartmental matrix as in (13).Then, all eigenvalues of L(σ, τ ) have non-positive real part, and any eigenvalue with real part zero is equal to zero.Moreover R n is the direct sum of the kernel and the image of L(σ, τ ).
Consider a mass-action reaction network G. Let G 1 , . . ., G r be the connected components of G and further order the set of complexes according to the connected component they belong to.If A i (κ) stands for the Laplacian matrix of G i , then A(κ) becomes a block diagonal matrix with r blocks, The form of the kernel of the Laplacian matrix A(κ) of a digraph with κ ∈ R m >0 is well known, in particular, in the context of reaction networks [13,Thm 16.4.2].It derives from the Matrix-Tree theorem [5,35,44].
• The dimension of the kernel of A(κ) agrees with the number of terminal strongly connected components and is independent of κ ∈ R m >0 .• If the digraph is strongly connected, then dim ker A(κ) = 1 and a generator of ker A(κ) is given by the sequence of signed principal minors (which are positive).
• If the digraph is not strongly connected, then any complex in the support of a vector in ker A(κ) belongs to a terminal strongly connected component.Furthermore, a basis of ker A(κ) can be chosen such that the support of each vector is exactly one terminal strongly connected component and the non-zero entries are positive.These entries arise as the signed principal minors of the restriction of the matrix to the nodes in the component.
• The vector e = (1, . . ., 1) belongs to the left-kernel of A(κ), and generates it when the digraph has one terminal strongly connected class.
These facts lead to the following lemma.
Lemma 3.9.Let G be a mass-action reaction network with labelling κ ∈ R m >0 .Let G 1 , . . ., G r be the connected components of G and assume that the set of complexes is ordered in accordance with the components.Let T be the number of terminal strongly connected components of G.
Then, the rank of A(κ) does not depend on the choice of κ ∈ R m >0 , and in particular (a) dim ker A(κ) = T .
(b) ker A(κ) has non-trivial intersection with the positive orthant R d >0 , if and only if G is weakly reversible.

Proof. (a-b) are direct consequences of the properties of the kernel of A(κ) discussed above. (c)
The column sums in each block A i (κ) are zero, as each submatrix is a Laplacian.The second part follows from (a), as T = r.
With the notation in Lemma 3.9, if a connected component of G has more than one terminal strongly connected component, then the vectors given in Lemma 3.9(c) do not form a basis of the left-kernel of A(κ).To obtain a basis, one has to augment them by vectors that might depend on the particular entries of A(κ), that is, on κ; see Example 2.4 for an illustration.

Complex-balancing and TFPVs
We now turn to an important class of reation networks called complex-balanced reaction networks and the existence of TFPVs for this class.Complex-balanced reaction networks are characterised by their equilibria, called complex-balanced equilibria.According to Horn and Jackson [29] (see also Feinberg [13,Ch.15ff.]), a positive equilibrium z ∈ R n >0 of ( 2) is complex-balanced for the parameter value κ * if A(κ * )z Y = 0.
By Lemma 3.9(b), the existence of a positive complex-balanced equilibrium implies that the reaction network is weakly reversible [ • Let e 1 , . . ., e d denote the standard basis of R d and let ∆ : The dimension δ of Ker Y ∩ span ∆ is called the deficiency, and satisfies δ = d − (n − s * ) − r.
• If system (2) admits a complex-balanced equilibrium in R n >0 , then every SCC contains precisely one positive equilibrium, which also is complex-balanced, and the Jacobian has n−s * eigenvalues with negative real part and s * zero eigenvalues (counted with multiplicity).As a consequence, the positive equilibria of the system form a manifold of dimension s * .
If G is weakly reversible and δ = 0, then all positive equilibria are complex-balanced, irrespective of the (positive) reaction rate constants.In general, there are δ algebraically independent relations on the rate parameters κ, characterizing when the reaction network admits positive complexbalanced equilibria.These relations are explicit [8,10,15].Complex-balanced equilibria form a manifold of dimension s * , and the rank of the Jacobian of system (2) evaluated at the equilibrium is n − s * .
In the following theorem we use the notation G( κ) introduced in Definition 3.5.
Theorem 3.11.Let G be a mass-action reaction network of codimension s * .Let κ ∈ R m ≥0 such that • G( κ) is weakly reversible of codimension s > s * .
• The non-zero coordinates of κ satisfy the relations for the existence of positive complex-balanced equilibria in G( κ).
Then, κ is a TFPV for dimension s of system (2).Furthermore, the minimal coordinate subspace containing κ is contained in W s .
Proof.We verify properties (vi)-(iv) of TFPVs.By Proposition 3.10, the dimension of the set of positive equilibria of G( κ) is s.The remaining properties of a TFPV follow from the properties of complex-balanced equilibria in Proposition 3.10.The last statement follows from Proposition 3.6(b).
An immediate consequence of Theorem 3.11 arises when G is weakly reversible and has deficiency zero.A key point is that the deficiency of any subnetwork obtained from G by removing reactions can only decrease [31,Prop. 8.2].In particular, if G has deficiency zero, then so does any subnetwork.Theorem 3.12.Let a weakly reversible reaction network G of deficiency zero be given, with dynamics governed by system (2).Let κ ∈ R m ≥0 be such that the induced subnetwork G( κ) is weakly reversible and has more connected components than G.
Then κ is a TFPV of system (2) for dimension n − d + r, with d and r the number of complexes, respectively, connected components of G( κ).This dimension equals the codimension of G( κ).
Proof.Let r * be the number of connected components of G.As the deficiencies of G and G( κ) are zero, the codimensions of G and G( κ) are s * = n − d + r * and s = n − d + r, respectively.As r > r * , we have s > s * .Furthermore, all parameter values κ yield complex-balanced equilibria for G( κ) as the deficiency is zero.The statement now follows from Theorem 3.11.
In particular, when the hypotheses of Theorem 3.12 hold, then κ lies in a coordinate subspace of the parameter space.Moreover, the connected components of G( κ) identify a coordinate subspace in W s for the appropriate s > s * .For some classes of reaction networks, including weakly reversible reaction networks of deficiency zero, an explicit formula for the singular perturbation reduction was derived in [16].
Example 3.14.The competitive inhibition reaction network with reversible product formation [32], has two additional reactions compared to the reversible Michaelis-Menten reaction network, see Example 2.10, namely, inhibition of the enzyme (X 2 ) by an inhibitor (X 5 ) via formation of an intermediate complex (X 6 ).The reaction network is weakly reversible with deficiency zero (five complexes, two linkage classes and stoichiometric subspace of dimension three).By Theorem 3.12, setting either κ 1 = κ −1 = 0, or κ 2 = κ −2 = 0, or κ 3 = κ −3 = 0, the number of connected components increases by one, and the resulting rate parameters are TFPVs for dimension 4. In addition, choosing two of the three pairs to be zero, one obtains TFPVs for dimension 5.
Example 3.15.Consider the following reaction network, which is the futile cycle with one phosphorylation site [46]: This reaction network is not weakly reversible and has codimension 3.An easy computation shows that the stationary set admits a parametrization with three free variables x 1 , x 2 , x 3 , and hence has dimension 3. Proposition 3.3(a) applies.Alternatively, Proposition 3.3(c) is applicable: For a specific choice of N ′ , the matrix We have Λ = R 3 >0 .The minor given by columns 1, 2, 5 is λ 1 λ 2 2 , which is non-zero.In conclusion, no TFPVs with positive entries and critical manifold intersecting the positive orthant exist.Upon setting κ 3 = κ 6 = 0, the resulting reaction network is weakly reversible and has deficiency 0 with codimension s = 4. Hence, by Theorem 3.11, any rate parameter of the form (κ 1 , κ 2 , 0, κ 4 , κ 5 , 0) with non-zero entries being positive, is a TFPV for dimension 4.
In this case, the positive part of the stationary variety of G( κ) always admits a parametrization.Using Hurwitz determinants, one confirms that the variety is linearly attracting.Therefore, the whole positive part of this particular coordinate subspace is formed by TFPVs.We conclude with an example of a weakly reversible reaction network admitting TFPVs that are not included in a proper coordinate subspace.
Example 3.17.This example is introduced in [4, Example 4.1], where the purpose is to show the existence of weakly reversible reaction networks with infinitely many equilibria in some SCC.The reaction network consists of four connected components, written in rows for convenience: This reaction network is weakly reversible of codimension s * = 0.When all parameters are set to 1 except for κ 3 = κ 8 = κ 10 = κ 13 = a, with a > 5, then the stationary variety has dimension 1: it consists of one unstable point (1, 1) and one attracting closed curve around (1,1).Hence, any such rate parameter is a TFPV for dimension 1, which does not belong to a proper coordinate subspace.One might note that the reaction network is of the form discussed in Example 3.2 with all negative terms in the ODE system being multiples of x 2 .

TFPVs for first order reaction networks
In this section, we consider the special case of a mass-action reaction network G = (Y, R, κ) containing only first order reactions; thus d = n or d = n + 1 and non-zero complexes may be identified with species.In the formulation the matrix Y is simply the identity matrix if 0 / ∈ Y and the identity matrix with an extra zero column otherwise.Hence either x Y = x or x Y = (x, 1) ⊤ .Remark 3.19.By Lemma 3.9, the rank of A(κ) does not depend on κ ∈ R m >0 .Let T be the number of terminal strongly connected components of G and s * be the codimension of G.We make the following observations: ∈ Y, then Y is the identity matrix, and the solution set to ( 14) in R n ≥0 is ker A(κ) and s * = T .
-If 0 belongs to a terminal strongly connected component of G, then the solution set to (14) in R n ≥0 is the linear affine subspace of ker A(κ) ∩ R n ≥0 with last coordinate equal to 1.By the description of ker A(κ) in Lemma 3.9, this subspace has dimension T − 1.
-If 0 does not belong to a terminal strongly connected component of G, then (14) has no solution.Indeed, the last entry of x Y is equal to 1, and hence positive, but any vector in ker A(κ) has last entry zero.
With this in mind, we obtain the following proposition.
Proposition 3.20.Let A(κ) be the Laplacian matrix of a mass-action reaction network G = (Y, R, κ) consisting only of first order reactions with dynamics governed by system (14) in R n ≥0 .Let T be the number of terminal strongly connected components of G, and s * the dimension of the solution set to (14).
≥0 be in a proper coordinate subspace of R m , and consider the subnetwork G( κ).Then κ is a TFPV if and only if G( κ) has more than T terminal strongly connected components, and additionally the complex 0 belongs to one such component, provided 0 is a complex of G.
Proof.(a) The proof is straightforward as the dimension of the solution set to (14) does not depend on κ, provided all entries are positive.(b) We first make a digression.Consider G( κ) and assume 0 ∈ Y. Then the last column of the matrix Y is zero and the last entry of v(x) = x Y is 1.Let A( κ) be the submatrix of A( κ) obtained by removing the last row and column.Let β ∈ R d−1 be the vector formed by the first d − 1 entries of the last column of A(κ).Then Y A( κ)v(x) = A( κ)x + β.To prove (b), we apply Lemma 3.8 to the compartmental matrices A( κ) or A( κ), depending on whether 0 is a complex of G. (c) is a direct consequence of (a) and (b), as Π s is a union of coordinate subspaces of R m ≥0 .
Rephrasing the statement of Proposition 3.20, all TFPVs are found by setting rate parameters to zero such that the number of terminal strongly connected components increases, and taking into consideration the role of the zero complex.We note that the irreducible components of any W s can be identified by inspecting the graph G.
If the considered first order reaction network G in addition is weakly reversible, then for this network Theorem 3.11 and Theorem 3.12 are both consequences of Proposition 3.20.For Theorem 3.11, note that if the subnetwork G( κ) of G is weakly reversible with codimension s > s * , then it must be that κ ∈ R m ≥0 belongs to a proper coordinate subspace of R m and the number of terminal strongly connected components of G( κ) exceeds the number of terminal strongly connected components of G. Hence, the conclusions of Theorem 3.11 follow from Proposition 3.20(b),(c).Note that the second condition of Theorem 3.11 is trivially fulfilled because G has deficiency zero, hence any subnetwork, in particular G( κ), has also deficiency zero [31,Prop. 8.2].For Theorem 3.12, we remark that it is a consequence of Theorem 3.11, hence also of Proposition 3.20.Alternatively, it follows directly from Proposition 3.20 by similar arguments to above.
Example 3.21.Consider a first order reaction network with three complexes and four reactions, This reaction network has one terminal strongly connected component.By Remark 3.19, s * = 1.There are three coordinate subspaces yielding TFPVs for dimension 2. These arise from the three ways to increase the number of terminal strongly connected components: Example 3.22.For the first order reaction network we have two connected components and s * = 1 (Remark 3.19), but this reaction network has no stationary points.Upon setting κ 2 = 0, we have three connected components and 0 belongs to a terminal strongly connected component.Hence, by Proposition 3.20, (κ 1 , κ −1 , 0) is a TFPV for dimension 2.

Scalings, stoichiometry and TFPVs
In this section, we start from LTC variable sets and the scaling approach to singular perturbation reductions of system (2), as initiated by Heineken et al. [28] (recall Subsection 2.3 on slow-fast systems).A priori, there are no TFPVs that correspond to scalings, but these may appear when the system is restricted to SCCs, as new parameters are introduced.For motivation, we look again at the reversible Michaelis-Menten system.
Example 4.1.We continue Example 2.10.The LTC variable set {x 2 , x 3 } corresponds to the stoichiometric first integral φ 1 = x 2 + x 3 , and the LTC variable set {x 1 , x 3 , x 4 } corresponds to the stoichiometric first integral φ 2 = x 1 + x 3 + x 4 .Moreover, on the SCC given by x 2 + x 3 = e 0 and x 1 + x 3 + x 4 = s 0 , one obtains the 2-dimensional system This system admits a TFPV with e 0 = 0, and all other parameters > 0, with a degenerate (one dimensional) SCC forming the critical manifold, and a subsequent singular perturbation reduction.
(For a TFPV with s 0 = 0, the SCC degenerates into a single point.) Quite generally, LTC variable sets point to bifurcation scenarios, and possibly interesting dynamics may appear for small perturbations.In general there is no perfect correspondence to stoichiometric first integrals, as shown by examples in [34].But stoichiometric first integrals which correspond to LTC variable sets may, in turn, yield TFPVs of the system on SCCs.
We start by characterizing LTC species sets.

A characterization of LTC species sets
A useful modification of system (2) is the following, when some complexes are non-reactant complexes, that is, they only appear as product complexes.Complex Y j is non-reactant if and only if column j of A(κ) is zero.Thus, one may form Y * from Y , respectively, A * (κ) from A(κ), by removing all columns that correspond to indices of non-reactant complexes, to rewrite ( Let d * be the number of reactant complexes, hence Y * ∈ N n×d * 0 .Note that A * (κ) is not a square matrix unless all complexes are reactant complexes and thus A * (κ) = A(κ).
We will first and foremost discuss sets that are LTC species sets for all parameter values κ ∈ R m >0 .The equations x Y = 0, respectively, x Y * = 0 define varieties with coordinate subspaces as irreducible components, and the corresponding variables are obviously LTC variables.We will first show that all LTC variable sets of system (15) (which is the same as system (1)) are of this type.
The following proposition characterizes the LTC species sets.Recall that reaction networks with inflow reaction do not admit any LTC species sets (remark below Definition 2.7).Proposition 4.2.Let system (1) be given.Then {i 1 , . . ., i u } with u < n and Proof.The non-trivial assertion is the "only if" part.The "if" part follows by definition.We need then also x Y * = 0. We may assume that the LTC index set is {1, . . ., u}, and that complexes are ordered such that the first d * are reactant complexes.We let y 1 , . . ., y d * denote the columns of Y * .
We argue by contradiction and assume that some = 0 and i = 1, . . ., d * , we have (16) x We may assume that ( 16) holds precisely for the indices d ′ ≤ i ≤ d * , for some d ′ ≤ d * .Thus, we aim to show d ′ = d * .
Let K κ be the (m × d * )-matrix with non-negative entries such that Each entry of K κ is one of the rate parameters: the (i, j)-th entry is κ ℓj if the i-th reaction is Y j → Y ℓ .As y i , i = 1, . . ., d * , are pairwise different, the monomials x y i , d ′ ≤ i ≤ d * , are linearly independent over R. Using N K κ x Y * = 0, we obtain with the last d * − d ′ + 1 ≥ 1 columns equal to zero.The equality tells us that the last d * − d ′ + 1 columns of K κ belong to ker(N ) for all κ ∈ R m >0 .The sum of these columns lies also in ker(N ).The entries of the sum are positive when they correspond to reactions with Y d ′ , . . ., Y d * among the reactant species, and zero otherwise.As κ varies, we thus obtain a relatively open and non-empty subset in some proper coordinate subspace C of R m .The row space of N is therefore orthogonal to C; hence N has at least one zero column, and we have reached a contradiction, as a reaction network does not have self-edges.
Corollary 4.3.LTC species sets are identifiable from the reactant complexes: A set of species {X i 1 , . . ., X iu } is an LTC species set if and only if in every reactant complex, one of the X i k appears with positive coefficient.
An enumeration of all LTC species sets may start from those reactant complexes that contain the fewest species (that is, the species appearing with positive stoichiometric coefficients).First, a species that appears alone in some reactant complex is necessarily contained in every LTC species set.Then proceed with complexes containing two species, and so on.From this observation, one also finds that LTC species sets for first order reaction networks (with every complex consisting of one species) are comprised of all species in reactant complexes.Hence, the notions of LTC species and LTC variables are of real interest only for non-linear systems.
Example 4.4.In Example 2.10, the reversible Michaelis-Menten reaction network, the reactant complexes are X 1 + X 2 , X 3 and X 4 + X 2 .Thus X 3 must lie in every LTC species set, and so must X 2 or X 1 .The first alternative yields the LTC species set {X 2 , X 3 }, while the second yields the LTC species set {X 1 , X 3 , X 4 }.These are the only two LTC species sets.In contrast, the standard irreversible Michaelis-Menten reaction network without the reaction X 4 + X 2 κ −2 −−→ X 3 has reactant complexes X 1 + X 2 and X 3 , with two LTC species sets, {X 1 , X 2 } and {X 2 , X 3 }.
Example 4.5.We consider again the futile cycle with one phosphorylation site, see Example 3.15: Here, X 5 and X 6 are contained in every LTC species set, and altogether one finds the following LCT species sets, Only the first and the last of these are also LTC species sets for the fully reversible system with the additional reactions

LTC species and first integrals
We proceed to study the relation between LTC species sets and linear first integrals.We first note a relation between LTC indices and the complex matrix.
Lemma 4.6.Let {i 1 , . . ., i u } with u < n and Then the following statements are equivalent.
(b) The support of every column of Y * contains some i k .
(c) There exists a non-negative row vector ω ∈ N n 0 with support {i 1 , . . ., i u } and such that every entry of ω • Y * is positive.
Proof.The equivalence of (a) and (b) is a restatement of Corollary 4.3.As for the equivalence of (b) and (c), note that Thus, the j-th entry of ω • Y * is positive if and only if y i ℓ ,j > 0 for some i ℓ .As the (i, j)-entry of Y * is the stoichiometric coefficient of X i in the complex Y j , we have that (ω • Y * ) j > 0 if and only if the support of the j-th column of Y * intersects {i 1 , . . ., i u }.The assertion follows.
The next example shows that non-negativity of the coefficients of the stoichiometric first integral in Proposition 4.7 cannot be discarded in general.
with κ 3 > 0. As in the system without degradation (Example 4.4), {X 2 , X 3 } is an LTC species set, but the only stoichiometric first integral (up to multiples) is φ = x 1 −x 2 +x 4 , due to (1, −1, 0, 1) Y = 0.The set {1, 2, 4} is not an LTC index set.As noted in [34], the example also shows that the scaling approach may yield singular perturbation scenarios which are not directly related to TFPVs (even after restricting to SCCs).One verifies that scaling x 2 and x 3 yields a system that admits a singular perturbation reduction to dimension two, with trivial reduced equation.
Remark 4.10.There remains the question under which conditions the existence of LTC species sets in turn implies the existence of stoichiometric first integrals with corresponding support.We give a characterization for reaction networks with one connected component and one terminal strongly connected component.Thus, let system (2) represent such a reaction network.Assume without loss of generality that {X 1 , . . ., X u } is an LTC species set, and denote by ȳ1 , . . ., ȳu the first rows of the complex matrix Y .By Lemma 3.9 the system admits a stoichiometric first integral if, and only if, e = (1, . . ., 1) is a multiple of some element in the closed convex hull of ȳ1 , . . ., ȳu .(Note that due to Lemma 4.6(c), there exist integers ω 1 > 0, . . ., ω u > 0 such that u i=1 ω i ȳi > 0 (coordinate-wise).)

Stoichiometry and TFPVs
We now address TFPVs of system (2) versus TFPVs of its restriction to stoichiometric compatibility classes.As seen in Example 4.1, the restricted system may admit additional TFPVs.We first fix some notation.We introduce the abbreviation In the following, we will assume that system (17) admits a maximal set of independent stoichiometric first integrals φ 1 , . . ., φ s * .Then every SCC is the intersection of R n ≥0 with the common level set which we abbreviate as S θ , θ = (θ 1 , . . ., θ s * ).One may choose x ∈ R n−s * with entries from x 1 , . . ., x n , such that the Jacobian of ( x, φ 1 (x), . . ., φ s * (x)) has full rank n.This yields an equivalent version which for given θ represents system (17) on S θ .We are interested in TFPVs of the (n − s * )dimensional system (18) for dimension s > 0. Possible candidates for TFPVs are as follows.
• TFPVs via "inheritance" from ( 17): If κ is a TFPV of (17) for dimension s > s * , then ( κ, θ) is a TFPV of ( 18) for some dimension > 0 and some θ whenever a transversality condition is satisfied.This condition is rather weak: It does not hold only if the sum of the stoichiometric subspace and the tangent space of V(h(•, κ)) at every x 0 ∈ V(h(•, κ)) has dimension < n−s * +s.
• TFPV candidates from stoichiometric first integrals: Let the setting of Proposition 4.7 be given and assume that the stoichiometric first integral φ ℓ with non-negative coefficients corresponds to an LTC variable set.If there exists a TFPV ( κ, θ) with θ ℓ = 0, then the critical variety will be a coordinate subspace, and consequently by [27] the singular perturbation reduction will agree with the "classical" QSS reduction (in the sense of Subsection 2.3) for the LTC variables.
(We restrict attention to a single first integral here, since we are interested in minimal LTC sets; cf.Section 2.3.) There remains to establish manageable criteria for TFPVs from stoichiometric first integrals.The next result yields conditions for parameter values that are "almost TFPV".Proposition 4.11.Let system (1) be given, and assume that every SCC of this system is compact (equivalently, the left-kernel of N in (1) has a vector with all entries positive [3]).Moreover assume that there exists a parameter θ ∈ R s * such that: (a) No stationary points in S θ are isolated relatively to S θ .
(b) For every ρ > 0, there exists some θ such that θ − θ < ρ and ˙ x = h( x, κ, θ) admits an isolated linearly attracting stationary point.(Here • denotes some norm.)Then, ˙ x = h( x, κ, θ) admits a non-isolated stationary point whose Jacobian has only eigenvalues with non-positive real part, and admits zero as an eigenvalue.
Proof.Given a compact subset K of the parameter space, the union of the SCCs S θ with θ ∈ K is compact.In the following, let K be a compact neighborhood of S θ .
For every positive integer L let θ L ∈ K be such that θ L − θ < 1/L and ˙ x = h( x, κ, θ L ) admits an isolated linearly attracting stationary point z L .By compactness, the sequence ( z L ) L in R n−s * has an accumulation point z, in S θ .Since z is not isolated, the Jacobian of D 1 h( z, κ, θ) has the eigenvalue zero.Moreover, the map which sends ( x, κ, θ) to the coefficients of the characteristic polynomial (19) χ of D 1 h( x, κ, θ) is continuous.Thus, if some eigenvalue of the Jacobian had positive real part, the same would hold for some eigenvalue of the Jacobian of D 1 h( z L , κ, θ L ) with L sufficiently large (see e.g. the reasoning in Gantmacher [22], Ch.V, section 3).
Corollary 4.12.Assume that system (18) describes the dynamics of a weakly reversible deficiency zero reaction network.Then, the conclusion of Proposition 4.11 holds for every θ, such that no stationary points are isolated in S θ .
Remark 4.13.We clarify here what is meant by "almost TFPV" prior to the statement of Proposition 4.11.For this we discuss the conditions for TFPV in Definition 2.5 in Subsection 2.2 for system (18) and parameter value θ.
• Condition (i) is always satisfied for some dimension > 0, due to condition (a) in Proposition 4.11.
• Condition (ii) requires equality of geometric and algebraic multiplicity for the eigenvalue 0. This holds automatically when the algebraic multiplicity is equal to one.Generally this property can be checked by algebraic methods: For x in the critical manifold, τ divides the characteristic polynomial in (19).Obtain a new polynomial η from χ ( x,κ,θ) by dividing out a power of τ such that a single factor τ remains.Then the multiplicity equals one if and only if η annihilates D 1 h( x, κ, θ).But (as mentioned e.g. in [25,Example 4]) there exist realistic reaction networks for which the direct sum condition on the kernel and the image does not hold.
• Finally, to guarantee condition (iii), one needs to verify that there exist no purely imaginary eigenvalues except 0. However, if (iii) is not satisfied, then the system may admit some interesting dynamics, like zero-Hopf bifurcations.
We note a sharper result for the case of a one dimensional critical variety.The non-trivial restrictions on the σ i in Corollary 4.14(b) suggest that there will be non-zero purely imaginary eigenvalues only in exceptional cases.There is an obvious (but less readily applicable) extension of Corollary 4.14 to TFPVs for dimension strictly larger than one, with an additional requirement that the geometric and the algebraic multiplicities of the zero eigenvalue are equal in Corollary 4.14(a), and that Corollary 4.14(b) is left unchanged except for the number of variables of Φ.
The polynomial Φ can be determined explicity.We recall some cases for SCCs of small generic dimension from [33,  (c) If the SCCs of system (2) generically have dimension four, and the hypotheses of Proposition 4.11 are satisfied, then θ is a TFPV for dimension one whenever σ 3 = 0 and σ 1 σ 2 = σ 3 at ( z, κ, θ).
The dynamics on an SCC is described by the ODE system, The Jacobian on the SCC with e 0 = 0 (thus, on the critical manifold with x 3 = x 6 = 0) is equal to: and the coefficients of its characteristic polynomial are Since both σ 1 and σ 2 are positive when 0 ≤ x 1 ≤ s 0 , the conditions in Remark 4.15 and Corollary 4.12 are satisfied, and θ is a TFPV for dimension one.This system was discussed by elementary means in [24], with no reference to reaction network theory.A comparison shows that the approach developed here saves substantial computational effort.Moreover, one verifies that Proposition 4.11 is also applicable to the system with irreversible product formation (that is, κ −2 = 0).

Partial scalings: An outlook
So far we considered on the one hand TFPVs that, in reaction network interpretation, arise from "switching off" certain reactions (Theorems 3.11 and 3.12).On the other hand, we introduced LTC species sets, with characterizing property that if their concentrations are zero, then all reactions of the reaction network are precluded from taking place, and discussed their relation to TFPVs (Proposition 4.11 and Corollary 4.14).It is suggestive to combine these approaches by switching off certain reactions and determining LTC species for the remaining reactant complexes, and this combination will be sketched next.In the setting of Subsection 2.3 and in particular expansion (7), we consider LTC variable sets for a specific choice of the parameter π * in (7), for reaction networks π * = κ.This yields a slow-fast system which may further admit a Tikhonov-Fenichel reduction.We will not attempt to establish necessary or sufficient conditions for a singular perturbation setting.We start again from (15), but now we consider some κ such that A( κ) has zero columns, thus there are additional non-reactant complexes in the reaction network G( κ).We may assume that the remaining reactant complexes correspond to columns y 1 , . . ., y d of Y * , thus The matrices of this type define some coordinate subspace of parameter space.We denote by Y 1 the matrix with columns y 1 , . . ., y d , and by Y 2 , the matrix with the remaining columns of Y * .LTC variable sets for G( κ) can be identified via Proposition 4.2 with the complex matrix Y 1 .
Upon relabelling, we may assume that x 1 , . . ., x u form an LTC variable set for G( κ).Considering a curve ε → κ + ερ + . . . in parameter space, we obtain with A 11 ∈ R u× d .Moreover, set x i = εx * i for 1 ≤ i ≤ u, then we have One would arrive at the same slow-fast system by starting from a different vantage point: First designate LTC variables and then switch off all reactions whose source complexes do not contain these variables.
Since the fast part of the scaled system involves slow reactions corresponding to A * 12 , the results from the previous subsections do not carry over to partial scalings.We will not discuss these matters any further here.
To close the present paper, we are satisfied to indicate by example that this heuristic is worth pursuing for general reaction networks.
additionally n − s * eigenvalues with negative real part), then the same holds for a norm-open neighborhood of κ * , and thus for a Zariski dense subset of R m >0 containing κ * .In particular, an irreducible component of V(h(•, κ * )) has dimension s * and intersects the positive orthant.(b) If κ * ∈ R m ≥0 is a TFPV for dimension s > s * with s the codimension of G(κ * ), then the minimal coordinate subspace containing κ * is contained in W s .Proof.(a) We will use that D 1 h(x * , κ * ) has n − s * eigenvalues with negative real part, if and only if the corresponding n − s * Hurwitz determinants of its characteristic polynomial, divided by τ s * , are positive, see Gantmacher [22, Ch.V, section 6].The rank of D 1 h(x * , κ * ) being n − s * implies that an irreducible component of V(h(•, κ * )) has dimension s * and intersects the positive orthant [7, §9.6 Thm 9].Let V be the real algebraic variety in the variables x, κ consisting of points where h(x, κ) = 0 and D 1 h(x, κ) has rank strictly smaller than n − s * , respectively, at least one of the Hurwitz determinants vanishes.By hypothesis, there exists (x * , κ * be an open Euclidean ball containing (x * , κ * ).The intersection of U and the zero set of h, which is non-empty, consists of points (x, κ) such that x ∈ V(h(•, κ)) and the Jacobian has maximal rank n − s * , respectively, all Hurwitz determinants are positive.The projection U of U onto R m >0 in the variable κ contains a non-empty open Euclidian ball of parameters κ 0 for which there exists x 0 ∈ R n >0 such that D 1 h(x 0 , κ 0 ) has rank n − s * , respectively, additional n − s * eigenvalues with negative real part.By the Implicit Function Theorem applied to h at (x * , κ * ), U contains an open ball centred at κ * such that R n >0 ∩ V(h(•, κ)) = ∅ for all κ in the ball.As any Euclidean ball is Zariski dense, this concludes the proof of (a).(b)Let C be the minimal coordinate subspace containing κ * .The parameter κ * and the reaction network G(κ * ) satisfy the hypotheses of (a), after restricting R m >0 to C. Therefore, there exists a norm-open and Zariski dense set U (relative to C) such that any κ

Lemma 3 .
18.A first order mass-action reaction network has deficiency zero.Proof.With the notation introduced in Proposition 3.10, one has Ker Y ∩ span ∆ = {0}, due to the form of Y .The assertion follows.

Example 4 . 9 . 3 −
Consider the reversible Michaelis-Menten mass-action reaction network in Example 2.10, with degradation of the intermediate complex (the reaction X 3 κ −→ 0), governed by the ODE system,

Corollary 4 . 14 .
In the setting of Proposition 4.11, consider system(18), with characteristic polynomial of the Jacobian given by(19).Let θ be such that σ n−s * = 0 and σ n−s * −1 = 0 for θ = θ, some κ ∈ R m >0 and some stationary point z ∈ R n−s * ≥0 .Then it holds: (a) The eigenvalue 0 of D 1 h( z, κ, θ) has multiplicity one.(b) There exists a polynomial Φ in n − s * − 1 variables with the following property: The Jacobian D 1 h( z, κ, θ) admits non-zero purely imaginary eigenvalues if and only if Φ( σ 1 , . . ., σ n−s * −1 ) = 0 at ( z, κ, θ).Proof.(a) is obvious.(b) There exists a polynomial Φ in the coefficients of the characteristic polynomial that vanishes if and only if a pair of (non-zero) eigenvalues adds up to zero; see e.g.[33, Lemma 4.1, Appendix B].Since all eigenvalues have real part ≤ 0, such a pair of eigenvalues must have zero real parts.