A generalization of Birchs theorem and vertex-balanced steady states for generalized mass-action systems

1.   Introduction

Reaction networks are commonly used to model natural phenomena in disciplines ranging from chemistry, biochemistry, epidemiology to population dynamics. In these systems, entities interact to form other entities as prescribed by a directed graph, the reaction network. For example, the reaction network

describes an enzymatic system, where a substrate $\rm{S}_0$ is converted into a product $\rm{S}_1$ by an enzyme ${\rm{E}}$ via an intermediate species $\rm{ES}_0$.

The concentrations of the chemical species in a network are often modelled by a system of ordinary differential equations. One of the most common assumptions in chemistry and biochemistry is that of mass-action kinetics, where the reaction rate is proportional to the concentrations of its reactants. According to mass-action kinetics, the reaction $\rm{E + S}_0 {\rightarrow} \rm{ES}_0$ proceeds at rate $\kappa_1 [{\rm{E}}][\rm{S}_0]$, where $\kappa_1 > 0$ is a rate constant, and $[\rm{X}]$ is the concentration of species $\rm{X}$ as a function of time $t$. The rates of change of the concentrations of ${\rm{E}}$, ${\rm{S_0}}$ and ${\rm{ES_0}}$ due to this single reaction are

The rates of change due to all reactions in the network is the sum over its individual reactions, e.g.

Mass-action systems have been studied extensively. Reaction network theory, as initially developed by Horn, Jackson and Feinberg ^[1,2,3], tries to conclude dynamical properties from simple characteristics of the underlying network. Moreover, as the reaction rate constant is usually obtained empirically and thus subjected to uncertainty, an ideal theoretical result does not depend on the precise values of the rate constants; indeed this is the case for many classical results in reaction network theory.

Mass-action kinetics assumes that the system is dilute (having low concentrations) and homogeneous (well-mixed). In the context of systems biology, that is not the typical environment; the cell is typically crowded and highly structured. Various models have been developed to account for this difference.

Biochemical systems theory ^[4,5] proposes power-law kinetics, where the exponents (or kinetic orders) in the reaction rate functions need not follow the stoichiometric coefficients. In the catalysis example above, we may want the concentration of ${\rm{E}}$ to be modelled by the equation

for some constants $\alpha$, $\beta$, $\gamma$, $\delta > 0$. This is an example of power-law kinetics. Classical mass-action kinetics and power-law kinetics can be incorporated into the framework of generalized mass-action kinetics as formulated in ^[6] (based on ^[7]).

Generalized mass-action systems can also be used to study mass-action systems that do not obviously admit nice dynamical properties. This is done by network translation, where a mass-action system is rewritten as a generalized mass-action system, where the underlying network has better properties (e.g., weakly reversible) ^[8].

Generalized mass-action systems are essentially dynamical systems of the form

where $\kappa_i \in {\mathbb{R}}_{ > 0}$, and ${\boldsymbol{u}_i}$, ${\boldsymbol{v}_i} \in {\mathbb{R}}^n$. (For ${\boldsymbol{x}} \in {\mathbb{R}}_{ > 0}^n$ and ${\boldsymbol{u}} \in {\mathbb{R}}^n$, we are using the notation ${\boldsymbol{x}}^{{\boldsymbol{u}}} = x_1^{u_1} x_2^{u_2} \cdots x_n^{u_n}$). For example, any polynomial dynamical system is of the form (1.1). Moreover, many classes of nonlinear ODEs can be recast as generalized mass-action systems ^[9,10]. For a complete definition of generalized mass-action systems, see Section 2. In this work, we are interested in the existence and uniqueness of steady states of these systems, as it relates to geometric properties of the vectors $\{{\boldsymbol{u}_i}, {\boldsymbol{v}_i}\}_{i \in I}$.

In classical mass-action systems, some classes of positive steady states enjoy certain algebraic and dynamical properties. Dating back to Boltzmann's kinetic theory, detailed-balanced equilibria can be regarded as thermodynamic equilibria. Their generalization, complex-balanced equilibria, are dynamically stable because of the existence of an associated Lyapunov function ^[2,11,12], and admit monomial parametrizations ^[13].

For generalized mass-action systems, the analogue of complex-balanced equilibria are the vertex-balanced steady states. Unsurprisingly, the theory of vertex-balanced steady states is quite complicated. Some necessary conditions for stability have been found recently ^[14]. Also, they admit monomial parametrizations that may be very useful in applications ^[6].

In this paper, we are interested in how many (if any) vertex-balanced steady states there are within each invariant affine subspace of a generalized mass-action system. In particular, we aim to understand which reaction networks admit vertex-balanced steady states, and whether they are unique. Interestingly, this question can be reformulated in two different ways, one related to a generalization of Birch's theorem in statistics ^[15], and the other to the bijectivity of generalized polynomial maps, similar to ones which appear in geometric modelling ^[7,16]. Indeed, the following questions are essentially equivalent:

1. When does a generalized mass-action system have exactly one vertex-balanced steady state within each invariant affine subspace, for any choice of rate constants?

2. Given vector subspaces $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$, when does the intersection^* $({\boldsymbol{x}_0} + S)\cap ({\boldsymbol{x}}^* \circ \exp\widetilde{S}^\perp)$ consist of exactly one point, for any ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^*\in {\mathbb{R}}_{ > 0}^n$?

3. Given vectors ${\boldsymbol{w}}^1, \dots, {\boldsymbol{w}}^n, \widetilde{{{\boldsymbol{w}}}}^1, \dots, \widetilde{{{\boldsymbol{w}}}}^n \in {\mathbb{R}}^d$, when is the generalized polynomial map on ${\mathbb{R}}_{ > 0}^d$ defined by

^*Thereby, ${\boldsymbol{x}} \circ {\boldsymbol{y}}$ denotes the component-wise product of the vectors ${\boldsymbol{x}}$ and ${\boldsymbol{y}}$; see Section 1.1.

bijective onto the relative interior of the polyhedral cone generated by ${\boldsymbol{w}}^1, \dots, {\boldsymbol{w}}^n$, for any ${\boldsymbol{x}}^*\in {\mathbb{R}}_{ > 0}^n$?

These questions will be expanded upon and explained in detail in Section 4.

Among the questions above, we initially focus on question 2, which is strongly related to Birch's theorem. One way to state Birch's theorem is: given a vector subspace $S \subseteq {\mathbb{R}}^n$, the intersection $({\boldsymbol{x}_0} + S) \cap ({\boldsymbol{x}}^* \circ \exp S^\perp)$ consists of exactly one point, for any ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$. In question 2, we have two vector subspaces $S$, $\widetilde{S}$, so it should not come as a surprise that an additional hypothesis is needed, in order for this intersection to consist of exactly one point.

This additional hypothesis is given in terms of sign vectors. For a subset $S\in {\mathbb{R}}^n$, its set of sign vectors $\sigma(S)$ is the image of vectors in $S$ under the coordinate-wise sign function. Its closure $\overline{{\sigma(S)}}$ contains $\sigma(S)$ and all sign vectors where one or more coordinates may be replaced by zeros (see Definition 5.1).

One of our main results is the following generalization of Birch's theorem:

