GitHub

Sensitivity Analysis of the Carcione2020 Epidemic Model

cd(@__DIR__)
using SBMLToolkit, ModelingToolkit, EasyModelAnalysis, UnPack

xmlfile = "../assets/Carcione2020.xml"

SBMLToolkit.checksupport_file(xmlfile)
mdl = readSBML(xmlfile, doc -> begin
    set_level_and_version(3, 2)(doc)
    convert_simplify_math(doc)
end)

rs = ReactionSystem(mdl)  # If you want to create a reaction system
odesys = convert(ODESystem, rs)  # Alternatively: ODESystem(mdl)

\[ \begin{align} \frac{\mathrm{d} \mathrm{Infected}\left( t \right)}{\mathrm{d}t} =& - \alpha \mathrm{Infected}\left( t \right) + \epsilon \mathrm{Exposed}\left( t \right) - \Gamma \mathrm{Infected}\left( t \right) \\ \frac{\mathrm{d} \mathrm{Exposed}\left( t \right)}{\mathrm{d}t} =& \frac{\beta \mathrm{Susceptible}\left( t \right) \mathrm{Infected}\left( t \right)}{\mathrm{Total}_{population}\left( t \right)} - \epsilon \mathrm{Exposed}\left( t \right) \\ \frac{\mathrm{d} \mathrm{Deceased}\left( t \right)}{\mathrm{d}t} =& \alpha \mathrm{Infected}\left( t \right) \\ \frac{\mathrm{d} \mathrm{Recovered}\left( t \right)}{\mathrm{d}t} =& \Gamma \mathrm{Infected}\left( t \right) \\ \frac{\mathrm{d} \mathrm{Total}_{population}\left( t \right)}{\mathrm{d}t} =& 0 \\ \frac{\mathrm{d} \mathrm{Susceptible}\left( t \right)}{\mathrm{d}t} =& \frac{ - \beta \mathrm{Susceptible}\left( t \right) \mathrm{Infected}\left( t \right)}{\mathrm{Total}_{population}\left( t \right)} \\ 0 =& - \lambda\left( t \right) + \frac{\mu \mathrm{Total}_{population}\left( t \right)}{City} \end{align} \]

sys = structural_simplify(odesys)

\[ \begin{align} \frac{\mathrm{d} \mathrm{Infected}\left( t \right)}{\mathrm{d}t} =& - \alpha \mathrm{Infected}\left( t \right) + \epsilon \mathrm{Exposed}\left( t \right) - \Gamma \mathrm{Infected}\left( t \right) \\ \frac{\mathrm{d} \mathrm{Exposed}\left( t \right)}{\mathrm{d}t} =& \frac{\beta \mathrm{Susceptible}\left( t \right) \mathrm{Infected}\left( t \right)}{\mathrm{Total}_{population}\left( t \right)} - \epsilon \mathrm{Exposed}\left( t \right) \\ \frac{\mathrm{d} \mathrm{Deceased}\left( t \right)}{\mathrm{d}t} =& \alpha \mathrm{Infected}\left( t \right) \\ \frac{\mathrm{d} \mathrm{Recovered}\left( t \right)}{\mathrm{d}t} =& \Gamma \mathrm{Infected}\left( t \right) \\ \frac{\mathrm{d} \mathrm{Total}_{population}\left( t \right)}{\mathrm{d}t} =& -0 \\ \frac{\mathrm{d} \mathrm{Susceptible}\left( t \right)}{\mathrm{d}t} =& \frac{ - \beta \mathrm{Susceptible}\left( t \right) \mathrm{Infected}\left( t \right)}{\mathrm{Total}_{population}\left( t \right)} \end{align} \]

@unpack Infected, Exposed, Deceased, Recovered, Total_population, Susceptible = sys
@unpack alpha, epsilon, gamma, mu, beta, City = sys
tspan = (0.0, 1.0)
prob = ODEProblem(odesys, [], tspan, [])
sol = solve(prob, Rodas5())
plot(sol, idxs = Deceased)
Example block output
pbounds = [
    alpha => [0.003, 0.006],
    epsilon => [1 / 6, 1 / 2],
    gamma => [0.1, 0.2],
    mu => [0.01, 0.02],
    beta => [0.7, 0.9]
]
create_sensitivity_plot(prob, 100.0, Deceased, pbounds; samples = 2000)
Example block output