Semistability of complex balanced kinetic systems with arbitrary time delays

TL;DR

This paper introduces a class of complex balanced kinetic systems with arbitrary time delays, proving their semistability and local asymptotic stability using Lyapunov-Krasovskii functionals.

Contribution

It extends the theory of complex balanced kinetic systems to include arbitrary time delays and establishes their semistability and stability properties.

Findings

Proves uniqueness of equilibrium solutions within delayed stoichiometric classes.

Establishes semistability of equilibria using Lyapunov-Krasovskii functionals.

Shows local asymptotic stability of positive equilibria.

Abstract

In this letter we introduce a class of delayed kinetic systems derived from mass action type reaction network models. We define the time delayed positive stoichiometric compatibility classes and the notion of complex balanced time delayed kinetic systems. We prove the uniqueness of equilibrium solutions within the time delayed positive stoichiometric compatibility classes for such models. In our main result we prove the semistability of the equilibrium solutions for complex balanced systems with arbitrary time delays using an appropriate Lyapunov-Krasovskii functional and LaSalle's invariance principle. As a consequence, we obtain that every positive complex balanced equilibrium solution is locally asymptotically stable relative to its positive stoichiometric compatibility class.

Equations

$$\dot{x}(t)=\sum _{k=1}^M{\kappa }_k(x(t){)}^{y_k}\left[y_k^{\prime }-y_k\right],t\ge 0,$$

$$\dot{x}(t)=\sum _{k=1}^M{\kappa }_k\left[(x(t-{\tau }_k){)}^{y_k}{y}_k^{\prime }-(x(t){)}^{y_k}{y}_k\right],t\ge 0,$$

$$\dot{x}(t)=F(x_t),$$

$$F(\varphi )=\sum _{k=1}^M{\kappa }_k\left[(\varphi (-{\tau }_k){)}^{y_k}{y}_k^{\prime }-(\varphi (0){)}^{y_k}{y}_k\right],\varphi \in {\overline{\mathcal{C}}}_{+}.$$

$$F_i(\varphi )=\sum _{k=1}^M{\kappa }_k(\varphi (-{\tau }_k){)}^{y_k}{(y_k^{\prime })}_i\ge 0.$$

$$\mathcal{S}=\text{span}\left\{y_k^{\prime }-y_k\mid k=1,\dots ,M\right\},$$

$${\mathcal{S}}_p=\left\{x\in {\overline{\mathbb{R}}}_{+}^N\mid x-p\in \mathcal{S}\right\}.$$

$$c_v(\psi )=v^T\left[\psi (0)+\sum _{k=1}^M\left({\kappa }_k{\int }_{-{\tau }_k}^0(\psi (s){)}^{y_k}ds\right)y_k\right],\psi \in {\overline{\mathcal{C}}}_{+}.$$

$${\mathcal{D}}_{\theta }=\{\psi \in {\overline{\mathcal{C}}}_{+}\mid c_v(\psi )=c_v(\theta )\text{for all v∈S⊥}\}.$$

$$\psi (0)-\theta (0)+\sum _{k=1}^M\left({\kappa }_k{\int }_{-{\tau }_k}^0[(\psi (s){)}^{y_k}-(\theta (s){)}^{y_k}]ds\right)y_k\in \mathcal{S}.$$

$$\begin{array}{rl}\frac{d}{dt}(c_v(x_t))& =v^T\sum _{k=1}^M{\kappa }_k(x(t-{\tau }_k){)}^{y_k}\left(y_k^{\prime }-y_k\right)\\ & =\sum _{k=1}^M{\kappa }_k(x(t-{\tau }_k){)}^{y_k}{v}^T(y_k^{\prime }-y_k)=0,\end{array}$$

$$c_v(x_t^{\psi })=c_v(x_0^{\psi })=c_v(\psi )=c_v(\theta )$$

$$\sum _{k=1}^M{\kappa }_k(\overline{x}{)}^{y_k}\left[y_k^{\prime }-y_k\right]=0.$$

$$\sum _{k:\eta =y_k}{\kappa }_k(\overline{x}{)}^{y_k}=\sum _{k:\eta =y_k^{\prime }}{\kappa }_k(\overline{x}{)}^{y_k},$$

$$\mathcal{E}=\{\tilde{x}\in {\mathbb{R}}_{+}^N\mid \text{Ln}(\tilde{x})-\text{Ln}(\overline{x})\in {\mathcal{S}}^{\perp }\}.$$

$$\tilde{x}-\overline{x}+\sum _{k=1}^M\left({\kappa }_k{\int }_{-{\tau }_k}^0[{\tilde{x}}^{y_k}-{\overline{x}}^{y_k}]ds\right)y_k\in \mathcal{S}.$$

