Models with many conditions specific parameters
This tutorial demonstrates how to create a PEtabODEproblem
for a small ODE-model with many parameters to estimate. Specifically, the ODE-system has $\leq 20$ states and $\leq 20$ parameters, but there are approximately 70 parameters to estimate, since most parameters are specific to a subset of simulation conditions. For example, cond1 has a parameter τcond1, and cond2 has τcond2, which maps to the ODE-system parameter τ, respectively.
To run the code, you need the Beer PEtab files, which you can find here. You can also find a fully runnable example of this tutorial here.
First, we load the necessary libraries and read the model.
using PEtab
using OrdinaryDiffEq
using Printf
pathYaml = joinpath(@__DIR__, "Beer", "Beer_MolBioSystems2014.yaml")
petabModel = readPEtabModel(pathYaml, verbose=true)
PEtabModel for model Beer. ODE-system has 4 states and 9 parameters.
Generated Julia files are at ...
Handling condition-specific parameters
When dealing with small ODE-systems like Beer, the most efficient gradient method is gradientMethod=:ForwardDiff
. Additionally, we can compute the hessian via hessianMethod=:ForwardDiff
. However, since there are several condition-specific parameters in Beer, we have to perform as many forward-passes (solve the ODE model) as there are model-parameters when we compute the gradient and hessian via a single call to ForwardDiff.jl. To address this issue, we can use the option splitOverConditions=true
to force one ForwardDiff.jl call per simulation condition. This is most efficient for models where a majority of parameters are specific to a subset of simulation conditions.
For a model like Beer, the following options are thus recommended:
odeSolverOptions
- Rodas5P() (works well for smaller models with up to 15 states) and we use the defaultabstol, reltol .= 1e-8
.gradientMethod
- For small models like Beer, forward mode automatic differentiation (AD) is the fastest, so we choose:ForwardDiff
.hessianMethod
- For small models like Boehm with up to 20 parameters, it is computationally feasible to compute the full Hessian via forward-mode AD. Thus, we choose:ForwardDiff
.splitOverConditions=true
- This forces a call to ForwardDiff.jl per simulation condition.
odeSolverOptions = ODESolverOptions(Rodas5P(), abstol=1e-8, reltol=1e-8)
petabProblem = createPEtabODEProblem(petabModel, odeSolverOptions,
gradientMethod=:ForwardDiff,
hessianMethod=:ForwardDiff,
splitOverConditions=true)
p = petabProblem.θ_nominalT
gradient = zeros(length(p))
hessian = zeros(length(p), length(p))
cost = petabProblem.computeCost(p)
petabProblem.computeGradient!(gradient, p)
petabProblem.computeHessian!(hessian, p)
@printf("Cost = %.2f\n", cost)
@printf("First element in the gradient = %.2e\n", gradient[1])
@printf("First element in the hessian = %.2f\n", hessian[1, 1])
Cost = -58622.91
First element in the gradient = 7.17e-02
First element in the hessian = 755266.33