Chemical reaction network theory

Chemical reaction network theory is a framework for modeling the evolution of chemical concentrations resulting from simultaneously occurring chemical reactions. A key feature of the theory is the relationship between the graphical structure of the reaction network and the resulting dynamics. A strong emphasis, consequently, is placed on results which hold regardless of the parameter values of the network, i.e. results which depend on the network structure alone.

The foundations of chemical reaction network theory were laid down in a series of seminal papers by Fritz Horn, Roy Jackson and Martin Feinberg in the early 1970's . In these papers, the authors were primarily focused on developing conditions sufficient for uniqueness and stability of equilibrium concentrations, but their foundation has since between adapted to questions of multistability, injectivity, monotonicity, persistence, equivalence of mass-action systems, model reduction, oscillations, and applications. The models have also been adapted to the stochastic framework, reaction-diffusion equations, and kinetic schemes other than mass-action (e.g. Michaelis-Menten, Hill kinetics, etc.).

Chemical reaction networks

A chemical reaction network is given by sets of chemical species (or reactants) which react to form different sets of species. For example, consider the set of reactions corresponding to the Michaelis-Menten model for enzyme-facilitated transformation of a single substrate into a single product:

${\begin{array}{rcl}S+E&\rightarrow &SE\\SE&\rightarrow &S+E\\SE&\rightarrow &P+E.\end{array}}$

In this model S denotes the generic substrate, E denotes the enzyme facilitating the transfer, SE denotes the intermediate substrate-enzyme quantity, and P denotes the end product of the process. We define the species set ${\mathcal {S}}=\left\{S,E,C,P\right\}$ of the network to be the collection of chemical species capable of undergoing basic chemical change. Another fundamental feature are the combined terms which appear on the left and right of each individual reaction (the net input and net output of a reaction). These terms are called complexes. If we count repeated complexes only once, we can define the complex set of the network as ${\mathcal {C}}=\left\{S+E,SE,P+E\right\}$ (since SE and S + E are repeated). Lastly, we define the reaction set to be the set of ordered pairs between elements of the complex set given by ${\mathcal {R}}=\left\{(S+E,SE),(SE,S+E),(SE,P+E)\right\}$. The chemical reaction network is given by the triplet ${\mathcal {N}}=({\mathcal {S}},{\mathcal {C}},{\mathcal {R}})$.

This formulation allows us to restructure the list of elementary reactions given above with each stoichiometrically distinct complex appearing only once. The resulting condensed graph is called the reaction graph of the network and is given by

$S+E\rightleftarrows SE\rightarrow P+E.$

In general, a chemical reaction network is defined by ${\mathcal {N}}=({\mathcal {S}},{\mathcal {C}},{\mathcal {R}})$ where ${\mathcal {S}}$ is the species set, ${\mathcal {C}}$ is the complex set, and ${\mathcal {R}}$ is the reaction set.

See also:

• Chemical reaction networks
• Indexing

The reaction graph

A chemical reaction network ${\mathcal {N}}=({\mathcal {S}},{\mathcal {C}},{\mathcal {R}})$ can be naturally represented as a directed graph G(V,E) where the set of vertices V is given by the complex set ${\mathcal {C}}$ and the set of directed edges E is given by the reaction set ${\mathcal {R}}$. This reaction graph keeps track of the connections (through reactions) between the nodes of a graph (the stoichiometric distinct complexes). Many topics from graph theory are relevant to the study of the reaction graph of a chemical reaction network, and consequently to the study of the dynamics of the dynamical systems corresponding to reaction networks.

We can express any reaction network via a reaction graph by writing each stoichiometric distinct complex exactly once as a node in a graph, and representing the reactions as connections between these nodes. For example, we can reformulate the following list of reactions (left) as the reaction graph (right) by writing each stoichiometrically distinct complex only once

${\begin{array}{rclccc}{\mathcal {A}}_{1}&\;\longrightarrow \;&{\mathcal {A}}_{2}&{\mathcal {A}}_{1}\;\longrightarrow \;{\mathcal {A}}_{2}\\[0.1cm]{\mathcal {A}}_{2}&\;\longrightarrow \;&{\mathcal {A}}_{3}+{\mathcal {A}}_{4}\;\;\;\;\;\;\;\;\;\;\Longrightarrow \;\;\;\;\;\;\;\;\;\;&\nwarrow \;\;\;\;\;\swarrow \\[0.1cm]{\mathcal {A}}_{3}+{\mathcal {A}}_{4}&\;\longrightarrow \;&{\mathcal {A}}_{1}&{\mathcal {A}}_{3}+{\mathcal {A}}_{4}.\end{array}}$

See also:

• Reaction graph
• Stoichiometry
• SR graph
• SCL graph

General kinetic framework

In order to translate a chemical reaction network into a dynamical framework, several assumptions need to be made. It is common to assume that the chemical species are in sufficient quantity to ignore stochastic effects (so that we may model the reaction continuously and deterministically over chemical concentrations) and that the medium where the reactions are occurring is well-mixed (so that there are no pockets of high concentrations of one species or another). (See the respective pages for information on stochastic and diffusive models.)