Theorem 5.7. Let $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$ be vector subspaces of equal dimension with $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$. Then for any positive vectors ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$, the intersection $({\boldsymbol{x}_0} + S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$ consists of exactly one point.

By using this theorem, we obtain a sufficient condition for question 1 in Theorem 5.8. More precisely, provided that certain conditions hold, we show that if a generalized mass-action system has at least one vertex-balanced steady state, then there is exactly one vertex-balanced steady state within every invariant affine subspace.

We introduce generalized mass-action systems and vertex-balanced steady states in Section 2 and 3 respectively. We prove Theorem 5.7 and Theorem 5.8 in Section 5, and conclude with an example in Section 6.

1.1.   Notation

There are several component-wise operations on vectors and matrices that will appear frequently. In the list below, let ${\boldsymbol{x}}$, ${\boldsymbol{z}} \in {\mathbb{R}}^n$ with ${\boldsymbol{x}} = (x_1, x_2, \dots, x_n)^T$ and ${\boldsymbol{z}} = (z_1, z_2, \dots, z_n)^T$. Let $Y = ({\boldsymbol{y}_1}, {\boldsymbol{y}_2}, \cdots, {\boldsymbol{y}_m})$ be a $n\times m$ matrix.

We write ${\boldsymbol{x}} \geq 0$ to mean that every component of the vector is non-negative. Similarly, ${\boldsymbol{x}} > 0$ means that every component of the vector is positive. We let ${\mathbb{R}}_{\geq 0}^n = \{ {\boldsymbol{x}} \in {\mathbb{R}}^n : {\boldsymbol{x}} \geq 0\}$, and ${\mathbb{R}}_{ > 0}^n = \{ {\boldsymbol{x}} \in {\mathbb{R}}^n : {\boldsymbol{x}} > 0\}$. We denote the cardinality of a set $M$ as $|M|$.

The vector and matrix operations we will use are:

When the above operations are applied to a subset of ${\mathbb{R}}^n$, they are applied to elements of the set. For example, given a set $S \subseteq {\mathbb{R}}^n$, we have $\exp(S) = \{ \exp({\boldsymbol{x}}) : {\boldsymbol{x}} \in S\}$, and ${\boldsymbol{x}} \circ S = \{ {\boldsymbol{x}} \circ {\boldsymbol{z}} : {\boldsymbol{z}} \in S\}$.

2.   Generalized mass-action systems

Consider a simple directed graph $G = (V, E)$ and the corresponding weighted digraph $G_{\boldsymbol \kappa} = (V, E, \boldsymbol \kappa)$ with $\boldsymbol \kappa \in {\mathbb{R}}_{ > 0}^{E}$ providing a positive weight for each edge in $E$. Let $V = \{v_1, v_2, \dots, v_m\}$ be the set of vertices. Given an edge $e = v_i {\rightarrow} v_j \in E$, we call $v_i$ the source of $e$, and $v_j$ its target. Let us denote by $V_s \subseteq V$ the set of source vertices, that is, the set of vertices that are sources of some edges. The weight $\kappa_e > 0$ on the edge $e = v_i {\rightarrow} v_j$ is called a rate constant, and we refer to the vector $\boldsymbol \kappa \in {\mathbb{R}}_{ > 0}^{E}$ as the vector of rate constants, or more simply as the rate constants. Often, we use the indices of the source and target vertices as edge label, i.e., $\kappa_{v_i {\rightarrow} v_j} = \kappa_{ij}$.

Let $\Phi: V \to {\mathbb{R}}^n$ be a map assigning to each vertex $v \in V$ a stoichiometric complex $\Phi(v) \in {\mathbb{R}}^n$, and let $\widetilde{\Phi}: V_s \to {\mathbb{R}}^n$ be another map that assigns to each source vertex $v \in V_s$ a kinetic-order complex $\widetilde{\Phi}(v) \in {\mathbb{R}}^n$. An edge $v_i {\rightarrow} v_j$ is called a reaction, and the vector $\Phi(v_j) - \Phi(v_i)$ is the reaction vector associated to the edge $v_i {\rightarrow} v_j$. For convenience, we often write ${\boldsymbol{y}_i}$ instead of $\Phi(v_i)$, and $\widetilde{{{\boldsymbol{y}}}}_i$ instead of $\widetilde{\Phi}(v_i)$. The graph $G$ and the two maps $\Phi$, $\widetilde{\Phi}$ on $G$ provide all the ingredients needed to define a generalized reaction network, while the weighted digraph $G_{\boldsymbol \kappa}$ and the maps $\Phi$, $\widetilde{\Phi}$ are all that is needed to define a generalized mass-action system.

Definition 2.1. A generalized reaction network is given by $(G, \Phi, \widetilde{\Phi})$, where $G = (V, E)$ is a simple directed graph, and $\Phi: V \to {\mathbb{R}}^n$, $\widetilde{\Phi}: V_s \to {\mathbb{R}}^n$ respectively assign to each vertex a stoichiometric complex and to each source vertex a kinetic-order complex.

Remark. We follow the definition of a generalized reaction network given by Müller and Regensburger in ^[6], rather than the one given in ^[7]. In particular, we do not assume that the maps $\Phi$ and $\widetilde{\Phi}$ are injective.

Remark. Throughout this paper, we are concerned with generalized reaction networks where $V_s = V$. The digraphs $\Phi(G)$ and $\widetilde{\Phi}(G)$ are two Euclidean embedded graphs ^[17,18,19]. One of the equivalent definitions of a (classical) reaction network is a directed graph $\mathcal{G} = (\mathcal{V}, \mathcal{E})$, where the set $\mathcal{V}$ of vertices (complexes) is a subset of ${\mathbb{R}}^n$. Using the notation above, a reaction network is given by $\mathcal{G} = \Phi(G)$, where $\Phi$ is injective ^[6].

Example 2.2. To illustrate the terminology above, we consider a directed graph $G = (V, E)$ and the corresponding weighted digraph $G_{\boldsymbol \kappa} = (V, E, \boldsymbol \kappa)$:

The set of vertices is $V = \{ v_1, v_2, v_3, v_4, v_5\}$, which coincides with the set of source vertices $V_s$. The set of edges is $E = \{v_1 {\rightarrow} v_2, \, v_2 {\rightarrow} v_1, \, v_3 {\rightarrow} v_4, \, v_4 {\rightarrow} v_3, \, v_4 {\rightarrow} v_5, \, v_5 {\rightarrow} v_3 \}$. The maps $\Phi$ and $\widetilde{\Phi}$ (both from $V$ to ${\mathbb{R}}^2$) are given in Figure 1.

Note that the vertex $v_3$ is mapped differently by $\Phi$ and $\widetilde{\Phi}$. Indeed, $v_3$ is mapped by $\Phi$ to the stoichiometric complex $(1, 2)^T$ and by $\widetilde{\Phi}$ to the kinetic-order complex $(3, 1)^T$. Further, note that $\Phi$ maps the vertices $v_2$ and $v_3$ to the same stoichiometric complex, whereas $\widetilde{\Phi}$ maps $v_2$ and $v_3$ to different kinetic-order complexes. Hence, the number of connected components and the number of vertices are different in the graphs $\Phi(G)$ and $\widetilde{\Phi}(G)$.

Now we are in a position to define generalized mass-action systems and the associated dynamical systems.

Definition 2.3. A generalized mass-action system is given by $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$, where $(G, \Phi, \widetilde{\Phi})$ is a generalized reaction network, with directed graph $G = (V, E)$, and $\boldsymbol \kappa \in {\mathbb{R}}_{ > 0}^{E}$ is a vector of rate constants.

Definition 2.4. For a generalized mass-action system $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$, the associated dynamical system on ${\mathbb{R}}_{ > 0}^n$ is given by

As the ODE system (2.1) is our main object of interest, we pause to make two observations. First, the rate of change $\frac{d{\boldsymbol{x}}}{dt}$ is restricted to the stoichiometric subspace $S = \text{span}_ {\mathbb{R}} \{ {\boldsymbol{y}_j} - {\boldsymbol{y}_i} : v_i {\rightarrow} v_j \in E\}$. Consequently, every trajectory ${\boldsymbol{x}}(t)$ of this dynamical system is restricted to a stoichiometric compatibility class ${\boldsymbol{x}}(0) + S$. Second, if $v_i {\rightarrow} v_j$ is a reaction and $\Phi(v_i) = \Phi(v_j)$, then this particular reaction does not contribute to the dynamics.

It is sometimes more convenient to write the ODE system (2.1) in matrix form. Let $Y \in {\mathbb{R}}^{n\times m}$ be the stoichiometric complex matrix, the $j$-th column of which is the stoichiometric complex ${\boldsymbol{y}_j}$. Let the kinetic-order complex matrix $\widetilde{Y} \in {\mathbb{R}}^{n\times m}$ be defined analogously; in particular, its $j$-th column is the kinetic-order complex ${\boldsymbol{{\widetilde{y}}}}_j$ if $v_j \in V_s$ and ${\boldsymbol{0}}$ if $v_j \not\in V_s$.^† Let $A_{\boldsymbol \kappa} \in {\mathbb{R}}^{m\times m}$ be the negative transpose of the Laplacian of the weighted directed graph $G_{\boldsymbol \kappa}$, i.e.,

^† The choice of $\widetilde{{{\boldsymbol{y}}}}_j = {\boldsymbol{0}}$ when $v_j\not\in V_s$ is arbitrary, since the $j$-th column of $\widetilde{Y}$ does not appear in the equations that are of interest to us ^[6]. In particular, it does not affect the vector $A_{\boldsymbol \kappa} {\boldsymbol{x}}^{\widetilde{Y}}$ and hence does not contribute to the right-hand side $YA_{\boldsymbol \kappa} {\boldsymbol{x}}^{\widetilde{Y}}$ of the system of differential equations (2.3).

The dynamical system (2.1) can be rewritten as

Example 2.5. Returning to Example 2.2, the dynamical system associated to $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$ is

where each term corresponds to an edge in the graph $G$. Expanding the equations, we recognize it to be a polynomial (more generally, a power-law) dynamical system:

Its stoichiometric subspace is $S = {\mathbb{R}}^2$. The stoichiometric complex matrix and kinetic-order complex matrix are

respectively. The matrix

is the negative transpose of the Laplacian of the weighted digraph $G_{\boldsymbol \kappa}$.

The definitions above and Example 2.2 are relatively abstract; one may wonder how generalized mass-action systems show up in applications. Suppose we are interested in modelling the following chemical system:

Let us assume that based on experimental data, the reaction rate functions are as shown below, with rate constants $\kappa_{ij} > 0$:

Note that the third reaction follows power-law kinetics, but not classical mass-action kinetics. Moreover, note that the second and third reactions are an example of branching reactions, that is, two reactions with the same source complex, different target complexes, and (most importantly) with different kinetics: mass-action for the second reaction and power-law for the third. Exactly this chemical system is specified as a generalized mass-action system $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$ in Example 2.2, see Figure 1(a) for the reaction network and Figure 1(b) for the reaction rate functions. The system of ordinary differential equations modelling this system is precisely (2.4).

Remark. We defined a generalized reaction network as a triple $(G, \Phi, \widetilde{\Phi})$. As pointed out in an earlier remark, if $\Phi$ is injective, then $\Phi(G)$ is a (classical) reaction network. A classical mass-action system can be obtained as a special case of a generalized mass-action system $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$, where $\Phi$ is injective and $\widetilde{\Phi} = \left. \Phi \right|_{V_s}$ ^[7]. It is thus natural to extend some of the standard definitions for classical mass-action systems to generalized mass-action systems.

We say the underlying graph $G$ is weakly reversible if every connected component of $G$ is strongly connected, i.e., every edge is part of a directed cycle. We have already defined the stoichiometric subspace $S$ as the span of reaction vectors. Whenever $V_s = V$ (in particular, when $G$ is weakly reversible), we define its kinetic analogue, the kinetic-order subspace $\widetilde{S} = \text{span}_ {\mathbb{R}} \{ \widetilde{{{\boldsymbol{y}}}}_j - \widetilde{{{\boldsymbol{y}}}}_i : v_i {\rightarrow} v_j \in E\}$.

The stoichiometric deficiency of the generalized reaction network $(G, \Phi, \widetilde{\Phi})$ is the non-negative integer

where $|V|$ is the number of vertices in $G$, $\ell_G$ is the number of connected components of $G$, and $S$ is the stoichiometric subspace. From the equivalent definition $\delta_G = \dim(\text{ker} Y \cap \text{im} I_E)$, where $I_E$ is the incidence matrix of $G$, it follows that $\delta_G$ is a non-negative integer ^[8]. In the case when $G$ is weakly reversible, we also have the formula ^[12]

Whenever $V_s = V$, the kinetic-order deficiency is defined as the non-negative integer

where $\widetilde{S}$ is the kinetic-order subspace.

Remark. In the definitions above, $|V|$ is the number of vertices in the underlying abstract graph $G$, not necessarily the number of distinct stoichiometric complexes or the number of kinetic-order complexes; $\ell_G$ is the number of connected components of $G$, not necessarily the number of connected components in $\Phi(G)$ or $\widetilde{\Phi}(G)$.

In Example 2.2, we have a weakly reversible network with $|V| = 5$ vertices and $\ell_G = 2$ connected components. We already observed that the stoichiometric subspace $S$ is all of ${\mathbb{R}}^2$. However, the kinetic-order subspace is $\widetilde{S} = \text{span}_ {\mathbb{R}} (1, 1)^T$. The stoichiometric deficiency in this example is $\delta_G = 5 - 2 -2 = 1$, but the kinetic-order deficiency is $\widetilde{\delta}_G = 5 - 2 - 1 = 2$.

Example 2.6. We have seen earlier that generalized mass-action systems arise naturally from power-law kinetics. This example illustrates how generalized mass-action systems also arise naturally in the study of mass-action systems, via a process called network translation ^[8,20]. Network translation produces a generalized mass-action system that has the same dynamics as the original mass-action system. We look at the $n$-site distributive phosphorylation-dephosphorylation system under mass-action kinetics.

This example first appeared in ^[8]; below we consider a different translation of the same mass-action system. Under the original definition of generalized mass-action system in ^[7], which requires the stoichiometric complex map $\Phi$ and the kinetic-order complex map $\widetilde{\Phi}$ to be injective, the translated network presented below would not have been a well-defined generalized reaction network. However, the later definition in ^[6] removes the requirements that $\Phi$ and $\widetilde{\Phi}$ are injective. As a result, many more dynamical systems can be written as a generalized mass-action system, and for this example, a more natural translation exists for the $n$-site distributive phosphorylation-dephosphorylation system.

Let ${\rm{E}}$, $\rm{F}$ be enzymes that catalyze the phosphorylation and dephosphorlyation processes, by forming intermediates $\rm{ES_j}$ and $\rm{FS_j}$ respectively. The $n$-site distributive phosphorylation-dephosphorylation system consists of the following reactions:

We assume that the reaction rates follow classical mass-action kinetics. There are $3n+3$ species involved, so the system of differential equations modelling their concentrations is defined on ${\mathbb{R}}_{ > 0}^{3n+3}$. We create a generalized mass-action system with the same differential equations by network translation. The main step involves changing the stoichiometric complexes: adding enzyme $\rm{F}$ to the series of reactions for phosphorylation by ${\rm{E}}$; and adding enzyme ${\rm{E}}$ to the series of reactions for dephosphorylation by $\rm{F}$. This process produces a weakly reversible network:

To define the generalized mass-action system, we take a more top-down approach, starting from a graph $G$ with $n$ components and $4n$ vertices:

Although $\Phi$ and $\widetilde{\Phi}$ map vertices to vectors in ${\mathbb{R}}^{3n+3}$, to make this example more readable, we will specify the images of $\Phi$ and $\widetilde{\Phi}$ in terms of formal linear combination of species.

The stoichiometric complexes are

and the kinetic-order complexes are

Note that the map $\Phi$ is not injective, as $\Phi(v_j) = \Phi(v_j')$ for $1 \leq j \leq n-1$. The image of the graph $G$ under $\Phi$ is connected:

One can check that $\dim S = \dim \widetilde{S} = 3n$; therefore, the stoichiometric deficiency and kinetic-order deficiency are $\delta_G = \widetilde{\delta}_G = (4n) - (n) - (3n) = 0$.

3.   Vertex-balanced steady states

Given the dynamical system associated to a generalized mass-action system $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$, written either as $\frac{d{\boldsymbol{x}}}{dt} = \sum_{v_i {\rightarrow} v_j\in E} \kappa_{ij} {\boldsymbol{{x}}}^{\widetilde{{{\boldsymbol{y}}}}_i} ({\boldsymbol{y}_j} - {\boldsymbol{y}_i})$ or in matrix form $\frac{d{\boldsymbol{x}}}{dt} = YA_{\boldsymbol \kappa} {\boldsymbol{x}}^{\widetilde{Y}}$, it is natural to ask how many steady states there are. We define the set of positive steady states as

For a classical mass-action system, an important subset of positive steady states is the set of complex-balanced equilibria ^[2], also known as complex-balancing equilibria or vertex-balanced equilibria ^[27]. Horn and Jackson introduced the idea of complex balancing at equilibrium to generalize the physical assumption of detailed balancing at thermodynamic equilibrium ^[2].

We illustrate the intuition behind the definition of such a steady state before introducing its analogue for a generalized mass-action system. Consider the graph $G$ of the reaction network, and associate to each edge $v_i {\rightarrow} v_j$ a reaction rate function $\kappa_{ij} {\boldsymbol{x}}^{{\boldsymbol{y}_i}}$. A concentration vector ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$ is a complex-balanced equilibrium of the classical mass-action system if at every vertex $v_i \in V$ of the graph, the sum of incoming fluxes balances the sum of outgoing fluxes, that is, for all $v_i \in V$,

This occurs if and only if $A_{\boldsymbol \kappa} ({\boldsymbol{x}}^*)^{Y} = {\boldsymbol{0}}$ ^[2]. Clearly, a complex-balanced equilibrium is a positive solution to a system of polynomial equations. Surprisingly, it is also a positive solution to a system of binomial equations ^[13].

For a generalized mass-action system, one can define a vertex-balanced steady state analogously: it is a positive steady state at which the net flux is zero across every vertex of the graph, where the flux is given by generalized mass-action kinetics.

Definition 3.1. The set of vertex-balanced steady states for a generalized mass-action system $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$ is the set

Note that ${\boldsymbol{x}}^* \in Z_{\boldsymbol \kappa}$ if and only if for all $v_i \in V$,

Remark. What we call vertex-balanced steady state here, is also called complex balancing equilibrium ^[6,7] or generalized complex-balanced steady state ^[8].

We call such a steady state vertex-balanced instead of complex-balanced to avoid a subtle point of confusion. In the case when $\Phi$ is not injective, the balancing of in-fluxes and out-fluxes occurs at each vertex $v \in V$ of the underlying abstract graph $G$. This in turn implies the balancing of fluxes at each stoichiometric complex $\Phi(v) \in \Phi(V)$; however, the converse is generally false. For example, consider the generalized mass-action system given by the weighted digraph $G_{\boldsymbol \kappa}$,

and the maps $\Phi(v_1) = \Phi(v_4) = 0 \in {\mathbb{R}}^1$, $\Phi(v_2) = \Phi(v_3) = 1$, $\widetilde{\Phi}(v_1) = 0$, and $\widetilde{\Phi}(v_3) = 1$. The associated dynamical system is $\frac{dx}{dt} = \kappa - \kappa x$, which also arises from the classical mass-action system $\Phi(G_{\boldsymbol \kappa})$:

However, the equilibrium $x^* = 1$ is not vertex-balanced for the generalized mass-action system $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$, but is complex-balanced for the classical mass-action system $\Phi(G_{\boldsymbol \kappa})$.

Example 3.2. To illustrate the definition of vertex-balanced steady states, we consider Example 2.2 again. A vertex-balanced steady state is a point ${\boldsymbol{x}} = (x_1, x_2)^T \in {\mathbb{R}}_{ > 0}^2$ satisfying five polynomial equations, one equation for each vertex of the graph $G$:

However, a positive steady state ${\boldsymbol{x}} = (x_1, x_2)^T \in {\mathbb{R}}_{ > 0}^2$ has to satisfy only two polynomial equations:

The two polynomial equations (3.6) are linear combinations of the five polynomial equations (3.5); thus $Z_{\boldsymbol \kappa} \subseteq E_{\boldsymbol \kappa}$. This follows from the matrix expression of the associated dynamical system $\frac{d{\boldsymbol{x}}}{dt} = Y(A_{\boldsymbol \kappa} {\boldsymbol{x}}^{\widetilde{Y}})$.

Complex-balanced equilibria of classical mass-action systems have been studied extensively. Some of the classical results extend directly to the case of generalized mass-action systems, even when the maps $\Phi$ and $\widetilde{\Phi}$, assigning stoichiometric complexes and kinetic-order complexes respectively, are not injective. For example, it is known that ^[6,7]: ^‡

^‡Some of these results were first proved in ^[7], under the assumption that $\Phi$ and $\widetilde{\Phi}$ are injective, but the same proof goes through without these hypotheses.

ⅰ) If $Z_{\boldsymbol \kappa} \neq \varnothing$ for some ${\boldsymbol \kappa} > 0$, then the underlying graph $G$ is weakly reversible.

ⅱ) If $Z_{\boldsymbol \kappa} \neq \varnothing$ and ${\boldsymbol{x}}^* \in Z_{\boldsymbol \kappa}$, then $Z_{\boldsymbol \kappa} = \{ {\boldsymbol{x}} \in {\mathbb{R}}_{ > 0}^n : \ln {\boldsymbol{x}}- \ln {\boldsymbol{x}}^* \in \widetilde{S}^\perp\} = {\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp$.

ⅲ) For a weakly reversible generalized reaction network, $\widetilde{\delta}_G = 0$ if and only if $Z_{\boldsymbol \kappa} \neq \varnothing$ for any choices of rate constants ${\boldsymbol \kappa} > 0$.