$$\begin{array}{rl}0& =(\text{Ln}(\tilde{x})-\text{Ln}(\overline{x}){)}^T\left[\tilde{x}-\overline{x}+\sum _{k=1}^M\left({\kappa }_k{\int }_{-{\tau }_k}^0[{\tilde{x}}^{y_k}-{\overline{x}}^{y_k}]ds\right)y_k\right]\\ & =\sum _{i=1}^N\left(\ln ({\tilde{x}}_i)-\ln ({\overline{x}}_i)\right)({\tilde{x}}_i-{\overline{x}}_i)+\sum _{k=1}^M{\kappa }_k{\tau }_k\left(\ln ({\tilde{x}}^{y_k})-\ln ({\overline{x}}^{y_k})\right)\left({\tilde{x}}^{y_k}-{\overline{x}}^{y_k}\right).\end{array}$$

$$e^a(b-a)\le e^b-e^a,$$

$$x\bigl{[}\ln(x)-\ln(b)-1\bigr{]}+b\geq c\ln\bigl{[}\,1+(x-b)^{2}\,\bigr{]}\geq 0.$$

$$\begin{split}V(\psi)=&\sum_{i=1}^{N}\bigl{(}\psi_{i}(0)(\ln(\psi_{i}(0))-\ln(\overline{x}_{i})-1)+\overline{x}_{i}\bigr{)}\\ &+\sum_{k=1}^{M}\kappa_{k}\int_{-\tau_{k}}^{0}\left\{(\psi(s))^{y_{k}}\left[\ln((\psi(s))^{y_{k}})-\ln(\overline{x}^{y_{k}})-1\right]+\overline{x}^{y_{k}}\right\}ds\end{split}$$

$$V(\psi )\ge \alpha (|\psi (0)-\overline{x}|),\psi \in {\mathcal{C}}_{+},$$

$$\begin{split}V(\psi)\geq&\sum_{i=1}^{N}c_{i}\ln\left[1+(\psi_{i}(0)-\overline{x}_{i})^{2}\right]\geq\gamma\sum_{i=1}^{N}\ln\left[1+(\psi_{i}(0)-\overline{x}_{i})^{2}\right]\\ &=\gamma\ln\prod_{i=1}^{N}\left[1+(\psi_{i}(0)-\overline{x}_{i})^{2}\right]\geq\gamma\ln\biggl{(}1+\sum_{i=1}^{N}(\psi_{i}(0)-\overline{x}_{i})^{2}\biggr{)}\\ &=\gamma\ln\left(1+|\psi(0)-\overline{x}|^{2}\right),\end{split}$$

$$\alpha (r)=\gamma \ln (1+r^2),r\ge 0.$$

$$\begin{array}{rl}\dot{V}(x_t)& =\sum _{k=1}^M{\kappa }_k\text{Ln}{\left(\frac{x(t)}{\overline{x}}\right)}^T\left[(x(t-{\tau }_k){)}^{y_k}{y}_k^{\prime }-(x(t){)}^{y_k}{y}_k\right]\\ & +\sum _{k=1}^M{\kappa }_k(x(t){)}^{y_k}\left(\ln \left({\left\{\frac{x(t)}{\overline{x}}\right\}}^{y_k}\right)-1\right)\\ & -\sum _{k=1}^M{\kappa }_k(x(t-{\tau }_k){)}^{y_k}\left(\ln \left({\left\{\frac{x(t-{\tau }_k)}{\overline{x}}\right\}}^{y_k}\right)-1\right)\\ & =\sum _{k=1}^M{\kappa }_k\left[(x(t-{\tau }_k){)}^{y_k}\ln \left({\left\{\frac{x(t)}{\overline{x}}\right\}}^{y_k^{\prime }}\right)-(x(t){)}^{y_k}\ln \left({\left\{\frac{x(t)}{\overline{x}}\right\}}^{y_k}\right)\right]\\ & +\sum _{k=1}^M{\kappa }_k\left[(x(t){)}^{y_k}\ln \left({\left\{\frac{x(t)}{\overline{x}}\right\}}^{y_k}\right)-(x(t-{\tau }_k){)}^{y_k}\ln \left({\left\{\frac{x(t-{\tau }_k)}{\overline{x}}\right\}}^{y_k}\right)\right]\\ & +\sum _{k=1}^M{\kappa }_k\left[(x(t-{\tau }_k){)}^{y_k}-(x(t){)}^{y_k}\right]\\ & =\sum _{k=1}^M{\kappa }_k{\overline{x}}^{y_k}{\left(\frac{x(t-{\tau }_k)}{\overline{x}}\right)}^{y_k}\left[\ln \left({\left\{\frac{x(t)}{\overline{x}}\right\}}^{y_k^{\prime }}\right)-\ln \left({\left\{\frac{x(t-{\tau }_k)}{\overline{x}}\right\}}^{y_k}\right)\right]\\ & +\sum _{k=1}^M{\kappa }_k{\overline{x}}^{y_k}\left[{\left(\frac{x(t-{\tau }_k)}{\overline{x}}\right)}^{y_k}-{\left(\frac{x(t)}{\overline{x}}\right)}^{y_k}\right].\end{array}$$

