Learn practical skills, build real-world projects, and advance your career
Updated 4 years ago
#Load all the packages.
using DifferentialEquations,Plots,DiffEqBiological
using Latexify
#constants
global growthrate = 1. #growth rate
global N_0 = 10^20 #population cap
@reaction_func growth(x) = growthrate*x*(1 - (A+B+B1+C)/N_0)
model = @reaction_network begin
#model equations
growth(A), 0 → A
growth(B), 0 → B
growth(B1), 0 → B1
growth(C), 0 → C
0.01, A → B
0.02, B → A
0.01, A → B1
0.02, B1 → A
0.001,A → 0
0.001,B → G
0.001,C → 0
0.001,B1 → G
0.002,G → 0
0.05, C + G → B1
end
┌ Warning: The RegularJump interface has changed to be matrix-free. See the documentation for more details.
└ @ DiffEqJump /home/ved/.julia/packages/DiffEqJump/sEB5p/src/jumps.jl:43
(::reaction_network) (generic function with 2 methods)
latexify(model,env=:chemical)
\begin{align}
\require{mhchem}
\ce{ \varnothing &->[growthrate \cdot A \cdot \left( 1 - \frac{A + B + B1 + C}{N_{0}} \right)] A}\\
\ce{ \varnothing &->[growthrate \cdot B \cdot \left( 1 - \frac{A + B + B1 + C}{N_{0}} \right)] B}\\
\ce{ \varnothing &->[growthrate \cdot B1 \cdot \left( 1 - \frac{A + B + B1 + C}{N_{0}} \right)] B1}\\
\ce{ \varnothing &->[growthrate \cdot C \cdot \left( 1 - \frac{A + B + B1 + C}{N_{0}} \right)] C}\\
\ce{ A &<=>[0.01][0.02] B}\\
\ce{ A &<=>[0.01][0.02] B1}\\
\ce{ A &->[0.001] \varnothing}\\
\ce{ B &->[0.001] G}\\
\ce{ C &->[0.001] \varnothing}\\
\ce{ B1 &->[0.001] G}\\
\ce{ G &->[0.002] \varnothing}\\
\ce{ C + G &->[0.05] B1}
\end{align}
latexify(model)
\begin{align}
\frac{dA(t)}{dt} =& growthrate \cdot A \cdot \left( 1 - \frac{A + B + B1 + C}{N_{0}} \right) - 0.01 \cdot A + 0.02 \cdot B - 0.01 \cdot A + 0.02 \cdot B1 - 0.001 \cdot A \\
\frac{dB(t)}{dt} =& growthrate \cdot B \cdot \left( 1 - \frac{A + B + B1 + C}{N_{0}} \right) + 0.01 \cdot A - 0.02 \cdot B - 0.001 \cdot B \\
\frac{dB1(t)}{dt} =& growthrate \cdot B1 \cdot \left( 1 - \frac{A + B + B1 + C}{N_{0}} \right) + 0.01 \cdot A - 0.02 \cdot B1 - 0.001 \cdot B1 + 0.05 \cdot C \cdot G \\
\frac{dC(t)}{dt} =& growthrate \cdot C \cdot \left( 1 - \frac{A + B + B1 + C}{N_{0}} \right) - 0.001 \cdot C - 0.05 \cdot C \cdot G \\
\frac{dG(t)}{dt} =& 0.001 \cdot B + 0.001 \cdot B1 - 0.002 \cdot G - 0.05 \cdot C \cdot G
\end{align}
tspan = (0.,500.)
initial = [0.,1.0 * 10^6,0.,1.0*10^4,0]
prob = ODEProblem(model,initial,tspan)
sol = solve(prob,maxiters = 1e10,saveat = 0.1);
timepoints = [0.0:0.1:500.0]
a = [sol[i][1] for i = 1:length(sol)];
b = [sol[i][2] for i = 1:length(sol)];
b1 = [sol[i][3] for i = 1:length(sol)];
c = [sol[i][4] for i = 1:length(sol)];
g = [sol[i][5] for i = 1:length(sol)];