Fig 3.03#
Michaelis-Menten kinetics
using OrdinaryDiffEq
using ModelingToolkit
using Catalyst
using Plots
Plots.default(linewidth=2)
Reaction network
rn = @reaction_network begin
(k1, km1), S + E <--> ES
k2, ES --> E + P
end
\[\begin{split} \begin{align*}
\mathrm{S} + \mathrm{E} &\xrightleftharpoons[\mathtt{km1}]{\mathtt{k1}} \mathrm{\mathtt{ES}} \\
\mathrm{\mathtt{ES}} &\xrightarrow{\mathtt{k2}} \mathrm{E} + \mathrm{P}
\end{align*}
\end{split}\]
setdefaults!(rn, [
:S=>5.,
:ES=>0.,
:E=>1.,
:P=>0.,
:k1 => 30.,
:km1 => 1.,
:k2 => 10.,
])
osys = convert(ODESystem, rn; remove_conserved = true) |> structural_simplify
\[\begin{split} \begin{align}
\frac{\mathrm{d} E\left( t \right)}{\mathrm{d}t} &= \mathtt{k2} \left( - E\left( t \right) + \Gamma_{1} \right) + \mathtt{km1} \left( - E\left( t \right) + \Gamma_{1} \right) - \mathtt{k1} E\left( t \right) S\left( t \right) \\
\frac{\mathrm{d} S\left( t \right)}{\mathrm{d}t} &= \mathtt{km1} \left( - E\left( t \right) + \Gamma_{1} \right) - \mathtt{k1} E\left( t \right) S\left( t \right)
\end{align}
\end{split}\]
observed(osys)
\[\begin{split} \begin{align}
\mathtt{ES}\left( t \right) &= - E\left( t \right) + \Gamma_{1} \\
P\left( t \right) &= E\left( t \right) - S\left( t \right) + \Gamma_{2}
\end{align}
\end{split}\]
tend = 1.0
1.0
prob = ODEProblem(osys, [], tend)
sol = solve(prob)
retcode: Success
Interpolation: 3rd order Hermite
t: 37-element Vector{Float64}:
0.0
0.006600703726699514
0.01022277056012863
0.015543467478140675
0.02051947168533509
0.02641018786018821
0.03261917715421904
0.03974230175678756
0.04773168899544981
0.0570841908779507
⋮
0.6470324023587559
0.6931031408477125
0.7367683827322631
0.7791059865493958
0.8212466198548
0.8642333668324144
0.9090993660547448
0.9568878660418696
1.0
u: 37-element Vector{Vector{Float64}}:
[1.0, 5.0]
[0.4180476197959369, 4.395302952712939]
[0.28277703335226406, 4.2362793319463865]
[0.1791117446616879, 4.091342551218704]
[0.13285191002575414, 4.002950476488607]
[0.10715692365888231, 3.9253138062440045]
[0.09568932433965065, 3.8580005265583313]
[0.09092109954691632, 3.788614104823389]
[0.0898206023516353, 3.7148221391871403]
[0.09058864991025081, 3.6304918096088135]
⋮
[0.7252277328768328, 0.04840338534337177]
[0.8069728851541074, 0.02310068905541521]
[0.866384281219644, 0.011833488316888211]
[0.9083028367930726, 0.006550597043571017]
[0.9376852665106246, 0.0038559086789705377]
[0.9582641290934262, 0.002355825567417385]
[0.9726435047202581, 0.0014581179992702214]
[0.9825982936657854, 0.0008953036977517524]
[0.9884434120994311, 0.0005834282632163945]
@unpack S, ES, E, P = osys
plot(sol, idxs=[S, ES, E, P], xlabel="Time (AU)", ylabel="Concentration (AU)", legend=:right, title="Fig 3.03")

rn303mm = @reaction_network begin
mm(S, k2 * ET, (km1 + k2) / k1), S ⇒ P ## using \Rightarrow
end
\[ \begin{align*}
\mathrm{S} &\xrightarrow{\frac{\mathtt{ET} S \mathtt{k2}}{S + \frac{\mathtt{k2} + \mathtt{km1}}{\mathtt{k1}}}} \mathrm{P}
\end{align*}
\]
setdefaults!(rn303mm, [
:S=>5.,
:ET=>1.,
:P=>0.,
:k1 => 30.,
:km1 => 1.,
:k2 => 10.,
])
osysmm = convert(ODESystem, rn303mm; remove_conserved = true) |> structural_simplify
\[ \begin{align}
\frac{\mathrm{d} S\left( t \right)}{\mathrm{d}t} &= \frac{ - \mathtt{ET} \mathtt{k2} S\left( t \right)}{\frac{\mathtt{k2} + \mathtt{km1}}{\mathtt{k1}} + S\left( t \right)}
\end{align}
\]
tend = 1.0
probmm = ODEProblem(osysmm, [], tend)
solmm = solve(probmm)
retcode: Success
Interpolation: 3rd order Hermite
t: 15-element Vector{Float64}:
0.0
0.08829951015219592
0.2967725537625472
0.5000004875811703
0.6196473972837556
0.6670104453577671
0.7099275100893456
0.7465737578395125
0.7827003202641042
0.8170029171212273
0.8515324938674043
0.8866728375768105
0.9240249128338719
0.9646817216791721
1.0
u: 15-element Vector{Vector{Float64}}:
[5.0]
[4.182468228248958]
[2.3146741114200724]
[0.7140966653045732]
[0.13296502509373836]
[0.04630250749684772]
[0.01563176034250213]
[0.005910624536905458]
[0.0022299117086426626]
[0.0008784935038914966]
[0.00034320139778785724]
[0.00013175527102757075]
[4.761595514034663e-5]
[1.5733007241346686e-5]
[6.0074984467007344e-6]
@unpack S, P = osys
fig = plot(sol, idxs=[S, P], line=(:dash), label=["S (full)" "P (full)"])
plot!(fig, solmm, idxs=[S, P], label=["S (MM)" "P (MM)"])
plot!(fig, title="Fig. 3.03",
xlabel="Time (AU)", ylabel="Concentration (AU)",
xlims=(0., tend), ylims=(0., 5.), legend=:right
)

This notebook was generated using Literate.jl.