$$\begin{array}{rl}& {\left(\frac{x(t-{\tau }_k)}{\overline{x}}\right)}^{y_k}\left[\ln \left({\left\{\frac{x(t)}{\overline{x}}\right\}}^{y_k^{\prime }}\right)-\ln \left({\left\{\frac{x(t-{\tau }_k)}{\overline{x}}\right\}}^{y_k}\right)\right]\\ & \le {\left(\frac{x(t)}{\overline{x}}\right)}^{y_k^{\prime }}-{\left(\frac{x(t-{\tau }_k)}{\overline{x}}\right)}^{y_k}\end{array}$$

$${\left(\frac{x(t)}{\overline{x}}\right)}^{y_k^{\prime }}={\left(\frac{x(t-{\tau }_k)}{\overline{x}}\right)}^{y_k}.$$

$$\begin{array}{rl}\dot{V}(x_t)& \le \sum _{k=1}^M{\kappa }_k{\overline{x}}^{y_k}\left[{\left(\frac{x(t)}{\overline{x}}\right)}^{y_k^{\prime }}-{\left(\frac{x(t)}{\overline{x}}\right)}^{y_k}\right]\\ & =\sum _{\eta \in \mathcal{K}}{\left(\frac{x(t)}{\overline{x}}\right)}^{\eta }\left[\sum _{k:\eta =y_k^{\prime }}{\kappa }_k{\overline{x}}^{y_k}-\sum _{k:\eta =y_k}{\kappa }_k{\overline{x}}^{y_k}\right]=0,\end{array}$$

$$\mathcal{R}=\{\psi \in {\hat{\mathcal{B}}}_{ϵ}(\overline{x})\mid \dot{V}(\psi )=0\}.$$

$$\mathcal{R}=\biggl{\{}\,\psi\in\hat{\mathcal{B}}_{\epsilon}(\overline{x})\mid\left(\frac{\psi(0)}{\overline{x}}\right)^{y_{k}^{\prime}}=\left(\frac{\psi(-\tau_{k})}{\overline{x}}\right)^{y_{k}}\quad\text{for $k=1,\dots,M$}\biggr{\}}.$$

$$\dot{x}(t)=\sum _{\eta \in \mathcal{K}}\left[\sum _{k:\eta =y_k^{\prime }}{\kappa }_k{\overline{x}}^{y_k}{\left(\frac{x(t-{\tau }_k)}{\overline{x}}\right)}^{y_k}-\sum _{k:\eta =y_k}{\kappa }_k{\overline{x}}^{y_k}{\left(\frac{x(t)}{\overline{x}}\right)}^{y_k}\right]\eta .$$

Full text

Keywords: Nonnegative systems; Kinetic systems; Chemical reaction networks; Stability Theory; Time delay; Logarithmic Lyapunov–Krasovskii functionals

1 Introduction

The class of kinetic systems has proven to be a useful representation of nonnegative system models not only in biochemistry, but also in other areas like population or disease dynamics, process systems, and even transportation networks [1, 2, 3]. A network-based description is often advantageous to describe key properties of potentially large, complex systems with many components [4, 5]. Kinetic systems are naturally equipped with a network (i.e., directed graph) structure called the reaction graph, which is the abstraction of a set of chemical reactions, where the chemical complexes and reactions can be represented by vertices and directed edges, respectively. One of the primary aims of chemical reaction network theory (CRNT) is to discover relations between the dynamical behaviour and the graph structure of kinetic systems [6, 7, 8]. Probably the most widely known results of general importance in this field are the Deficiency Zero and Deficiency One Theorems [9] and more recently the notion of absolute concentration robustness [10].

The detailed balance property of a thermodynamic system, defined originally by Boltzmann in the 19th century, means that at equilibrium, each elementary reaction step is equilibriated by the corresponding reverse reaction. A more general condition is complex balance, which requires that the signed sum of incoming and outgoing reaction rates at equilibrium is zero for each complex in a chemical reaction network [11, 12, 13]. It is worth remarking that complex balance does not depend on a particular equilibrium (if there exist multiple equilibria in a system), but it is a property of a chemical reaction network itself [6]. For a historical review of the notions of detailed and complex balance, see [14]. Generally, complex balance is related both to the structure and to the parameters of chemical reaction networks. Firstly, complex balance implies that each component of the reaction graph is strongly connected (i.e., the reaction network is weakly reversible) [6]. It is also important that deficiency zero weakly reversible reaction networks are complex balanced for any positive values of the reaction rate coefficients [7, 9]. However, complex balance becomes a parameter-dependent property when the deficiency of the network is higher than zero. The main significance of complex balance in systems and control theory stands in its stability implications [15]. According to the Global Attractor Conjecture, complex balanced kinetic systems are globally stable in the positive orthant with a logarithmic Lyapunov function that does not depend on the model parameters. The conjecture was proved for several special cases such as one linkage class networks [16], and a possible proof for the general problem has recently appeared in [17]. Using the non-uniqueness of reaction graphs corresponding to kinetic models [18], a feedback design method was proposed in [19] that transforms a polynomial model into a complex balanced closed loop system via nonlinear state feedback.