ⅳ) For a weakly reversible generalized reaction network, if $\delta_G = 0$, then for any choice of rate constants ${\boldsymbol \kappa} > 0$, any positive steady state is a vertex-balanced steady state, i.e., $E_{\boldsymbol \kappa} = Z_{\boldsymbol \kappa}$.

In the example of the $n$-site phosphorlyation-dephosphorlyation system (Example 2.6), we noted that $\delta_G = \widetilde{\delta}_G = 0$. By statements (ⅲ) and (ⅳ) above, we conclude for any rate constants ${\boldsymbol \kappa}$, the set of vertex-balanced steady states $Z_{\boldsymbol \kappa}$ is non-empty, and all positive steady states are vertex-balanced. Moreover, the set of positive steady states is given by $E_{\boldsymbol \kappa} = Z_{\boldsymbol \kappa} = {\boldsymbol{x}}^* \circ \exp \widetilde{{S}}^\perp$, where ${\boldsymbol{x}}^*$ is any positive steady state and $\widetilde{{S}}$ is the kinetic-order subspace, i.e., the vector space spanned by the differences of kinetic-order complexes according to the edges in the graph. It should be noted that the $n$-site phosphorlyation-dephosphorlyation system is multistationary when $n \geq 2$ ^[21,22], i.e., the system admits multiple steady states within the same stoichiometric compatibility class. In other words, for some choices of rate constants, there are multiple vertex-balanced steady states within some stoichiometric compatibility class. This contrasts with a classical complex-balanced mass-action system, where $Z_{\boldsymbol \kappa}$ meets every stoichiometric compatibility class at most once.

