Skip to article frontmatterSkip to article content
Site not loading correctly?

This may be due to an incorrect BASE_URL configuration. See the MyST Documentation for reference.

Model of G-protein signalling pathway

using ComponentArrays: ComponentArray
using SimpleUnPack
using OrdinaryDiffEq
using DiffEqCallbacks
using SteadyStateDiffEq
using CairoMakie
function model605!(D, u, p, t)
    @unpack kRL, kRLm, kGa, kGd0, kG1, L, Rtotal, Gtotal = p
    @unpack RL, Ga, Gd = u
    R = Rtotal - RL
    G = Gtotal - Ga - Gd
    Gbg = Ga + Gd
    v1 = kRL * R * L - kRLm * RL
    v2 = kGa * G * RL
    v3 = kGd0 * Ga
    v4 = kG1 * Gd * Gbg
    D.RL = v1
    D.Ga = v2 - v3
    D.Gd = v3 - v4
    nothing
end
model605! (generic function with 1 method)
ps605 = ComponentArray(
    kRL = 2e6,
    kRLm = 0.01,
    kGa = 1e-5,
    kGd0 = 0.11,
    kG1 = 1.0,
    Rtotal = 4e3,
    Gtotal = 1e4,
    L = 0.0
)

u0605 = ComponentArray(
    RL = 0.0,
    Ga = 0.0,
    Gd = 0.0,
)
ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0)

Events

affect_L1!(integrator) = integrator.p.L = 1e-9
affect_L2!(integrator) = integrator.p.L = 0.0
event_L1 = PresetTimeCallback([200.0], affect_L1!)
event_L2 = PresetTimeCallback([800.0], affect_L2!)
cbs = CallbackSet(event_L1, event_L2)
SciMLBase.CallbackSet{Tuple{}, Tuple{SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L1!)}, typeof(Main.var"##248".affect_L1!), DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L1!)}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}, SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L2!)}, typeof(Main.var"##248".affect_L2!), DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L2!)}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}}}((), (SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L1!)}, typeof(Main.var"##248".affect_L1!), DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L1!)}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}(DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L1!)}([200.0], true, SciMLBase.INITIALIZE_DEFAULT, Main.var"##248".affect_L1!), Main.var"##248".affect_L1!, DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L1!)}([200.0], true, SciMLBase.INITIALIZE_DEFAULT, Main.var"##248".affect_L1!), SciMLBase.FINALIZE_DEFAULT, Bool[1, 1], nothing, ()), SciMLBase.DiscreteCallback{DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L2!)}, typeof(Main.var"##248".affect_L2!), DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L2!)}, typeof(SciMLBase.FINALIZE_DEFAULT), Nothing, Tuple{}}(DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L2!)}([800.0], true, SciMLBase.INITIALIZE_DEFAULT, Main.var"##248".affect_L2!), Main.var"##248".affect_L2!, DiffEqCallbacks.PresetTimeFunction{Vector{Float64}, typeof(SciMLBase.INITIALIZE_DEFAULT), typeof(Main.var"##248".affect_L2!)}([800.0], true, SciMLBase.INITIALIZE_DEFAULT, Main.var"##248".affect_L2!), SciMLBase.FINALIZE_DEFAULT, Bool[1, 1], nothing, ())))

Fig 6.5 A

tend = 1200.0
prob605 = ODEProblem(model605!, u0605, tend, ps605)
ODEProblem with uType ComponentArrays.ComponentVector{Float64, Vector{Float64}, Tuple{ComponentArrays.Axis{(RL = 1, Ga = 2, Gd = 3)}}} and tType Float64. In-place: true Non-trivial mass matrix: false timespan: (0.0, 1200.0) u0: ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0)
@time sol = solve(prob605, KenCarp47(), callback=cbs)
  3.775149 seconds (13.93 M allocations: 914.755 MiB, 16.32% gc time, 99.97% compilation time)
retcode: Success Interpolation: 3rd order Hermite t: 40-element Vector{Float64}: 0.0 9.999999999999999e-5 0.0010999999999999998 0.011099999999999997 0.11109999999999996 1.1110999999999995 11.111099999999993 111.11109999999994 200.0 200.0 ⋮ 800.0 800.0 847.9373930554989 870.6078939443786 932.2846541128085 993.9614142812385 1064.8232097437851 1137.6531912327312 1200.0 u: 40-element Vector{ComponentArrays.ComponentVector{Float64, Vector{Float64}, Tuple{ComponentArrays.Axis{(RL = 1, Ga = 2, Gd = 3)}}}}: ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ComponentVector{Float64}(RL = 0.0, Ga = 0.0, Gd = 0.0) ⋮ ComponentVector{Float64}(RL = 666.1695113001948, Ga = 570.9744183731231, Gd = 0.10997883843719922) ComponentVector{Float64}(RL = 666.1695113001948, Ga = 570.9744183731231, Gd = 0.10997883843719922) ComponentVector{Float64}(RL = 412.474565011872, Ga = 393.7393586748587, Gd = 0.10996920131464943) ComponentVector{Float64}(RL = 328.80656703058014, Ga = 317.20351895585287, Gd = 0.10996217899119244) ComponentVector{Float64}(RL = 177.45505023790008, Ga = 174.01237304857798, Gd = 0.1099354043025367) ComponentVector{Float64}(RL = 95.77149002625285, Ga = 94.75768279985633, Gd = 0.10987287516491091) ComponentVector{Float64}(RL = 47.15206849509245, Ga = 46.903573253928194, Gd = 0.10974312915168877) ComponentVector{Float64}(RL = 22.762475136013283, Ga = 22.70415890135508, Gd = 0.10947139979007932) ComponentVector{Float64}(RL = 12.202781720217143, Ga = 12.18602384669326, Gd = 0.10902503584768246)
fig = Figure()
ax = Axis(fig[1, 1], xlabel="Time", ylabel="Concentration", title="Fig 6.05 (A)")
lines!(ax, 0 .. tend, t -> sol(t).RL, label="RL")
lines!(ax, 0 .. tend, t -> sol(t).Ga, label="Ga")
axislegend(ax, position=:rc)
fig
Figure()

Fig 6.5 B

lrange = range(0, 20 * 1e-9, 101)

prob_func = (prob, i, repeat) -> remake(prob, p=ComponentArray(ps605; L=lrange[i]))
prob605b = SteadyStateProblem(model605!, u0605, ps605)
trajectories = length(lrange)
alg = DynamicSS(KenCarp47())
eprob = EnsembleProblem(prob605b; prob_func)
@time sim = solve(eprob, alg; trajectories);
  2.120891 seconds (10.99 M allocations: 750.412 MiB, 4.87% gc time, 189.79% compilation time)

ga = map(s->s.u.Ga, sim)
rl = map(s->s.u.RL, sim)
fig = Figure()
ax = Axis(
    fig[1, 1],
    xlabel = "Ligand (nM)",
    ylabel = "Steady-state abundance",
    title = "Fig. 6.5 (B)"
)

lines!(ax, lrange .* 1e9, ga, label="Ga")
lines!(ax, lrange .* 1e9, rl, label="RL")
axislegend(ax, position=:rc)
fig
Figure()

This notebook was generated using Literate.jl.

Textbook figures
Fig 6.3
Textbook figures
Fig 6.07