Time-delays are often present in natural and technological processes, and the detailed mathematical treatment of such delays is sometimes necessary to model and understand important observed dynamical phenomena [20, 21]. An excellent summary of the fundamental results on nonnegative and compartmental systems with time-delay can be found in Chapter 3 of [2], where simple algebraic necessary and sufficient conditions are given for the asymptotic stability of delayed linear nonnegative models including linear compartmental systems. Among other results, the semistability of an important special class of nonlinear compartmental systems for arbitrary time-delays was shown in [22].

Motivated by the above results, the purpose of this paper is to introduce the complex balance condition for kinetic systems with delayed reactions, and to study the stability properties of such systems using logarithmic Lyapunov-Krasovskii functionals and LaSalle’s invariance principle.

Throughout the paper, we will use the following notations. If $N$ is a positive integer, ${\mathbb{R}}^{N}$ denotes the $N$-dimensional space of real column vectors. The symbols $\mathbb{R}^{N}_{+}$ and $\overline{\mathbb{R}}^{N}_{+}$ denote the set of (element-wise) positive and nonnegative vectors in $\mathbb{R}^{N}$, respectively. For $x,y\in\mathbb{R}^{N}_{+}$, the vector $\frac{x}{y}\in{\mathbb{R}}^{N}_{+}$ is defined by ${(\frac{x}{y})}_{i}=\frac{x_{i}}{y_{i}}$ for $i=1,\dots,N$. For $x,y\in\overline{\mathbb{R}}^{N}_{+}$, the vector exponential $x^{y}$ is defined as $x^{y}=\prod_{i=1}^{N}x_{i}^{y_{i}}$. The mapping $\text{Ln}:~{}\mathbb{R}^{N}_{+}\rightarrow\mathbb{R}^{N}$ is the element-wise logarithmic mapping defined by ${(\text{Ln}(x)})_{i}=\ln(x_{i})$ for $x\in\mathbb{R}^{N}_{+}$ and $i=1,\dots,N$. Recall that $\ln(x^{y})=y^{T}\text{Ln}(x)$ for $x\in\mathbb{R}^{N}_{+}$ and $y\in\mathbb{R}^{N}$, where $(\cdot)^{T}$ denotes the transpose, and $\text{Ln}(\frac{x}{y})=\text{Ln}(x)-\text{Ln}(y)$ whenever $x,y\in\mathbb{R}^{N}_{+}$. For every $\tau\geq 0$, the symbol $\mathcal{C}=C([-\tau,0],\mathbb{R}^{N})$ denotes the Banach space of continuous functions mapping the interval $\left[-\tau,0\right]$ into $\mathbb{R}^{N}$ with the norm $\|\psi\|=\sup_{-\tau\leq s\leq 0}|\psi(s)|$ for $\psi\in\mathcal{C}$, where $|\cdot|$ denotes the Eucledian norm in $\mathbb{R}^{N}$. Finally, let $\mathcal{C}_{+}=C([-\tau,0],\mathbb{R}^{N}_{+})$ and $\overline{\mathcal{C}}_{+}=C([-\tau,0],\overline{\mathbb{R}}^{N}_{+})$ denote the set of positive and nonnegative functions in $\mathcal{C}$.

2 Kinetic systems with time delays

In this section, we introduce mass-action kinetic systems with time delays and show that they generate a nonnegative semiflow.

Consider the ordinary mass-action kinetic system [7]

[TABLE]

where $x(t)\in\overline{\mathbb{R}}^{N}_{+}$ is the state vector. We have a set of complexes $\mathcal{K}\subset\overline{\mathbb{Z}}^{N}_{+}$ and there are $M$ reactions between the complexes. As usual, $\overline{\mathbb{Z}}^{N}_{+}$ denotes the set on nonnegative integers. Each reaction has a source and product complex $y_{k},y_{k}^{\prime}\in\mathcal{K}$, respectively, with a reaction rate constant $\kappa_{k}>0$, $k=1,\dots,M$. Solutions of (1) are determined by nonnegative initial vectors $x(0)=\eta\in\overline{\mathbb{R}}^{N}_{+}$. In this paper, we will consider the mass-action kinetic system with time delays

[TABLE]