Example 3.3. There is another surprising way in which vertex balancing differs from classical complex-balanced mass-action systems. While it is clear that $Z_{\boldsymbol \kappa} \subseteq E_{\boldsymbol \kappa}$, in generalized mass-action systems it is possible that $\varnothing \neq Z_{\boldsymbol \kappa} \subsetneq E_{\boldsymbol \kappa}$^§. For example, consider the 1-dimensional generalized mass-action system given by the weighted digraph $G_{\boldsymbol \kappa}$

^§For classical complex-balanced mass-action systems, it is always the case that $Z_{\boldsymbol \kappa} = E_{\boldsymbol \kappa}$ ^[2].

and the maps $\Phi$, $\widetilde{{\Phi}}$ given by

The associated dynamical system is

The point $x^* = 1$ is a vertex-balanced steady state. The system also has two other steady states $x^* \approx 0.27$ and $x^* \approx 3.72$, neither of which satisfy the vertex-balanced condition:

In applications, the vector of rate constants ${\boldsymbol \kappa} \in {\mathbb{R}}_{ > 0}^E$ is often not known precisely. Surprisingly, some important results for complex-balanced equilibria in classical mass-action systems hold irrespective of the precise values of the rate constants. We are interested in results for vertex-balanced equilibria of generalized mass-action systems that are in this sense independent of the choice of rate constants. We have observed that the solution trajectories are confined to a stoichiometric compatibility class ${\boldsymbol{x}_0} + S$, where ${\boldsymbol{x}_0} \in {\mathbb{R}}_{ > 0}^n$ is an initial state and $S$ is the stoichiometric subspace. Therefore, our object of study is the intersection $({\boldsymbol{x}_0} + S) \cap Z_{\boldsymbol \kappa}$ for any ${\boldsymbol{x}_0} \in {\mathbb{R}}_{ > 0}^n$ and any ${\boldsymbol \kappa} \in {\mathbb{R}}_{ > 0}^E$.

4.   Problem reformulations

In the introduction, we have mentioned that the following questions are essentially equivalent:

1. When does a generalized mass-action system have exactly one vertex-balanced steady state within each stoichiometric compatibility class, for any choice of rate constants?

2. Given vector subspaces $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$, when does the intersection $({\boldsymbol{x}_0} + S)\cap ({\boldsymbol{x}}^* \circ \exp\widetilde{S}^\perp)$ consist of exactly one point, for any ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^*\in {\mathbb{R}}_{ > 0}^n$?

