Generating reaction systems programmatically

Generating reaction systems programmatically#

There are two ways to create a reaction for a ReactionSystem:

  • Reaction() function

  • @reaction macro

using Catalyst
using ModelingToolkit

@parameters α K n δ γ β μ
@variables t
@species m₁(t) m₂(t) m₃(t) P₁(t) P₂(t) P₃(t)
\[\begin{split} \begin{equation} \left[ \begin{array}{c} m_1\left( t \right) \\ m_2\left( t \right) \\ m_3\left( t \right) \\ P_1\left( t \right) \\ P_2\left( t \right) \\ P_3\left( t \right) \\ \end{array} \right] \end{equation} \end{split}\]

Reaction function#

The Reaction(rate, substrates, products) function builds reactions. To allow for other stoichiometric coefficients we also provide a five argument form: Reaction(rate, substrates, products, substrate_stoichiometries, product_stoichiometries)

rxs = [
    Reaction(hillr(P₃, α, K, n), nothing, [m₁]),
    Reaction(hillr(P₁, α, K, n), nothing, [m₂]),
    Reaction(hillr(P₂, α, K, n), nothing, [m₃]),
    Reaction(δ, [m₁], nothing),
    Reaction(γ, nothing, [m₁]),
    Reaction(δ, [m₂], nothing),
    Reaction(γ, nothing, [m₂]),
    Reaction(δ, [m₃], nothing),
    Reaction(γ, nothing, [m₃]),
    Reaction(β, [m₁], [m₁, P₁]),
    Reaction(β, [m₂], [m₂, P₂]),
    Reaction(β, [m₃], [m₃, P₃]),
    Reaction(μ, [P₁], nothing),
    Reaction(μ, [P₂], nothing),
    Reaction(μ, [P₃], nothing)
]
15-element Vector{Reaction{Any, Int64}}:
 Catalyst.hillr(P₃(t), α, K, n), ∅ --> m₁
 Catalyst.hillr(P₁(t), α, K, n), ∅ --> m₂
 Catalyst.hillr(P₂(t), α, K, n), ∅ --> m₃
 δ, m₁ --> ∅
 γ, ∅ --> m₁
 δ, m₂ --> ∅
 γ, ∅ --> m₂
 δ, m₃ --> ∅
 γ, ∅ --> m₃
 β, m₁ --> m₁ + P₁
 β, m₂ --> m₂ + P₂
 β, m₃ --> m₃ + P₃
 μ, P₁ --> ∅
 μ, P₂ --> ∅
 μ, P₃ --> ∅
WARNING: both Symbolics and ModelingToolkit export "infimum"; uses of it in module Catalyst must be qualified
WARNING: both Symbolics and ModelingToolkit export "supremum"; uses of it in module Catalyst must be qualified

Use ReactionSystem(reactions, indenpendeent_variable) to collect these reactions. @named macro is used because every system in ModelingToolkit needs a name. @named x = System(...) is a short hand for x = System(...; name=:x)

@named repressilator = ReactionSystem(rxs, t)
\[\begin{split} \begin{align*} \require{mhchem} \varnothing &\xrightarrow{\frac{\alpha K^{n}}{K^{n} + P_3^{n}}} \mathrm{m_1} \\ \varnothing &\xrightarrow{\frac{\alpha K^{n}}{K^{n} + P_1^{n}}} \mathrm{m_2} \\ \varnothing &\xrightarrow{\frac{\alpha K^{n}}{K^{n} + P_2^{n}}} \mathrm{m_3} \\ \mathrm{m_1} &\xrightleftharpoons[\gamma]{\delta} \varnothing \\ \mathrm{m_2} &\xrightleftharpoons[\gamma]{\delta} \varnothing \\ \mathrm{m_3} &\xrightleftharpoons[\gamma]{\delta} \varnothing \\ \mathrm{m_1} &\xrightarrow{\beta} \mathrm{m_1} + \mathrm{P_1} \\ \mathrm{m_2} &\xrightarrow{\beta} \mathrm{m_2} + \mathrm{P_2} \\ \mathrm{m_3} &\xrightarrow{\beta} \mathrm{m_3} + \mathrm{P_3} \\ \mathrm{P_1} &\xrightarrow{\mu} \varnothing \\ \mathrm{P_2} &\xrightarrow{\mu} \varnothing \\ \mathrm{P_3} &\xrightarrow{\mu} \varnothing \end{align*} \end{split}\]

The @reaction macro provides the same syntax in the @reaction_network to build reactions. Note that @reaction macro only allows one-way reaction; reversible arrows are not allowed.

@variables t
@species P₁(t) P₂(t) P₃(t)

rxs = [(@reaction hillr($P₃, α, K, n),  --> m₁),
    (@reaction hillr($P₁, α, K, n),  --> m₂),
    (@reaction hillr($P₂, α, K, n),  --> m₃),
    (@reaction δ, m₁ --> ),
    (@reaction γ,  --> m₁),
    (@reaction δ, m₂ --> ),
    (@reaction γ,  --> m₂),
    (@reaction δ, m₃ --> ),
    (@reaction γ,  --> m₃),
    (@reaction β, m₁ --> m₁ + P₁),
    (@reaction β, m₂ --> m₂ + P₂),
    (@reaction β, m₃ --> m₃ + P₃),
    (@reaction μ, P₁ --> ),
    (@reaction μ, P₂ --> ),
    (@reaction μ, P₃ --> )]

@named repressilator2 = ReactionSystem(rxs, t)
\[\begin{split} \begin{align*} \require{mhchem} \varnothing &\xrightarrow{\frac{\alpha K^{n}}{K^{n} + P_3^{n}}} \mathrm{m_1} \\ \varnothing &\xrightarrow{\frac{\alpha K^{n}}{K^{n} + P_1^{n}}} \mathrm{m_2} \\ \varnothing &\xrightarrow{\frac{\alpha K^{n}}{K^{n} + P_2^{n}}} \mathrm{m_3} \\ \mathrm{m_1} &\xrightleftharpoons[\gamma]{\delta} \varnothing \\ \mathrm{m_2} &\xrightleftharpoons[\gamma]{\delta} \varnothing \\ \mathrm{m_3} &\xrightleftharpoons[\gamma]{\delta} \varnothing \\ \mathrm{m_1} &\xrightarrow{\beta} \mathrm{m_1} + \mathrm{P_1} \\ \mathrm{m_2} &\xrightarrow{\beta} \mathrm{m_2} + \mathrm{P_2} \\ \mathrm{m_3} &\xrightarrow{\beta} \mathrm{m_3} + \mathrm{P_3} \\ \mathrm{P_1} &\xrightarrow{\mu} \varnothing \\ \mathrm{P_2} &\xrightarrow{\mu} \varnothing \\ \mathrm{P_3} &\xrightarrow{\mu} \varnothing \end{align*} \end{split}\]