where $\tau_{k}\geq 0$, $k=1,\dots,M$. In the special case $\tau_{k}=0$, $k=1,\dots,M$, Eq. (2) reduces to the ordinary mass kinetic system (1). Solutions of (2) are generated by initial data $x(t)=\theta(t)$ for $-\tau\leq t\leq 0$, where $\tau=\max_{1\leq k\leq M}\tau_{k}$ is the maximum delay and $\theta\in\overline{\mathcal{C}}_{+}$ is a nonnegative continuous initial function. Throughout the paper, the solution of (2) with initial function $\theta\in\overline{\mathcal{C}}_{+}$ will be denoted by $x=x^{\theta}$. Note that the solutions of delay differential equations are usually interpreted in $\mathcal{C}$. For every $t\geq 0$, $x_{t}\in\mathcal{C}$ is defined by $x_{t}(s)=x(t+s)$ for $-\tau\leq s\leq 0$.

In the following theorem, we show that the semiflow generated by the time delay kinetic system (2) is nonnegative.

Theorem 1**.**

For every initial function $\theta\in\overline{\mathcal{C}}_{+}$, the solution $x^{\theta}$ of (2) is nonnegative, i.e., $x^{\theta}_{t}\in\overline{\mathcal{C}}_{+}$ for all $t\geq 0$.

Proof.

Eq. (2) can be written in the form

[TABLE]

where $F:\overline{\mathcal{C}}_{+}\rightarrow\mathbb{R}^{N}$ is given by

[TABLE]

It follows from the definition of the vector exponential that if $\phi\in\overline{\mathcal{C}}_{+}$ and $\phi_{i}(0)=0$ for some $i\in\{1,\dots,N\}$, then

[TABLE]

Here $\phi_{i}$ and $F_{i}$ denote the $i$-th coordinate function of $\phi$ and $F$, respectively. The conclusion follows from Theorem 2.1 in Chap. 5 of [23]. Alternatively, we can use the generalization of Proposition 3.1 of [2] to equations with multiple delays. ∎

3 Stoichiometric compatibility classes for delayed kinetic systems

Recall [7] that the stoichiometric subspace $\mathcal{S}$ for the ordinary mass-action kinetic system (1) is defined by

[TABLE]

and for each $p\in\overline{\mathbb{R}}^{N}_{+}$ the corresponding positive stoichiometric compatibility class $\mathcal{S}_{p}$ is given by

[TABLE]

It is well know that the positive stoichiometric classes $\mathcal{S}_{p}$ are positively invariant under the mass-action kinetic system (1), i.e. $x(0)\in\mathcal{S}_{p}$ implies $x(t)\in\mathcal{S}_{p}$ for all $t\geq 0$.

In this section we will extend the definition of the positive stoichiometric classes to the time delayed kinetic system (2) and we prove their invariance property.

For each $v\in{\mathbb{R}}^{N}$, define the functional $c_{v}:\overline{\mathcal{C}}_{+}\rightarrow\mathbb{R}$ by

[TABLE]

Let $\mathcal{S}^{\perp}$ denote the orthogonal complement of the stoichiometric subspace $\mathcal{S}$ given by \mathcal{S}^{\perp}=\{\,v\in\mathbb{R}^{N}\mid v^{T}y=0\text{ for all y\in\mathcal{S}}\,\}. Now we can formulate the definition of the positive stoichiometric compatibility classes for the delayed kinetic system (2). For each $\theta\in\overline{\mathcal{C}}_{+}$, the positive stoichiometric compatibility class of (2) corresponding to $\theta$ is denoted by $\mathcal{D}_{\theta}$ and is defined by

[TABLE]

It is easily seen that $\psi\in\mathcal{D}_{\theta}$ if and only if $\psi\in\overline{\mathcal{C}}_{+}$ and

[TABLE]

Therefore if we ignore the delays in (2), i.e. $\tau_{k}=0$ for $k=1,\dots,M$, then the above delayed positive stoichiometric compatibility classes coincide with the positive stoichiometric compatibility classes of the ordinary kinetic system  (1).

In the next theorem, we establish the invariance property of the above delayed positive stoichiometric compatibility classes.

Theorem 2**.**

For every $\theta\in\overline{\mathcal{C}}_{+}$, the positive stoichiometric compatibility class $\mathcal{D}_{\theta}$ is a closed subset of $\overline{\mathcal{C}}_{+}$. Moreover, $\mathcal{D}_{\theta}$ is positively invariant under Eq. (2), i.e. if $\psi\in\mathcal{D}_{\theta}$, then $x^{\psi}_{t}\in\mathcal{D}_{\theta}$ for all $t\geq 0$.

Proof.