3. Given vectors ${\boldsymbol{w}}^1, \dots, {\boldsymbol{w}}^n, \widetilde{{{\boldsymbol{w}}}}^1, \dots, \widetilde{{{\boldsymbol{w}}}}^n \in {\mathbb{R}}^d$, when is the generalized polynomial map on ${\mathbb{R}}_{ > 0}^d$ defined by

bijective onto the relative interior of the polyhedral cone generated by ${\boldsymbol{w}}^1, \dots, {\boldsymbol{w}}^n$, for any ${\boldsymbol{x}}^*\in {\mathbb{R}}_{ > 0}^n$?

Before we discuss the relationship between these problems in detail, let us first make a historical note. When speaking of a weakly reversible classical mass-action system, Horn and Jackson ^[2] proved that if the system has at least one complex-balanced equilibrium, then every stoichiometric compatibility class has exactly one complex-balanced equilibrium. Indeed, they showed that every positive steady state of such a system is complex-balanced and locally asymptotically stable within its stoichiometric compatibility class. A complex-balanced equilibrium is globally stable within its stoichiometric compatibility class when the network has a single connected component ^[23], or is strongly endotactic ^[24], or when the system is in ${\mathbb{R}}^3$ ^[25,26]. A general proof of global stability of complex-balanced equilibrium within its stoichiometric compatibility class was proposed for all complex-balanced systems in ^[27].

The first of the three questions above is phrased in the context of reaction networks. We start with a generalized reaction network and suppose that for some rate constants ${\boldsymbol \kappa}$, there is a vertex-balanced steady state ${\boldsymbol{x}}^* \in Z_{\boldsymbol \kappa}$. What is a condition (E) on the network $(G, \Phi, \widetilde{\Phi})$ for the existence of a vertex-balanced steady state within every stoichiometric compatibility class? What is a condition (U) on $(G, \Phi, \widetilde{\Phi})$ so that a vertex-balanced steady state is unique within its stoichiometric compatibility class? We would like to obtain conditions for these to hold or fail that are independent of the rate constants ${\boldsymbol \kappa}$. More precisely:

Problem 1. Let $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$ be a generalized mass-action system. Suppose that $Z_{\boldsymbol \kappa} \neq \varnothing$. What are conditions (E) and (U) on $(G, \Phi, \widetilde{\Phi})$, so that the following statements are true?

1. If $(G, \Phi, \widetilde{\Phi})$ satisfies condition (E), then there is at least one vertex-balanced steady state in every stoichiometric compatibility class, i.e., $({\boldsymbol{x}_0} +S) \cap Z_{\boldsymbol \kappa}$ contains at least one point for any ${\boldsymbol{x}_0} \in {\mathbb{R}}_{ > 0}^n$.

2. If $(G, \Phi, \widetilde{\Phi})$ satisfies condition (U), then there is at most one vertex-balanced steady state in every stoichiometric compatibility class, i.e., $({\boldsymbol{x}_0} + S) \cap Z_{\boldsymbol \kappa}$ contains at most one point for any ${\boldsymbol{x}_0} \in {\mathbb{R}}_{ > 0}^n$.

