About
editThis page assists a talk at the Julia@BGI community event at 11 of May 2023 at the Max Planck Institute for Biogeochemistry.
Logistics
The talk's materials are accessible on-wiki via https://w.wiki/6geW. You can also explore the page in slide form
- specifying and solving Ordinary Differential Equations (ODE) and variants
- i.e. Box-models
- part of Sciml organization
- specify symbolically
- MTK simplifies and generates fast code for DifferentialEqations.jl
- possible to compose systems of smaller components
function sesam3C(;name, k_N=60.0) @parameters t D = Differential(t) @parameters ϵ ϵ_tvr κ_E a_E m τ ... @variables (begin B(t), L(t), R(t), cumresp(t), ... end) eqs = [ D(B) ~ dB, dB ~ syn_B - tvr_B, D(L) ~ dL, dL ~ -dec_L + i_L, D(R) ~ dR, dR ~ -dec_R + ϵ_tvr*tvr_B + (1-κ_E)*syn_Enz, syn_Enz ~ a_E*B, tvr_Enz ~ syn_Enz, ... ] ODESystem(eqs; name) end function sesam3N(;name, sC = sesam3C(name=:sC)) @unpack L, R, B, dec_L, dec_R, i_L, ϵ_tvr, tvr_B, tvr_B0, syn_B, syn_Enz, tvr_Enz, r_tvr, κ_E = sC eqs = [ β_NL ~ L/L_N, β_NR ~ R/R_N, D(L_N) ~ dL_N, dL_N ~ -dec_L/β_NL + i_L/β_Ni, ] extend(ODESystem(eqs, t, sts, ps; name), sC) end
SIR:Implementierung in Octave:
edit% Modifiziertes SIR Modell, Implementierung in Octave %------------------- %Parameter B=55000 %Kapazität: 2/3 der deutschen Bevölkerung (in Tausend) %Raten c=0.3238 w=1/15; t=0.003; %Perzentualanteil der erfassten Infizierten r=0.1 times=(0:0.1:180); %Anfangswerte yo=[B-16/1000;16/1000;0]; % Funktion f der rechten Seite des DGL-Systems y'=f(y,t) f=@(y,x) [-c*y(1)*y(2)/(r*B);c*y(1)*y(2)/B-w*y(2);w*y(2) -d*y(2)]; y = lsode (f, yo, times); plot (times, y (:,1),'.',times, y (:,2), '.', times, y (:,3),'.',times, y (:,3)+y (:,2),'.') legend ("Succeptible","Infected", "Removed", "cummulative Infected ","location", "east") xlabel ('Tage') ylabel ('Befoelkerung (in Tausend)') axis([0,180, 0 ,55600])
Experiences
edit- nice to work with observables
- syn_Enz unsed for simplifying the system but computed only when quering the solution object
- Iterative computations not possible any more within derivative function
- Hard to debug
- cannot look/step at generated derivative code
- Position of state variables and parameters only known after simplifying the system
- Hassle on optimizing a subset of parameters
- still some problems with Automatic differentiation
- cannot use all the solvers
- fast response on issue, but hard to understand the aid
- Bayesian analysis with flexible statistical model
- used for inferring parameters and their uncertainty given observations
usage in SesamFitSPP.jl.git
edit- util_sample.jl
@model function sppfit(obs, ::Type{T} = Float64) where {T}
tmp
edita
b
c