Let $\theta\in\overline{\mathcal{C}}_{+}$. The closedness of $\mathcal{D}_{\theta}$ is a simple consequence of the continuity of functionals $c_{v}$, $v\in\mathcal{S}^{\perp}$. We will show that for every $v\in\mathcal{S}^{\perp}$ the functional $c_{v}$ defined by (5) is constant along the solutions of Eq. (2). Indeed, if $x$ is a solution of (2), then we have for $t\geq 0$,

[TABLE]

the last equality being a consequence of the definition of $\mathcal{S}^{\perp}$. From this, we find that if $\psi\in\mathcal{D}_{\theta}$, then for every $v\in\mathcal{S}^{\perp}$ and $t\geq 0$,

[TABLE]

and hence $x^{\psi}_{t}\in\mathcal{D}_{\theta}$. This show that $\mathcal{D}_{\theta}$ invariant under Eq. (2). ∎

4 Semistability for delayed complex balanced kinetic systems

Before we formulate our main stability criterion, we recall some definitions.

By a positive equilibrium of (1) or (2), we mean a positive vector $\overline{x}\in\mathbb{R}^{N}_{+}$ such that $x(t)\equiv\overline{x}$ is a solution of (1) and (2), respectively. Note that Eqs. (1) and (2) share the same equilibria satisfying the algebraic equation

[TABLE]

A positive equilibrium $\overline{x}$ is called complex balanced if for every $\eta\in\mathcal{K}$,

[TABLE]

where the sum on the left is over all reactions for which $\eta$ is the source complex and the sum on the right is over all reactions for which $\eta$ is the product complex. Finally, an ordinary or delayed kinetic system is called complex balanced if it has a positive complexed balanced equilibrium.

It is well-known [24] that if Eq. (1) and hence (2) has a positive complex balanced equilibrium $\overline{x}$, then any other positive equilibrium is complex balanced and the set of all positive equilibria $\mathcal{E}$ can be characterized by

[TABLE]

Now we formulate the main result of the paper about the semistability of positive equilibria of delayed complex balanced systems in the sense of the following definition. A positive equilibrium $\overline{x}$ of Eq. (2) is called semistable if it is Lyapunov stable and there exists $\delta>0$ such that if $\theta\in\hat{\mathcal{B}}_{\delta}(\overline{x})$, then $x^{\theta}(t)$ converges to a Lyapunov stable equilibrium of (2) as $t\to\infty$. As usual, $\hat{\mathcal{B}}_{\delta}(\overline{x})=\{\,\psi\in\mathcal{C}\mid\|\psi-\overline{x}\|\leq\delta\,\}.$

Theorem 3**.**

Every positive complex balanced equilibrium of the delayed kinetic system (2) is semistable.

As a preparation for the proof of Theorem 3, we establish an auxiliary result about the uniqueness of positive equilibria in the positive stoichiometric classes of complex balanced systems.

Theorem 4**.**

Suppose that the delayed kinetic system (2) is complex balanced. Then for every $\theta\in\overline{\mathcal{C}}_{+}$ the corresponding delayed stoichiometric class $\mathcal{D}_{\theta}$ contains at most one positive equilibrium.

Proof.

Let $\theta\in\overline{\mathcal{C}}_{+}$. Suppose that $\tilde{x}$ and $\overline{x}$ are positive equilibria belonging to $\mathcal{D}_{\theta}$. From the characterization (7) of  $\mathcal{D}_{\theta}$, we find that

[TABLE]

This, together with (10), yields

[TABLE]

Since $(\ln(a)-\ln(b))(a-b)\geq 0$ whenever $a$, $b>0$ with equality if and only if $a=b$, this is possible only if $\tilde{x}_{i}=\overline{x}_{i}$ for all $i=1,\dots,N$. ∎

Now we are in a position to give a proof of Theorem 3. It will be based on the Lyapunov-Krasovskii method and LaSalle’s invariance principle [2], [25], [26], [27].

Proof of Theorem 3.

We will use the following two inequalities. For every $a,b\in\mathbb{R}$,

[TABLE]

with equality if and only if $a=b$. For every $b>0$ there exists $c>0$ such that for all $x>0$,

[TABLE]

Inequality (11) is not new. It is a simple consequence of the mean value theorem applied to the exponential function. Inequality (12) is less obvious. Its proof is given in the Appendix.

Consider the candidate Lyapunov–Krasovskii functional $V:\mathcal{C}_{+}\rightarrow\overline{\mathbb{R}}_{+}$ defined by

[TABLE]

for $\psi\in\mathcal{C}_{+}$. Clearly, $V(\overline{x})=0$. We will show that there exists a continuous strictly increasing function $\alpha:[0,\infty)\rightarrow[0,\infty)$ with $\alpha(0)=0$ such that

[TABLE]

where $|\cdot|$ is the Euclidean norm in $\mathbb{R}^{N}$. By virtue of (12), the second sum in (13) is nonnegative and the first sum in (13) can be estimated from below using the first inequality in (12). Thus, (12) implies the existence of positive constants  $c_{i}$, $1\leq i\leq N$, such that for all $\psi\in\mathcal{C}_{+}$,