Recall that $Z_{\boldsymbol \kappa} = {\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp$ for any ${\boldsymbol{x}}^* \in Z_{\boldsymbol \kappa}$. Thus, the vertex-balanced steady states within any stoichiometric compatibility class ${\boldsymbol{x}_0} +S$ belong to the intersection $({\boldsymbol{x}_0} + S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$. This leads us to the following reformulation of Problem 1:

Problem 2. Let $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$ be vector subspaces. What are conditions (E) and (U) on $S$, $\widetilde{S}$, so that the following statements are true?

1. If $S$, $\widetilde{S}$ satisfy condition (E), then $({\boldsymbol{x}_0} + S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$ contains at least one point, for any ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$.

2. If $S$, $\widetilde{S}$ satisfy condition (U), then $({\boldsymbol{x}_0} + S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$ contains at most one point, for any ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$.

If a generalized mass-action system happens to be a classical mass-action system, then its stoichiometric subspace $S$ is also the kinetic-order subspace $\widetilde{S}$. The existence and uniqueness of a point in the intersection $({\boldsymbol{x}_0} + S) \cap ({\boldsymbol{x}}^* \circ \exp S^\perp)$ for any ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$ is the content of Birch's theorem in algebraic statistics ^[15].

Another reformulation of the above problems was introduced by Müller and Regensburger ^[7], in terms of injectivity/surjectivity of an exponential map or a generalized polynomial map onto a polyhedral cone. Such polynomial maps appear in other applications; for example, a renormalized version of the generalized polynomial appears in computer graphics and geometric modelling, where the map being injective implies that the curve or surface does not self-intersect ^[16].

Let ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$ be an arbitrary vector, and $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$ be vector subspaces, with $d = \text{codim} S$, $\widetilde{d} = \text{codim} \widetilde{S}$. Choose a basis for $S^\perp$ and let the basis vectors be the rows of the matrix $W \in {\mathbb{R}}^{d\times n}$. Similarly, choose a basis for $\widetilde{S}^\perp$ and let the basis vectors be the rows of $\widetilde{W} \in {\mathbb{R}}^{\widetilde{d} \times n}$. Write the two matrices in terms of their columns: $W = ({\boldsymbol{w}}^1, {\boldsymbol{w}}^2, \cdots, {\boldsymbol{w}}^n)$ and $W = (\widetilde{{{\boldsymbol{w}}}}^1, \widetilde{{{\boldsymbol{w}}}}^2, \cdots, \widetilde{{{\boldsymbol{w}}}}^n)$. In this manner, $S^\perp = \text{im} (W^T)$, $S = \text{ker} W$, and $\widetilde{S}^\perp = \text{im} (\widetilde{W}^T)$, $\widetilde{S} = \text{ker} \widetilde{W}$. Finally, write ${{C}}^0(W)$ for the relative interior of the polyhedral cone $C(W)$, i.e., ${{C}}^0(W)$ is the set of all positive combinations of $\{{\boldsymbol{w}}^i \}_{i = 1}^n$. For any ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$, define the maps

and

Problem 2 is equivalent to the following (see ^[6,7] for details):

Problem 3. Let $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$ be vector subspaces. What are conditions (E) and (U) on $S$, $\widetilde{S}$, so that the following statements are true?

1. If $S$, $\widetilde{S}$ satisfy condition (E), then the map $f_{{\boldsymbol{x}}^*}$ (respectively $F_{{\boldsymbol{x}}^*}$) is surjective onto $C^0(W)$, for any ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$.

2. If $S$, $\widetilde{S}$ satisfy condition (U), then the map $f_{{\boldsymbol{x}}^*}$ (respectively $F_{{\boldsymbol{x}}^*}$) is injective, for any ${\boldsymbol{x}}^*\in {\mathbb{R}}_{ > 0}^n$.

Müller and Regensburger characterized when the maps $f_{{\boldsymbol{x}}^*}$, $F_{{\boldsymbol{x}}^*}$ are injective, namely, if and only if $\sigma(S) \cap \sigma(\widetilde{S}^\perp) = \{{\boldsymbol{0}}\}$ [7,Theorem 3.6]. Recall that, for a subset $S \subseteq {\mathbb{R}}^n$, its set of sign vectors $\sigma(S)$ is the image of vectors in $S$ under the coordinate-wise sign function (Definition 5.1). They also provided a sufficient condition for bijectivity: if $\sigma(S) = \sigma(\widetilde{S})$ and $(+, +, \cdots, +)^T \in \sigma(S^\perp)$, then $f_{{\boldsymbol{x}}^*}$, $F_{{\boldsymbol{x}}^*}$ are bijective (and indeed, real analytic isomorphisms) [7,Proposition 3.9]. Our main result (Theorem 5.8) can be regarded as a generalization of this result. Recently, Müller, Hofbauer, and Regensburger have used Hadamard's global inversion theorem to characterize when $f_{{\boldsymbol{x}}^*}$, $F_{{\boldsymbol{x}}^*}$ are bijective for arbitrary ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$ ^[28].

5.   Main result

In previous work as well as in ours, the conditions (E) and (U) are stated in terms of sign vectors. For a brief introduction to sign vectors of linear subspaces, we refer the reader to the appendix in ^[7].

Definition 5.1. Given a vector ${\boldsymbol{x}} \in {\mathbb{R}}^n$, we define its sign vector to be

The set of sign vectors for a subset $S \subseteq {\mathbb{R}}^n$ is the collection $\sigma(S) = \{ \sigma({\boldsymbol{x}}) : {\boldsymbol{x}} \in S\}$.

We introduce a partial order on $\{0, +, -\}$ by $0 < -$ and $0 < +$ (but no relation between $-$ and $+$). This induces a partial order on $\{0, +, -\}^n$: $\tau \leq \tau'$ if $\tau_j \leq \tau'_j$ for all $j$. The closure of a set of sign vectors $\Lambda \subseteq \{0, +, -\}^n$ is the set

We define an orthant of ${\mathbb{R}}^n$ to be a maximal subset of ${\mathbb{R}}^n$ on which $\sigma$ is constant. Geometrically, the sign vector $\sigma({\boldsymbol{x}})$ tells us which orthant $\mathcal{O}_{{\boldsymbol{x}}}$ the vector ${\boldsymbol{x}}$ lies in, while the closure $\overline{{\sigma({\boldsymbol{x}})}}$ refers to the union of $\mathcal{O}_{{\boldsymbol{x}}}$ and its boundary. Finally, we define an orthogonality relation on $\{0, +, - \}^n$; we say that two sign vectors $\tau$ and $\tau'$ are orthogonal (denoted $\tau \perp \tau'$) if

This differs from the typical definition of an orthant of ${\mathbb{R}}^n$, which is full dimensional.

either $\tau_j\cdot \tau_j' = 0$ for all $1 \leq j \leq n$

or there exist indices $i, j$ such that $\tau_i\cdot \tau_i' = +$ and $\tau_j\cdot \tau_j' = -$,

where the product operation on signs is as one would expect:

It is easy to see that if ${\boldsymbol{x}}$, ${\boldsymbol{y}} \in {\mathbb{R}}^n$ are orthogonal vectors, then $\sigma({\boldsymbol{x}}) \perp \sigma({\boldsymbol{y}})$.

We show in this section that if $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$ and $\dim S = \dim \widetilde{S}$, then for any ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$, the intersection $({\boldsymbol{x}_0} +S) \cap ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$ contains exactly one point. The intuitive idea is that the sign condition $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$ is related to a transversal intersection of the two manifolds $({\boldsymbol{x}_0} + S)$ and $({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$. If we have one intersection point, say ${\boldsymbol{x}}^* \in ({\boldsymbol{x}}^*+S) \cap ({\boldsymbol{x}}^*\circ \exp\widetilde{S}^\perp)$, we cannot lose the intersection point as we translate the affine plane from $({\boldsymbol{x}}^* +S)$ to $({\boldsymbol{x}_0} +S)$.

We first show in Lemma 5.2 that our sign condition $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$ implies the uniqueness condition $\sigma(S) \cap \sigma(\widetilde{S}^\perp) = \{ {\boldsymbol{0}}\}$ in ^[7]. In Lemma 5.4, we establish transversality of the manifolds $({\boldsymbol{x}} + S)$ and $({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$. Lemma 5.3 prevents our desired intersection point from escaping to the boundary of ${\mathbb{R}}_{ > 0}^n$ or to infinity. Finally, these results lead to Theorem 5.7, concluding the existence and uniqueness of a point in the intersection $({\boldsymbol{x}_0} + S) \cap ({\boldsymbol{x}}^* \circ \exp\widetilde{S}^\perp)$. In Theorem 5.8 and Corollary 5.9, we apply this result to generalized mass-action systems.

Lemma 5.2 (Uniqueness). Let $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$ be vector subspaces. If $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$, then $\sigma(S) \cap \sigma(\widetilde{S}^\perp) = \{{\boldsymbol{0}}\}$. In particular, for any ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$ the intersection $({\boldsymbol{x}_0} + S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$ contains at most one point.

Proof. By assumption, $\sigma(S) \cap \sigma(\widetilde{S}^\perp) \subseteq \overline{{\sigma(\widetilde{S})}} \cap \sigma(\widetilde{S}^\perp)$. We show that $\overline{{\sigma(\widetilde{S})}} \cap \sigma(\widetilde{S}^\perp) = \{ {\boldsymbol{0}} \}$. Let $\tau \in \overline{{\sigma(\widetilde{S})}} \cap \sigma(\widetilde{S}^\perp)$ be a sign vector. There exist vectors ${\boldsymbol{x}} \in \widetilde{S}$ and ${\boldsymbol{y}} \in \widetilde{S}^\perp$ such that $\tau \leq \sigma({\boldsymbol{x}})$ and $\tau = \sigma({\boldsymbol{y}})$. It is easy to see that if $\tau \leq \sigma({\boldsymbol{x}})$, and $\tau \perp \sigma({\boldsymbol{x}})$, then $\tau = {\boldsymbol{0}}$.

By ^[7], $\sigma(S) \cap \sigma(\widetilde{S}^\perp) = \{{\boldsymbol{0}}\}$ is necessary and sufficient for the intersection $({\boldsymbol{x}_0}+S)\cap ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$ to contain at most one point for any ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$.

Lemma 5.3 (Compactness). Let $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$ be vector subspaces, and let $K \subseteq {\mathbb{R}}_{ > 0}^n$ be a compact subset, and ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$. Suppose $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$. Then $(K+S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$ is a compact subset of ${\mathbb{R}}_{ > 0}^n$.

Proof. Let $\Gamma = (K+S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$. Since ${\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp \subseteq {\mathbb{R}}_{ > 0}^n$, the intersection $\Gamma$ also lies in the positive orthant. We first show that $\Gamma$ is bounded away from infinity and from the boundary of ${\mathbb{R}}_{ > 0}^n$.

Suppose that is not the case. Let ${\boldsymbol{x}}^k \in \Gamma$ be a sequence such that either $\limsup_{k\to\infty} x^k_i = \infty$ or $\liminf_{k\to\infty} x^k_i = 0$ for some index $1 \leq i \leq n$. Passing to a subsequence, we may assume that

where $I_1$, $I_2$, $I_3$ partition the index set $\{1, 2, \dots, n\}$, and $I_1 \cup I_2 \neq \varnothing$.

On one hand, ${\boldsymbol{x}}^k \in K+S$, so decompose it as ${\boldsymbol{x}}^k = {\boldsymbol{v}}^k + {\boldsymbol{s}}^k$, where ${\boldsymbol{v}}^k \in K$ and ${\boldsymbol{s}}^k \in S$. Since $K \subseteq {\mathbb{R}}_{ > 0}^n$ is compact, each component of ${\boldsymbol{v}}^k$ is uniformly bounded from above and below from zero. Thus for $i \in I_1$, we have $s^k_i \to \infty$; in particular, $s^k_i > 0$ for sufficiently large $k$. Similarly, if $i \in I_2$, then $s^k_i < 0$ for sufficiently large $k$, because $s^k_i + v^k_i \to 0$ and $v^k_i > 0$ is bounded away from zero. Hence the sign of $s^k_i$ is constant for any $i \in I_1 \cup I_2$ for sufficiently large $k$. Because $\sigma({\boldsymbol{s}}^k) \in \sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$, there is a vector ${\boldsymbol{u}} \in \widetilde{S}$ such that $u_i > 0$ if $i \in I_1$ and $u_i < 0$ if $i \in I_2$.

On the other hand, ${\boldsymbol{x}}^k \in {\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp$, that is, $\log \left(\frac{{\boldsymbol{x}}^k}{{\boldsymbol{x}}^*}\right) \in \widetilde{S}^\perp$, where the division is understood to be component-wise. Hence, ${\boldsymbol{u}} \perp \log \left(\frac{{\boldsymbol{x}}^k}{{\boldsymbol{x}}^*} \right)$ for all $k$, and we have

The sum over $I_3$ is uniformly bounded for all $k$. Now let $k\to \infty$. For $i \in I_1$, we know $u_i > 0$ and $x_i^k \to \infty$, so the sum over $I_1$ is positive and unbounded. For $i \in I_2$, we know $u_i < 0$ and $x_i^k \to 0$, so $\log\left(\frac{x_i^k}{x_i^*} \right) \to -\infty$, so the sum over $I_2$ is also positive and unbounded. Consequently, $0 = \lim_{k\to\infty} \langle{{{\boldsymbol{u}}}}, {{\log\left(\frac{{\boldsymbol{x}}^k}{{\boldsymbol{x}}^*} \right)}}\rangle = \infty$, a contradiction. Hence, $\Gamma \subseteq {\mathbb{R}}_{ > 0}^n$ is bounded away from infinity and away from the boundary of the positive orthant.

Next, we want to show that $\Gamma \subseteq {\mathbb{R}}_{ > 0}^n$ is a closed subset. Let us fix $\varepsilon > 0$ such that $\Gamma$ lies inside the hypercube $Q = [\varepsilon, \varepsilon^{-1}]^n \subseteq {\mathbb{R}}_{ > 0}^n$. Being the intersection of two closed sets, $Q \cap (K+S)$ is closed. The set $Q \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$ is diffeomorphic to $[\log \varepsilon, \log \varepsilon^{-1}]^n \cap (\ln {\boldsymbol{x}}^* + \widetilde{S}^\perp)$, which is again a closed set. Therefore, the set $(K+S)\cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp) = [Q \cap (K+S)] \cap [Q \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)]$ is the intersection of two closed sets, and thus it is closed in ${\mathbb{R}}_{ > 0}^n$.

Two manifolds $X$ and $Y$ of ${\mathbb{R}}^n$ intersect transversally if at each point ${\boldsymbol{p}} \in X \cap Y$, their tangent spaces span the entire Euclidean space, i.e., $T_{{\boldsymbol{p}}}(X) + T_{{\boldsymbol{p}}}(Y) = {\mathbb{R}}^n$. We refer the reader to ^[29,30] for the theory of transversality and intersection.

Again, let ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^*\in {\mathbb{R}}_{ > 0}^n$ be two arbitrary vectors in what follows. In Lemma 5.2, we showed that our sign condition $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$ implies $\sigma(S) \cap \sigma(\widetilde{S}^\perp) = \{{\boldsymbol{0}}\}$, which is equivalent to the intersection $({\boldsymbol{x}_0}+S) \cap ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$ containing at most one point. Indeed, this weaker sign condition together with $\dim S = \dim \widetilde{S}$ is enough to conclude that the two manifolds ${\boldsymbol{x}_0}+S$ and ${\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp$ intersect transversally. This is the content of the follow lemma.

Lemma 5.4 (Transversality). Let $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$ be vector subspaces. Assume $\sigma(S) \cap \sigma(\widetilde{S}^\perp) = \{ {\boldsymbol{0}}\}$. Let ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$ be any two positive vectors. Then the tangent spaces of ${\boldsymbol{x}_0} +S$ and ${\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp$ satisfy

for any point ${\boldsymbol{p}} \in ({\boldsymbol{x}_0} +S) \cap ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$.

If we further assume that $\dim S = \dim \widetilde{S}$, then $T_{\boldsymbol{p}}({\boldsymbol{x}_0} + S) + T_{\boldsymbol{p}}({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp) = {\mathbb{R}}^n$ for any intersection point ${\boldsymbol{p}} \in ({\boldsymbol{x}_0}+S) \cap ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$, i.e., ${\boldsymbol{x}_0} +S$ and ${\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp$ intersect transversally.

Proof. For any intersection point ${\boldsymbol{p}} \in ({\boldsymbol{x}_0} +S) \cap ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$, we note that $T_{\boldsymbol{p}}({\boldsymbol{x}_0} + S) = S$ and $T_{\boldsymbol{p}}({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp) = {\boldsymbol{p}} \circ \widetilde{S}^\perp$ and hence $\sigma(T_{\boldsymbol{p}}({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)) = \sigma(\widetilde{S}^\perp)$.

Now consider ${\boldsymbol{x}} \in T_{\boldsymbol{p}}({\boldsymbol{x}_0} + S) \cap T_{\boldsymbol{p}}({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$. Then $\sigma({\boldsymbol{x}}) \in \sigma(S) \cap \sigma(\widetilde{S}^\perp) = \{{\boldsymbol{0}}\}$, which implies ${\boldsymbol{x}} = {\boldsymbol{0}}$. Consequently, $T_{\boldsymbol{p}} ({\boldsymbol{x}_0} + S) \cap T_{\boldsymbol{p}}({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp) = \{ {\boldsymbol{0}}\}$.

If we further assume that $\dim S = \dim \widetilde{S}$, we note that $T_{\boldsymbol{p}} ({\boldsymbol{x}_0} +S) + T_{\boldsymbol{p}} ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$ is of dimension $n$. In other words, the manifolds ${\boldsymbol{x}_0}+S$ and ${\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp$ intersect transversally.

Now we are ready to state and prove our main result. The proof starts with a known intersection point, ${\boldsymbol{x}}^* \in ({\boldsymbol{x}}^* +S)\cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$. Next, we translate the affine space $({\boldsymbol{x}}^* + S)$ to $({\boldsymbol{x}_0} +S)$, creating a $(d+1)$-dimensional strip of the form $K+S$, where $d = \dim S$ and $K$ is a compact subset of ${\mathbb{R}}_{ > 0}^n$. This strip intersects ${\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp$ transversally, and we use Corollary 5.6 below to conclude that the intersection $(K+S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$ is a one-dimensional manifold, whose boundary lies on the boundary of the affine strip $K+S$. Finally, we conclude the existence of a boundary point on ${\boldsymbol{x}_0} +S$ by the uniqueness condition.

Consider the following differential topology result:

Theorem 5.5. ([30,Theorem 3.5.1]). Let $X$ and $Y$ be manifolds and $Z \subseteq Y$ a submanifold, where $Z$ and $Y$ are boundaryless. Let $f: X \to Y$ be a smooth map. Suppose $f$ intersects $Z$ transversally and $\left. f \right|_{ \partial X}$ also intersects $Z$ transversally. Then $f^{-1}(Z)$ is a submanifold of $X$ with boundary $\partial(f^{-1}(Z)) = \partial X \cap f^{-1}(Z)$, and $\text{codim}_{X}(f^{-1}(Z)) = \text{codim}_Y(Z)$.

Consider the setting where the ambient manifold is $Y = {\mathbb{R}}^n_{ > 0}$. If $f$ is the inclusion map of a submanifold $X$ into ${\mathbb{R}}^n_{ > 0}$, to say that the maps $f$ and $\left. f\right|_{ \partial X}$ intersect the manifold $Z$ transversally is equivalent to the manifolds $X$ and $\partial X$ intersect $Z$ transversally. The preimage $f^{-1}(Z)$ is the submanifold $X \cap Z$. Moreover, the dimension of the intersection $X \cap Z$ is given by the equation

In other words, $\dim(X\cap Z) = \dim X + \dim Z - n$. We arrive at the following corollary:

Corollary 5.6. Let $X$, $Z \subseteq {\mathbb{R}}^n_{ > 0}$ be submanifolds, where $Z$ is boundaryless. Suppose $X$ intersects $Z$ transversally and $\partial X$ also intersects $Z$ transversally. Then $X \cap Z$ is a manifold with boundary $\partial(X \cap Z) = \partial X \cap Z$ and of dimension $\dim (X \cap Z) = \dim X + \dim Z - n$.

Our main result is:

Theorem 5.7. Let $S$, $\widetilde{S} \subseteq {\mathbb{R}}^n$ be vector subspaces of equal dimension with $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$. Then for any positive vectors ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$, the intersection $({\boldsymbol{x}_0}+S) \cap ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$ consists of exactly one point.

Proof. Let ${\boldsymbol{x}_0}$, ${\boldsymbol{x}}^* \in {\mathbb{R}}_{ > 0}^n$ be arbitrary positive vectors. Lemma 5.2 implies that the intersection $({\boldsymbol{x}_0}+S) \cap ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp)$ contains at most one point. Consider first ${\boldsymbol{x}}^* \in {\boldsymbol{x}_0} +S$. Clearly, $({\boldsymbol{x}_0} +S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp) = \{ {\boldsymbol{x}}^*\}$.

Now consider the case when ${\boldsymbol{x}}^* \not\in {\boldsymbol{x}_0} +S$. Let $d = \dim S$. We define a $(d+1)$-dimensional affine strip, which we will intersect with $({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$. To define this affine strip, we consider the interpolation function

Since the line segment $K([0, 1]) \subseteq {\mathbb{R}}_{ > 0}^n$ is compact, the intersection $(K([0, 1]) + S) \cap ({\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp) \subseteq {\mathbb{R}}_{ > 0}^n$ is compact by Lemma 5.3. Moreover, the manifolds $K([0, 1])+S$ and ${\boldsymbol{x}}^*\circ \exp \widetilde{S}^\perp$ intersect transversally, as a consequence of Lemma 5.4, i.e.,

By Corollary 5.6, the intersection $\Gamma = (K([0, 1])+S)\cap ({\boldsymbol{x}}^* \circ \exp\widetilde{S}^\perp)$ is a manifold with boundary $\partial \Gamma \subseteq \partial (K([0, 1])+S) = ({\boldsymbol{x}}^* + S) \cup ({\boldsymbol{x}_0} +S)$. In addition, $\Gamma$ is 1-dimensional because

Consider the connected component $\Gamma^* \subseteq \Gamma$ containing the point ${\boldsymbol{x}}^*$. The point ${\boldsymbol{x}}^*$ must be an endpoint of $\Gamma^*$; otherwise uniqueness fails at $K(\delta_0) +S$ for some small $\delta_0 > 0$. Since $\Gamma^*$ is compact, it is a curve with two endpoints. As $\partial \Gamma^* \subseteq \partial \Gamma = ({\boldsymbol{x}}^* + S) \cup ({\boldsymbol{x}_0} + S)$, by uniqueness the other endpoint of $\Gamma^*$ must be in ${\boldsymbol{x}_0} +S$. Thus, a point exists in $({\boldsymbol{x}_0} +S) \cap ({\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp)$.

We apply Theorem 5.7 to show the existence and uniqueness of vertex-balanced steady state for a generalized mass-action system.

Theorem 5.8 (Vertex-balanced steady states of a generalized mass-action system). Let $(G, \Phi, \widetilde{\Phi})$ be a weakly reversible generalized reaction network, with stoichiometric subspace $S$ and kinetic-order subspace $\widetilde{S}$. Assume that $\dim S = \dim \widetilde{S}$ and $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$. Then the following statements hold:

i) Suppose for some rate constants ${\boldsymbol \kappa}$, the generalized mass-action system $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$ admits a vertex-balanced steady state ${\boldsymbol{x}}^*$. Then every stoichiometric compatibility class contains exactly one vertex-balanced steady state.

ii) $\widetilde{\delta}_G = 0$ if and only if the generalized mass-action system ($G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$ admits a vertex-balanced steady state ${\boldsymbol{x}}^*$ for all rate constants ${\boldsymbol \kappa}$. In this case, every stoichiometric compatibility class contains exactly one vertex-balanced steady state.

iii) Under the premises of i), additionally suppose $\delta_G = 0$. Then every stoichiometric compatibility class contains exactly one positive steady state, which is vertex-balanced.

Proof. As ${\boldsymbol{x}}^*$ is a vertex-balanced steady state for $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$, the set of vertex-balanced steady state is $Z_{\boldsymbol \kappa} = {\boldsymbol{x}}^* \circ \exp \widetilde{S}^\perp$. By Theorem 5.7, $Z_{\boldsymbol \kappa}$ intersects the stoichiometric compatibility class ${\boldsymbol{x}_0} +S$ exactly once for any ${\boldsymbol{x}_0} \in {\mathbb{R}}_{ > 0}^n$. This proves statement i).

The first part of statement ii) is the content of [6,Theorem 1(a)]. By statement i), we conclude that every stoichiometric compatibility class contains exactly one vertex-balanced steady state.

If in addition, $\delta_G = 0$, then $E_{\boldsymbol \kappa} = Z_{\boldsymbol \kappa}$, i.e., there are no positive steady states that are not vertex-balanced. Consequently, there exists a unique steady state within each stoichiometric compatibility class, which is vertex-balanced. This proves statement iii).

We state a simpler version of iii) in the theorem above.

Corollary 5.9. Let $(G, \Phi, \widetilde{\Phi})$ be a weakly reversible generalized reaction network, with stoichiometric subspace $S$ and kinetic-order subspace $\widetilde{S}$. Suppose that $\dim S = \dim \widetilde{S}$, $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$, and $\delta_G = \widetilde{\delta}_G = 0$. Then for any choice of rate constants, every stoichiometric compatibility class contains exactly one positive steady state, which is vertex-balanced.

We have focused almost exclusively on the existence and uniqueness of vertex-balanced steady states for generalized mass-action systems. For complex-balanced equilibria of classical mass-action systems, more is known. For example, complex-balanced equilibria are locally asymptotically stable within their stoichiometric compatibility classes. They are conjectured to be globally stable in their stoichiometric compatibility classes; this is known as the global attractor conjecture ^[13,27]. In particular, it has been shown that a complex-balanced equilibrium of a mass-action system is globally stable within its stoichiometric compatibility class if the network has a single connected component ^[23], or is strongly endotactic ^[25,26], or if the system is in ${\mathbb{R}}^3$ ^[25,26]. A proof of the global attractor conjecture in full generality has been proposed in ^[27].

For planar generalized mass-action systems (in particular, S-systems), local and even global stability of vertex-balanced steady states have been characterized in ^[31,32,33]. For generalized mass-action systems of arbitrary dimension, necessary conditions for linear stability have been given in ^[14]. Obviously, it is not true that a vertex-balanced steady state is always globally stable within its stoichiometric compatibility class, since it is possible for a generalized mass-action system to have multiple vertex-balanced steady states within the same stoichiometric compatibility class. Consider, for example, the following generalized mass-action system:

where each box is a vertex of the graph; the top entry in each box is the stoichiometric complex of that vertex ($0$ and $\rm{X_1 + X_2}$), and the bottom entry in the parentheses is the kinetic-order complex ($\rm{2X_1}$ and $\rm{X_1 + 2 X_2}$). The associated dynamical system of this generalized mass-action system is given by

One can check that the set of vertex-balanced steady states is $Z_{\boldsymbol \kappa} = \{ (t^2, t) : t > 0 \}$. If ${\boldsymbol{x}_0} = (0, \varepsilon)^T$ where $0 < \varepsilon < \frac{1}{4}$, then there are two vertex-balanced steady states in ${\boldsymbol{x}_0} + S = \{ (r, \varepsilon + r) : r \in {\mathbb{R}}\}$. In particular, this implies that these vertex-balanced steady states cannot be globally stable in their stoichiometric compatibility class.

Moreover, it is also possible for a unique vertex-balanced steady state (within its stoichiometric compatibility class) to be unstable. Consider the generalized mass-action system:

Its associated dynamical system is

This is an example of a reversible generalized mass-action system with $\delta_G = \widetilde{\delta}_G = 0$, and its stoichiometric subspace $S$ and its kinetic-order subspace $\widetilde{S}$ are ${\mathbb{R}}^2$. There is a unique positive steady state $\boldsymbol{x}^* = (1, 1)^T$, which is vertex-balanced; nonetheless, it can be shown that this steady state is a saddle point. Moreover, all solutions that start outside its stable manifold converge to the origin or infinity; in particular, the system is neither persistent nor permanent.

6.   An illustrative example

We conclude by applying Theorem 5.8 to the following example of a family of generalized mass-action systems. Let $a$, $b$, $\kappa_i > 0$. Consider the generalized mass-action system $(G_{\boldsymbol \kappa}, \Phi, \widetilde{\Phi})$

At each vertex (box), a stoichiometric complex (top entry) and a kinetic-order complex (second entry in parentheses) are assigned. Let $x_i$ be the concentration of species $\rm{X_i}$, for $1 \leq i \leq 4$, and ${\boldsymbol{x}} = (x_1, x_2, x_3, x_4)^T$. The stoichiometric complexes and kinetic-order complexes are

The associated dynamical system is

Another way to write the system of differential equations is

The stoichiometric subspace and the kinetic-order subspace are

respectively. Their sign vectors are

where the dots indicate the negatives of the listed sign vectors. By visual inspection, we find that $\sigma(S) \subseteq \overline{{\sigma(\widetilde{S})}}$. Moreover, one can check that the deficiency $\delta_G$ and the kinetic-order deficiency $\widetilde{\delta}_G$ are zero. Therefore, Corollary 5.9 applies and we conclude that, for any choice of rate constants, every stoichiometric compatibility class contains exactly one positive steady state, which is vertex-balanced.

Acknowledgments

This project started during the AIM SQuaRE workshop on "Dynamical properties of deterministic and stochastic models of reaction networks". The authors thank the American Institute of Mathematics for the productive environment and hospitality. GC was partially supported by NSF-DMS grants 1412643 and 1816238; SM was supported by the Austria Science Fund (FWF), project P28406; CP was partially supported by NSF-DMS grant 1517577; and PYY was partially supported by a NSERC PGS-D award and NSF-DMS grant 1412643.

Conflict of interest

The authors declare there is no conflicts of interest.