These assumptions lead to an ordinary differential equation model over the chemical concentrations. We let ${\vec {c}}=[c_{1},c_{2},\ldots ,c_{n}]^{T}$ denote the vector of chemical concentrations corresponding to the n chemical species $\left\{{\mathcal {A}}_{1},\ldots ,{\mathcal {A}}_{n}\right\}$. In order to determine the dynamics, we need to determine the net stoichiometric change resulting from each reaction and the corresponding rate of the reaction. We will number the reactions from $1,\ldots ,r$ and define ${\vec {y}}_{i}$ and ${\vec {y}}'_{i}$ to be the vectors of stoichiometric coefficients of the reactant and product complex of the $i^{th}$ reaction, respectively. We can then represent the dynamics of the chemical reaction network ${\mathcal {N}}$ as

${\frac {d{\vec {c}}}{dt}}=\sum _{i=1}^{r}({\vec {y}}_{i}'-{\vec {y}}_{i})\;R_{i}(t,{\vec {c}})$

where $R_{i}(t,{\vec {c}})\geq 0$ corresponds to the rate of the $i^{th}$ reaction. This can be more compactly written in matrix form as

${\frac {d{\vec {c}}}{dt}}=\Gamma \;{\vec {R}}(t,{\vec {c}})$

where the stoichiometric matrix $\Gamma \in \mathbb {Z} ^{n\times r}$ contains the reaction vectors ${\vec {y}}_{i}'-{\vec {y}}_{i}$, $i=1,\ldots ,r,$ as its columns and ${\vec {R}}(t,{\vec {c}})=[R_{1}(t,{\vec {c}}),\ldots ,R_{r}(t,{\vec {c}})]^{T}\in \mathbb {R} _{\geq 0}^{r}$.

It is important to note that this formulation does not depend on specific assumptions on the rate functions $R_{i}(t,{\vec {c}})$ of the individual reactions, and that those rates may be allowed to depend on time. They may be corresponded with a specific kinetic framework or be left general.

See also:

• Mass action kinetics
• k-variable kinetics
• Michaelis-Menten kinetics
• Hill kinetics

Mass action systems

The most commonly used form for the rate functions $R_{i}(t,{\vec {c}})$ is that of mass action kinetics. Under this assumption, each reaction is assigned a positive rate constant $k_{i}>0$, $i=1,\ldots ,r$, from a rate constant set ${\mathcal {K}}=\left\{k_{i}>0\;|\;i=1,\ldots ,r\right\}$. The rate functions are then defined according to $R_{i}(t,{\vec {c}})=k_{i}c_{1}^{y_{i1}}\cdots c_{n}^{y_{in}},i=1,\ldots ,r.$ If we introduce the shorthand notation ${\vec {c}}^{{\vec {y}}_{i}}=c_{1}^{y_{i1}}\cdots c_{n}^{y_{in}}$, we can write the system of ordinary differential equations governing the dynamics of a chemical reaction network under mass-action kinetics as

${\frac {d{\vec {c}}}{dt}}=\sum _{i=1}^{r}k_{i}({\vec {y}}_{i}'-{\vec {y}}_{i}){\vec {c}}^{{\vec {y}}_{i}}.$

The mass action system associated with the chemical reaction network ${\mathcal {N}}$ and the rate constant set ${\mathcal {K}}$ is given by the quadruple $({\mathcal {S}},{\mathcal {C}},{\mathcal {R}},{\mathcal {K}})$.

In order to most clearly emphasize the connection between the reaction graph of a chemical reaction network and the behaviour of solutions under mass-action kinetics, it is common to reformulate the dynamical representation of mass-action systems by breaking the stoichiometric matrix $\Gamma$ and rate vector ${\vec {R}}(t,{\vec {c}})$ into forms which emphasize the connections between stoichiometrically distinct complexes. To this end, we reindex the complexes according to their stoichiometric distinctiveness, ${\mathcal {C}}_{i}$, $i=1,\ldots ,m$, and define the complex matrix $Y\in \mathbb {Z} _{\geq 0}^{n\times m}$ to be the matrix where the $i^{th}$ column is given by the stoichiometric vector ${\vec {y}}_{i}$ corresponding to the $i^{th}$ complex and the mass-action vector $\Psi ({\vec {c}})\in \mathbb {R} _{\geq 0}^{m}$ to be the vector with entries $\Psi _{i}({\vec {c}})=c_{1}^{y_{i1}}\cdots c_{n}^{y_{in}}$ for $i=1,\ldots ,m$. We can write the mass-action system above as

${\frac {d{\vec {c}}}{dt}}=YA_{k}\Psi ({\vec {c}})$

where $A_{k}\in \mathbb {R} _{\geq 0}^{m\times m}$ is the kinetic or Kirchhoff matrix of the network.

Under the assumption of mass-action kinetics, chemical reaction networks give rise to a system of autonomous, polynomial ordinary differential equations (with no negative cross-effects). Despite this simple formulation, mass-action systems are known which exhibit varied and exotic behaviour, including stability, multistability, switching behaviour, finite escape time, oscillations, limit cycles, chaotic behaviour, and Hopf bifurcations. It is also worth noting that two mass-action systems with the same underlying reaction network ${\mathcal {N}}$ do not necessarily share the same qualitative behaviour (i.e. the dynamical behaviour of mass-action systems in general is rate-constant-dependent as well as network dependent).

See also:

• Detailed balanced systems
• Complex balanced systems
• Model reduction
• Polynomial ordinary differential equations
• Indexing