[TABLE]

where $\gamma=\min_{1\leq i\leq N}c_{i}$. Thus, (14) holds with

[TABLE]

Next, it follows that the Lyapunov-Krasovskii directional derivative along trajectories of (2) is given by

[TABLE]

By virtue of (11), we have for each $k=1,\dots,M$,

[TABLE]

with equality if and only if for each $k=1,\dots,M$,

[TABLE]

From this, we find that

[TABLE]

where the last equality follows from the complex balanced property (9). This implies that the complex balanced equilibrium $\overline{x}$ of (2) is Lyapunov stable.

Choose $\epsilon$ such that $0<\epsilon<\min_{1\leq i\leq N}\overline{x}_{i}$ so that $\hat{\mathcal{B}}_{\epsilon}(\overline{x})\subset\mathcal{C}_{+}$. The Lyapunov stability of the equilibrium $\overline{x}$ implies the existence of $\delta>0$ such that if $\theta\in\hat{\mathcal{B}}_{\delta}(\overline{x})$, then $x^{\theta}_{t}\in\hat{\mathcal{B}}_{\epsilon}(\overline{x})$ for all $t\geq 0$. We will show that for every $\theta\in\hat{\mathcal{B}}_{\delta}(\overline{x})$ the solution $x^{\theta}(t)$ converges to a Lyapunov stable equilibrium of (2). Let

[TABLE]

From the previous calculations, we find that

[TABLE]

Let $\mathcal{M}$ be the largest set in $\mathcal{R}$ which is invariant under Eq. (2). We will show that every element of $\mathcal{M}$ is a positive equilibrium of  (2). Let $\psi\in\mathcal{M}$ and write $x=x^{\psi}$ for brevity. Rewrite Eq. (2) in the form

[TABLE]

Since $\mathcal{M}\subset\mathcal{R}$ is invariant, we have that $x_{t}\in\mathcal{R}$ for all $t\geq 0$ and hence

[TABLE]

for $t\geq 0$, where the last equality is a consequence of the complex balanced property (9). Thus, $x=x^{\psi}$ is a constant solution of (9) and hence $\psi\equiv\tilde{x}$ is a positive equilibrium. Now suppose that $\theta\in\hat{\mathcal{B}}_{\delta}(\overline{x})$. As noted before, $x^{\theta}_{t}\in\hat{\mathcal{B}}_{\epsilon}(\overline{x})$ for all $t\geq 0$. By the application of LaSalle’s invariance principle [26], we conclude that $\omega(\theta)\subset\mathcal{M}$, where \omega(\theta)=\{\,\phi\in\mathcal{C}\mid\text{there exists t_{n}\rightarrow\infty$such that$x^{\theta}{t{n}}\rightarrow\phi}\,\} is the omega limit set. On the other hand, since $\theta\in\mathcal{D}_{\theta}$ and according to Theorem 2 the stoichiometric class $\mathcal{D}_{\theta}$ is closed and invariant, it follows that $\omega(\theta)\subset\mathcal{D}_{\theta}$. Thus, $\omega(\theta)\subset\mathcal{M}\cap\mathcal{D}_{\theta}$. As shown before, every element of $\mathcal{M}$ is a positive equilibrium of (2), while Theorem 4 implies that $\mathcal{D}_{\theta}$ contains at most one positive equilibrium. Hence $\omega(\theta)=\{\tilde{x}\}$ for some $\tilde{x}\in\mathcal{E}$ and $x^{\theta}(t)\rightarrow\tilde{x}$ as $t\to\infty$. The Lyapunov stability of the positive equilibrium $\tilde{x}$ follows from the first part of the proof. ∎

Remark**.**

In the previous proof we have shown that for every positive initial function $\theta$ from a neighborhood of the positive complex balanced equilibrium $\overline{x}$ of (2) the stoichiometric class $\mathcal{D}_{\theta}$ contains exactly one positive equilibrium. A simple modification of the above proof can be used to show that if system (2) is comlexed balanced then $\mathcal{D}_{\theta}$ contains exactly one positive equilibrium whenever the closure of the forward orbit $\mathcal{O}^{+}_{\theta}=\{\,x^{\theta}_{t}\mid t\geq 0\,\}$ remains in $\mathcal{C}_{+}$. This last condition certainly holds if the solution $x^{\theta}$ is persistent in the sense that $\liminf_{t\to\infty}x^{\theta}_{i}(t)>0$ for each $i=1,\dots,N$.**

Let $\overline{x}$ be a positive complex balanced equilibrium of Eq. (2). Theorem 4 implies that $\overline{x}$ is the only positive equilibrium in its positive stoichiometric compatibility class $\mathcal{D}_{\overline{x}}$. This, together with Theorem 3 yields the following analogue of a known result for ordinary kinetic systems.

Theorem 5**.**

Every positive complex balanced equilibrium $\overline{x}$ of the delayed kinetic system (2) is locally asymptotically stable relative to its positive stoichiometric compatibility class $\mathcal{D}_{\overline{x}}$.

5 Example

In this section, we will illustrate our results and notations on a simple example. The studied system is intentionally low dimensional in order to be able to simply illustrate the relations and differences between non-delayed and delayed kinetic systems.

Let the time delayed complex balanced kinetic system be given by a reversible reaction $2X_{1}\leftrightarrows X_{2}$ containing one undelayed and a delayed reaction as follows

[TABLE]

Then, the corresponding time-delay differential equation is

[TABLE]

where $x=[x_{1},x_{2}]^{T}\in\overline{\mathbb{R}}^{2}_{+}$ are the states and $\tau_{2}$ is the time delay of the second reaction. It is easily verified that $[2,2]^{T}$ is a positive complex balanced equilibrium of  (15). The stoichiometric subspace is

[TABLE]

The dimension of $S$ is one, therefore Eq. (15) has infinitely many positive equilibria given by the set

[TABLE]

For $\overline{x}\in\mathcal{E}$, consider the set $\mathcal{X}_{\overline{x}}$ of those positive constant functions which belong to $\mathcal{D}_{\overline{x}}$:

[TABLE]

According to Theorem 5, if $\theta\equiv\eta\in\mathcal{X}_{\overline{x}}$ is close to the equilibrium $\overline{x}\in\mathcal{E}$, then $x^{\theta}(t)\rightarrow\overline{x}$ as $t\to\infty$.

Figure 1 shows the phase portrait of the system (15) with $\tau_{2}=0.5$ and with different constant initial conditions. The initial conditions are chosen such that the corresponding solutions converge to three different equilibria. Figure 2 shows the time domain behavior of the system (15) when there are different time delays, but the initial conditions are same.

6 Conclusions

In this paper a class of delayed kinetic systems is introduced, where different constant time-delays can be assigned to the individual reactions of the network. The complex balance property is defined for this system in a straightforward way. It is shown that the equilibrium solutions of complex balanced kinetic systems can be directly obtained from the equilibria of the corresponding non-delayed kinetic system. Therefore, the complex balance property of a delayed network can be checked in the same way as in the non-delayed model. The notion of stoichiometric compatibility classes is extended to delayed networks. It is shown that contrary to the classical mass action case, these classes are no longer linear manifolds in the state space. The uniqueness of equilibrium solutions within a time delayed positive stoichiometric compatibility class is proved for delayed complex balanced models. By introducing a logarithmic Lyapunov-Krasovskii functional and using LaSalle’s invariance principle, the semistability of equilibrium solutions in complex balanced systems with arbitrary time delays is also proved. As a consequence, a positive complex balanced equilibrium is always locally asymptotically stable relative to its positive stoichiometric class. The obtained results further underline the significance of the complex balance principle in the theory of dynamical systems. In the light of [16] and [17], an interesting question is whether the asymptotic stability of a delayed complex balanced system is global relative to its positive stoichiometric class.

Acknowledgements

We acknowledge the fruitful discussions with Professor Ferenc Hartung and Dr. Tamás Péni. This research has been supported by the National Research, Development and Innovation Office through grants K115694 and K120186. G. Sz. acknowledges the support of the grant PPKE KAP-1.1-16-ITK.

Appendix

We give a proof of inequality (12). Let $b>0$ be fixed. For $x>0$, define

[TABLE]

and

[TABLE]

Since $f^{\prime}(x)=\ln(x)-\ln(b)$ for $x>0$, $f^{\prime}<0$ on $(0,b)$ and $f^{\prime}>0$ on $(b,\infty)$. This implies that $f$ has a strict minimum at $x=b$. Hence $f(x)>f(b)=0$ for $x\in(0,b)\cup(b,\infty)$. Clearly, the same inequality holds for $g$. A repeated application of l’Hospital’s rule yields

[TABLE]

Therefore the function $h:(0,\infty)\rightarrow\mathbb{R}$ defined by

[TABLE]

is positive and continuous. Since $x\ln(x)\rightarrow 0$ as $x\to 0+$, $h$ can be extended continuously to the interval $[0,\infty)$ by

[TABLE]

Since $\lim_{x\to\infty}h(x)=\infty$, there esists $T>0$ such that $h(x)>h(0)$ for all $x>T$. The continuity of $h$ implies the existence of $c=\min_{0\leq x\leq T}h(x)$. Since $h(x)>h(0)\geq c>0$ for $x>T$, we have that $h(x)\geq c$ for all $x\geq 0$ which implies the desired inequality (12).