PEtab.jl API

PEtabModel

A PEtabModel for parameter estimation/inference can be created by importing a PEtab parameter estimation problem in the standard format, or it can be directly defined in Julia. For the latter, observables that link the model to measurement data are provided by PEtabObservable, parameters to estimate are defined by PEtabParameter, and any potential events (callbacks) are specified as PEtabEvent.

PEtab.PEtabObservableType
PEtabObservable(obs_formula, noise_formula; kwargs...)

Formulas defining the likelihood that links the model output to the measurement data.

obs_formula describes how the model output relates to the measurement data, while noise_formula describes the standard deviation (measurement error) and can be an equation or a numerical value. Both the observable and noise formulas can be a valid Julia equation. Variables used in these formulas must be either model species, model parameters, or parameters defined as PEtabParameter. The formulas can also include time-point-specific noise and observable parameters; for more information, see the documentation.

Keyword Argument

  • transformation: The transformation applied to the observable and its corresponding measurements. Valid options are :lin (normal measurement noise), :log, or :log10 (log-normal measurement noise). See below for more details.

Description

For a measurement y, an observable h = obs_formula, and a standard deviation σ = noise_formula, the PEtabObservable defines the likelihood that links the model output to the measurement data: $\pi(y \mid h, \sigma)$. For transformation = :lin, the measurement noise is assumed to be normally distributed, and the likelihood is given by:

\[\pi(y|h, \sigma) = \frac{1}{\sqrt{2\pi \sigma^2}}\mathrm{exp}\bigg( -\frac{(y - h)^2}{2\sigma^2} \bigg)\]

As a special case, if $\sigma = 1$, this likelihood reduces to the least-squares objective function. For transformation = :log, the measurement noise is assumed to be log-normally distributed, and the likelihood is given by:

\[\pi(y|h, \sigma) = \frac{1}{\sqrt{2\pi \sigma^2 y^2}}\mathrm{exp}\bigg( -\frac{\big(\mathrm{log}(y) - \mathrm{log}(h)\big)^2}{2\sigma^2} \bigg)\]

For transformation = :log10, the measurement noise is assumed to be log10-normally distributed, and the likelihood is given by:

\[\pi(y|h, \sigma) = \frac{1}{\sqrt{2\pi \sigma^2 y^2}\mathrm{log}(10) }\mathrm{exp}\bigg( -\frac{\big(\mathrm{log}_{10}(y) - \mathrm{log}_{10}(h)\big)^2}{2\sigma^2} \bigg)\]

For numerical stabillity, PEtab.jl works with the log-likelihood in practice.

source
PEtab.PEtabParameterType
PEtabParameter(x; kwargs...)

Parameter estimation information for parameter x.

All parameters to be estimated in a PEtabODEProblem must be declared as a PEtabParameter, and x must be the name of a parameter that appears in the model, observable formula, or noise formula.

Keyword Arguments

  • lb::Float64 = 1e-3: The lower parameter bound for parameter estimation.
  • ub::Float64 = 1e3: The upper parameter bound for parameter estimation.
  • scale::Symbol = :log10: The scale on which to estimate the parameter. Allowed options are :log10 (default), :log, and :lin. Estimating pest on the log10 scale typically improves performance and is recommended.
  • prior = nothing: An optional continuous univariate parameter prior distribution from Distributions.jl.
  • prior_on_linear_scale = true: Whether the prior is on the linear scale (default) or on the transformed scale. For example, if scale = :log10 and prior_on_linear_scale = false, the prior acts on the transformed value; log10(x).
  • estimate = true: Whether the parameter should be estimated (default) or treated as a constant.
  • value = nothing: Value to use if estimate = false, and value retreived by the get_x function. Defaults to the midpoint between lb and ub.

Description

If a prior $\pi(x_i)$ is provided, the parameter estimation problem becomes a maximum a posteriori problem instead of a maximum likelihood problem. Practically, instead of minimizing the negative log-likelihood,$-\ell(x)$, the negative posterior is minimized:

\[\min_{\mathbf{x}} -\ell(\mathbf{x}) - \sum_{i} \pi(x_i)\]

For all parameters $i$ with a prior.

source
PEtab.PEtabEventType
PEtabEvent(condition, affects, targets)

A model event triggered by condition that sets the value of targets to that of affects.

For a collection of examples with corresponding plots, see the documentation.

Arguments

  • condition: A Boolean expression that triggers the event when it transitions from false to true. For example, if t == c1, the event is triggered when the model time t equals c1. For S > 2.0, the event triggers when specie S passes 2.0 from below.
  • affects: An equation of of model species and parameters that describes the effect of the event. It can be a single expression or a vector if there are multiple targets.
  • targets: Model species or parameters that the event acts on. Must match the dimension of affects.
source

Then, given a dynamic model (as ReactionSystem or ODESystem), measurement data as a DataFrame, and potential simulation conditions as a Dict (see this tutorial), a PEtabModel can be created:

PEtab.PEtabModelType
PEtabModel(sys, observables::Dict{String, PEtabObservable}, measurements::DataFrame,
           parameters::Vector{PEtabParameter}; kwargs...)

From a ReactionSystem or an ODESystem model, observables that link the model to measurements and parameters to estimate, create a PEtabModel for parameter estimation.

For examples on how to create a PEtabModel, see the documentation.

See also PEtabObservable, PEtabParameter, and PEtabEvent.

Keyword Arguments

  • simulation_conditions = nothing: An optional dictionary specifying initial specie values and/or model parameters for each simulation condition. Only required if the model has multiple simulation conditions.
  • events = nothing: Optional model events (callbacks) provided as PEtabEvent. Multiple events should be provided as a Vector of PEtabEvent.
  • verbose::Bool = false: Whether to print progress while building the PEtabModel.
PEtabModel(path_yaml; kwargs...)

Import a PEtab problem in the standard format with YAML file at path_yaml into a PEtabModel for parameter estimation.

For examples on how to import a PEtab problem, see the documentation.

Keyword Arguments

  • ifelse_to_callback::Bool = true: Whether to rewrite ifelse (SBML piecewise) expressions to callbacks. This improves simulation runtime. It is strongly recommended to set this to true.
  • verbose::Bool = false: Whether to print progress while building the PEtabModel.
  • write_to_file::Bool = false: Whether to write the generated Julia functions to files in the same directory as the PEtab problem. Useful for debugging.
source

PEtabODEProblem

From a PEtabModel, a PEtabODEProblem can:

PEtab.PEtabODEProblemType
PEtabODEProblem(model::PEtabModel; kwargs...)

From model create a PEtabODEProblem for parameter estimation/inference.

If no options (kwargs) are provided, default settings are used for the ODESolver, gradient method, and Hessian method. For more information on the default options, see the documentation.

See also ODESolver, SteadyStateSolver, and PEtabModel.

Keyword Arguments

  • odesolver::ODESolver: ODE solver options for computing the likelihood (objective function).
  • odesolver_gradient::ODESolver: ODE solver options for computing the gradient. Defaults to odesolver if not provided.
  • ss_solver::SteadyStateSolver: Steady-state solver options for computing the likelihood. Only applicable for models with steady-state simulations.
  • ss_solver_gradient::SteadyStateSolver: Steady-state solver options for computing the gradient. Defaults to ss_solver if not provided.
  • gradient_method::Symbol: Method for computing the gradient. Available options and defaults are listed in the documentation.
  • hessian_method::Symbol: Method for computing the Hessian. As with the gradient, available options and defaults can be found in the documentation.
  • FIM_method=nothing: Method for computing the empirical Fisher Information Matrix (FIM). Accepts the same options as hessian_method.
  • sparse_jacobian=false: Whether to use a sparse Jacobian when solving the ODE model. This can greatly improve performance for large models.
  • sensealg: Sensitivity algorithm for gradient computations. Available and default options depend on gradient_method. See the documentation for details.
  • chunksize=nothing: Chunk size for ForwardDiff.jl when using forward-mode automatic differentiation for gradients and Hessians. If not provided, a default value is used. Tuning chunksize can improve performance but is non-trivial.
  • split_over_conditions::Bool=false: Whether to split ForwardDiff.jl calls across simulation conditions when computing the gradient and/or Hessian. This improves performance for models with many condition-specific parameters, otherwise it increases runtime.
  • reuse_sensitivities::Bool=false: Whether to reuse forward sensitivities from the gradient computation for the Gauss-Newton Hessian approximation. Only applies when hessian_method=:GaussNewton and gradient_method=:ForwardEquations. This can greatly improve performance when the optimization algorithm computes both the gradient and Hessian simultaneously.
  • verbose::Bool = false: Whether to print progress while building the PEtabODEProblem.

Returns

A PEtabODEProblem, where the key fields are:

  • nllh: Compute the negative log-likelihood function for an input vector x; nllh(x). For this function and the ones below listed below, the input x can be either a Vector or a ComponentArray.
  • grad!: Compute the in-place gradient of the nllh; grad!(g, x).
  • grad: Compute the out-of-place gradient of the nllh; g = grad(x).
  • hess!: Compute the in-place Hessian of the nllh; hess!(H, x).
  • hess: Compute the out-of-place Hessian of the nllh; H = hess(x).
  • FIM: Compute the out-of-place empirical Fisher Information Matrix (FIM) of the nllh; F = FIM(x).
  • chi2: Compute the chi-squared test statistic for the nllh (see mathematical definition below); χ² = chi2(x).
  • residuals: Computes the residuals between the measurement data and model output (see the mathematical definition below); r = residuals(x).
  • simulated_values: Computes the corresponding model values for each measurement in the measurements table, in the same order as the measurements appear.
  • lower_bounds: Lower parameter bounds for parameter estimation, as specified when creating the model.
  • upper_bounds: Upper parameter bounds for parameter estimation, as specified when creating the model.

Description

Following the PEtab standard, the objective function to be used for parameter estimation created by the PEtabODEProblem is a likelihood function, or, if priors are provided, a posterior function. The characteristics of the objective is defined in the PEtabModel. In practice, for numerical stability, a PEtabODEProblem works with the negative log-likelihood:

\[-\ell(\mathbf{x}) = - \sum_{i = 1}^N \ell_i(\mathbf{x}),\]

where $\ell_i$ is the likelihood for each measurement $i$. In addition, to accommodate numerical optimization packages, the PEtabODEProblem provides efficient functions for computing the gradient $-\nabla \ell(\mathbf{x})$ and the second derivative Hessian matrix $-\nabla^2 \ell(\mathbf{x})$, using the specified or default gradient_method and hessian_method.

In addition to $-\ell$ and its derivatives, the PEtabODEProblem provides functions for diagnostics and model selection. The χ² (chi-squared) value can be used for model selection [1], and it is computed as the sum of the χ² values for each measurement. For a measurement y, an observable h = obs_formula, and a standard deviation σ = noise_formula, the χ² is computed as:

\[\chi^2 = \frac{(y - h)^2}{\sigma^2}\]

The residuals $r$ can be used to assess the measurement error model and are computed as:

\[r = \frac{(y - h)}{\sigma}\]

Lastly, the empirical Fisher Information Matrix (FIM) can be used for identifiability analysis [2]. It should ideally be computed with an exact Hessian method. The inverse of the FIM provides a lower bound on the covariance matrix. While the FIM can be useful, the profile-likelihood approach generally yields better results for identifiability analysis [2].

  1. Cedersund et al, The FEBS journal, pp 903-922 (2009).
  2. Raue et al, Bioinformatics, pp 1923-1929 (2009).
source

A PEtabODEProblem has numerous configurable options. Two of the most important options are the ODESolver and, for models with steady-state simulations, the SteadyStateSolver:

PEtab.ODESolverType
ODESolver(solver, kwargs..)

ODE-solver options (solver, tolerances, etc.) used for solving the ODE model in a PEtabODEProblem.

Any solver from OrdinaryDiffEq.jl and Sundials.jl is supported. For solver recommendations and default options used when an ODESolver is not provided when creating a PEtabODEProblem, see the documentation.

More information on the available options and solvers can also be found in the documentation for DifferentialEquations.jl.

Keyword Arguments

  • abstol=1e-8: Absolute tolerance when solving the ODE model. It is not recommended to increase above 1e-6 for gradients in order to obtain accurate gradients.
  • reltol=1e-8: Absolute tolerance when solving the ODE model. It is not recommended to increase above 1e-6 for gradients in order to obtain accurate gradients.
  • dtmin=nothing: Minimum acceptable step size when solving the ODE model.
  • maxiters=10000: Maximum number of iterations when solving the ODE model. Increasing above the default value can cause parameter estimation to take substantial longer time.
  • verbose::Bool=true: Whether or not warnings are displayed if the solver exits early. true is recommended to detect if a suboptimal choice of solver.
  • adj_solver=solver: Solver to use when solving the adjoint ODE. Only applicable if gradient_method=:Adjoint when creating the PEtabODEProblem. Defaults to solver.
  • abstol_adj=abstol: Absolute tolerance when solving the adjoint ODE model. Only applicable if gradient_method=:Adjoint when creating the PEtabODEProblem. Defaults to abstol.
  • reltol_adj=abstol: Relative tolerance when solving the adjoint ODE model. Only applicable if gradient_method=:Adjoint when creating the PEtabODEProblem. Defaults to reltol.
source
PEtab.SteadyStateSolverType
SteadyStateSolver(method::Symbol; kwargs...)

Steady-state solver options (method, tolerances, etc.) for computing the steady state, where the ODE right-hand side f fulfills du = f(u, p, t) ≈ 0.

The steady state can be computed in two ways. If method=:Simulate, by simulating the ODE model until du = f(u, p, t) ≈ 0 using ODE solver options defined in the provied ODESolver to the PEtabODEProblem. This approach is **strongly** recommended. Ifmethod=:Rootfinding, by directly finding the root of the RHSf(u, p, t) = 0` using any algorithm from NonlinearSolve.jl. The root-finding approach is far less reliable than the simulation approach (see description below).

Keyword Arguments

  • termination_check = :wrms: Approach to check if the model has reached steady-state for method = :Simulate. Two approaches are supported:
    • :wrms: Weighted root-mean-square. Terminate when: $\sqrt{\frac{1}{N} \sum_i^N \Big( \frac{du_i}{reltol * u_i + abstol}\Big)^2} < 1$
    • :Newton: Terminate if the step for Newton's method Δu is sufficiently small: $\sqrt{\frac{1}{N} \sum_i^N \Big(\frac{\Delta u_i}{reltol * u_i + abstol}\Big)^2} < 1$ The :Newton approach requires that the Newton step Δu can be computed, which is only possible if the Jacobian of the RHS of the ODE model is invertible. If this is not the case, a pseudo-inverse is used if pseudoinverse = true, else wrms is used. The :Newton termination should perform better than :wrms as it accounts for the scale of the ODEs, but until benchmarks confirm this, we recommend :wrms.
  • pseudoinverse = true: Whether to use a pseudo-inverse if inverting the Jacobian fails for termination_check = :Newton.
  • rootfinding_alg = nothing: Root-finding algorithm from NonlinearSolve.jl to use if method = :Rootfinding. If left empty, the default NonlinearSolve.jl algorithm is used.
  • abstol: Absolute tolerance for the NonlinearSolve.jl root-finding problem or the termination_check formula. Defaults to 100 * abstol_ode, where abstol_ode is the tolerance for the ODESolver in the PEtabODEProblem.
  • reltol: Relative tolerance for the NonlinearSolve.jl root-finding problem or the termination_check formula. Defaults to 100 * reltol_ode, where reltol_ode is the tolerance for the ODESolver in the PEtabODEProblem.
  • maxiters: Maximum number of iterations to use if method = :Rootfinding.

Description

For an ODE model of the form:

\[\frac{\mathrm{d}u}{\mathrm{d}t} = du = f(u, p, t)\]

The steady state is defined by:

\[f(u, p, t) = 0.0\]

The steady state can be computed in two ways: either by simulating the ODE model from a set of initial values u₀ until du ≈ 0, or by using a root-finding algorithm such as Newton's method to directly solve f(u, p, t) = 0.0. While the root-finding approach is computationally more efficient (since it does not require solving the entire ODE), it is far less reliable and can converge to the wrong root. For example, in mass-action models with positive initial values, the feasible root should be positive in u. This is generally fulfilled when computing the steady state via simulation (a negative root typically only occurs if the ODE solver fails). However, with root-finding, there is no such guarantee, as any root that satisfies f(u, p, t) = 0.0 may be returned. Consistent with this, benchmarks have shown that simulation is the most reliable method [1].

Another alternative is to solve for the steady state symbolically. If feasible, this is the most computationally efficient approach [1].

  1. Fiedler et al, BMC system biology, pp 1-19 (2016)
source

PEtab.jl provides several functions for interacting with a PEtabODEProblem:

PEtab.get_xFunction
get_x(prob::PEtabODEProblem; linear_scale = false)::ComponentArray

Get the nominal parameter vector with parameters in the correct order expected by prob for parameter estimation/inference. Nominal values can optionally be specified when creating a PEtabParameter, or in the parameters table if the problem is provided in the PEtab standard format.

For ease of interaction (e.g., changing values), the parameter vector is returned as a ComponentArray. For how to interact with a ComponentArray, see the documentation and the ComponentArrays.jl documentation.

See also PEtabParameter.

Keyword argument

  • linear_scale: Whether to return parameters on the linear scale. By default, parameters are returned on the scale they are estimated, which by default is log10 as this often improves parameter estimation performance.
source
SciMLBase.remakeMethod
remake(prob::PEtabODEProblem, xchange::Dict)::PEtabODEProblem

Fix and consequently remove a set of parameters from being estimated without rebuilding prob.

xchange should be provided as a Dict with parameter names and their new values. For example, to fix k1 to 1.0, provide xchange = Dict(:k1 => 1.0). The parameters to be fixed can only be those that were originally set to be estimated.

remake is the most efficient approach for fixing parameters, as it does not re-compile the problem.

source

And additionally, functions for interacting with the underlying dynamic model (ODEProblem) within a PEtabODEProblem:

PEtab.get_u0Function
get_u0(res, prob::PEtabODEProblem; kwargs...)

Retrieve the ODEProblem initial values for simulating the ODE model in prob. res can be a parameter estimation result (e.g., PEtabMultistartResult) or a Vector with parameters in the order expected by prob (see get_x).

For information on keyword arguements see get_ps.

See also get_odeproblem and get_odesol.

source
PEtab.get_psFunction
get_ps(res, prob::PEtabODEProblem; kwargs...)

Retrieve the ODEProblem parameter values for simulating the ODE model in prob. res can be a parameter estimation result (e.g., PEtabMultistartResult) or a Vector with parameters in the order expected by prob (see get_x).

See also: get_u0, get_odeproblem, get_odesol.

Keyword Arguments

  • retmap=true: Whether to return the values as a map in the form [k1 => val1, ...]. Such a map can be directly used when building an ODEProblem. If false, a Vector is returned. This keyword is only applicable for get_u0 and get_ps.
  • cid::Symbol: Which simulation condition to return parameters for. If not provided, defaults to the first simulation condition. For other get functions, the ODEProblem, u0, or ODESolution for the specified cid is returned.
  • preeq_id: Which potential pre-equilibration (steady-state) simulation id to use. If a valid preeq_id is provided, the ODE is first simulated to steady state for preeq_id. Then the model shifts to cid, and the parameters for cid are returned. For other get functions, the ODEProblem, u0, or ODESolution for the specified cid is returned
source
PEtab.get_odeproblemFunction
get_odeproblem(res, prob::PEtabODEProblem; kwargs...) -> (sol, callbacks)

Retrieve the ODEProblem and callbacks (CallbackSet) for simulating the ODE model in prob. res can be a parameter estimation result (e.g., PEtabMultistartResult) or a Vector with parameters in the order expected by prob (see get_x).

For information on keyword arguements see get_ps.

See also: get_u0 and get_odesol.

source
PEtab.get_odesolFunction
get_odesol(res, prob::PEtabODEProblem; kwargs...)::ODESolution

Retrieve the ODESolution from simulating the ODE model in prob. res can be a parameter estimation result (e.g., PEtabMultistartResult) or a Vector with parameters in the order expected by prob (see get_x).

For information on keyword arguements see get_ps.

See also: get_u0 and get_odeproblem.

source
PEtab.solve_all_conditionsFunction
solve_all_conditions(x, prob::PEtabODEProblem, solver; kwargs)

Solve the ODE model in prob for all simulation conditions with the provided ODE-solver.

x should be a Vector or ComponentArray with parameters in the order expected by prob (see get_x).

Keyword Arguments

  • abstol=1e-8: Absolute tolerance for the ODE solver.
  • reltol=1e-8: Relative tolerance for the ODE solver.
  • maxiters=1e4: Maximum iterations for the ODE solver.
  • ntimepoints_save=0: The number of time points at which to save the ODE solution for each condition. A value of 0 means the solution is saved at the solvers default time points.
  • save_observed_t=false: When set to true, this option overrides ntimepoints_save and saves the ODE solution only at the time points where measurement data is available.

Returns

  • odesols: A dictionary containing the ODESolution for each simulation condition.
source

Parameter Estimation

A PEtabODEProblem contains all the necessary information for wrapping a suitable numerical optimization library, but for convenience, PEtab.jl provides wrappers for several available optimizers. In particular, single-start parameter estimation is provided via calibrate:

PEtab.calibrateFunction
calibrate(prob::PEtabODEProblem, x0, alg; kwargs...)::PEtabOptimisationResult

From starting point x0 using optimization algorithm alg, estimate unknown model parameters for prob, and get results as a PEtabOptimisationResult.

x0 can be a Vector or a ComponentArray, where the individual parameters must be in the order expected by prob. To get a vector in the correct order, see get_x.

A list of available and recommended optimization algorithms (alg) can be found in the documentation. Briefly, supported algorithms are from:

  • Optim.jl: LBFGS(), BFGS(), or IPNewton() methods.
  • Ipopt.jl: IpoptOptimizer() interior-point Newton method.
  • Fides.py: Fides() Newton trust region method.

Different ways to visualize the parameter estimation result can be found in the documentation.

See also PEtabOptimisationResult, Fides and IpoptOptimizer

Keyword Arguments

  • save_trace::Bool = false: Whether to save the optimization trace of the objective function and parameter vector. Only applicable for some algorithms; see the documentation for details.
  • options = DEFAULT_OPTIONS: Configurable options for alg. The type and available options depend on which package alg belongs to. For example, if alg = IPNewton() from Optim.jl, options should be provided as an Optim.Options() struct. A list of configurable options can be found in the documentation.
source
PEtab.PEtabOptimisationResultType
PEtabOptimisationResult

Parameter estimation statistics from single-start optimization with calibrate.

See also: calibrate

Fields

  • xmin: Minimizing parameter vector found by the optimization.
  • fmin: Minimum objective value found by the optimization.
  • x0: Starting parameter vector.
  • alg: Parameter estimation algorithm used.
  • niterations: Number of iterations for the optimization.
  • runtime: Runtime in seconds for the optimization.
  • xtrace: Parameter vector optimization trace. Empty if save_trace = false was provided to calibrate.
  • ftrace: Objective function optimization trace. Empty if save_trace = false was provided to calibrate.
  • converged: Convergence flag from alg.
  • original: Original result struct returned by alg. For example, if alg = IPNewton() from Optim.jl, original is the Optim return struct.
source

Multi-start (recommended method) parameter estimation, is provided via calibrate_multistart:

PEtab.calibrate_multistartFunction
calibrate_multistart(prob::PEtabODEProblem, alg, nmultistarts::Integer;
                     dirsave=nothing, kwargs...)::PEtabMultistartResult

Perform nmultistarts parameter estimation runs from randomly sampled starting points using the optimization algorithm alg to estimate the unknown model parameters in prob.

A list of available and recommended optimization algorithms (alg) can be found in the package documentation and calibrate documentation. If dirsave is provided, intermediate results for each run are saved in dirsave. It is strongly recommended to provide dirsave for larger models, as parameter estimation can take hours (or even days!), and without dirsave, all intermediate results will be lost if something goes wrong.

Different ways to visualize the parameter estimation result can be found in the documentation.

See also PEtabMultistartResult, get_startguesses, and calibrate.

Keyword Arguments

  • sampling_method = LatinHypercubeSample(): Method for sampling a diverse (spread out) set of starting points. See the documentation for get_startguesses, which is the function used for sampling.
  • sample_prior::Bool = true: See the documentation for get_startguesses.
  • seed = nothing: Seed used when generating starting points.
  • options = DEFAULT_OPTIONS: Configurable options for alg. See the documentation for calibrate.
source
PEtab.get_startguessesFunction
get_startguesses(prob::PEtabODEProblem, n::Integer; kwargs...)

Generate n random parameter vectors within the parameter bounds in prob.

If n = 1, a single random vector is returned. For n > 1, a vector of random parameter vectors is returned. In both cases, parameter vectors are returned as a ComponentArray. For details on how to interact with a ComponentArray, see the documentation and the ComponentArrays.jl documentation.

See also calibrate and calibrate_multistart.

Keyword Arguments

  • sampling_method = LatinHypercubeSample(): Method for sampling a diverse (spread out) set of parameter vectors. Any algorithm from QuasiMonteCarlo is allowed, but the default LatinHypercubeSample is recommended as it usually performs well.
  • sample_prior::Bool = true: Whether to sample random parameter values from the prior distribution if a parameter has a prior.
  • allow_inf::Bool = true: Whether to return parameter vectors for which the likelihood cannot be computed (typically happens because the ODEProblem cannot be solved). Often it only makes sense to use starting points with a computable likelihood for parameter estimation, hence it typically does not make sense to change this option.
source
PEtab.PEtabMultistartResultType
PEtabMultistartResult

Parameter estimation statistics from multi-start optimization with calibrate_multistart.

See also calibrate_multistart and PEtabOptimisationResult.

Fields

  • xmin: Best minimizer across all runs.
  • fmin: Best minimum across all runs.
  • alg: Parameter estimation algorithm.
  • nmultistarts: Number of parameter estimation runs.
  • sampling_method: Sampling method used for generating starting points.
  • dirsave: Path of directory where parameter estimation run statistics are saved if dirsave was provided to calibrate_multistart.
  • runs: Vector of PEtabOptimisationResult with the parameter estimation results for each run.
PEtabMultistartResult(dirres::String; which_run::String="1")

Import multistart parameter estimation results saved at dirres.

Each time a new optimization run is performed, results are saved with unique numerical endings. Results from a specific run can be retrieved by specifying the numerical ending with which_run.

source

Lastly, model selection is provided via petab_select:

PEtab.petab_selectFunction
petab_select(path_yaml, alg; nmultistarts = 100, kwargs...) -> path_res

For a PEtab-select problem, perform model selection with the method specified in the PEtab select problem files. Returns the path (path_res) to a YAML-file with model selection results.

The general model selection (e.g. to use forward-search) options are specified in the PEtab-select files. For details on how to set this up, see the PEtab-select documentation.

For each round of model selection, the candidate models are parameter estimated using multi-start parameter estimation, with nmultistarts performed for each model. The objective values obtained from parameter estimation are then used for the next round of model evaluation.

A list of available and recommended optimization algorithms (alg) can be found in the package documentation and calibrate documentation.

See also calibrate_multistart.

Keyword Arguments

source

For each case case, PEtab.jl supports the usage of optimization algorithms from Optim.jl, Ipopt.jl, and Fides.py:

PEtab.FidesType
Fides(hessian_method; verbose::Bool=false)

Setup the Fides box-constrained Newton-trust region optimizer for parameter estimation.

See also calibrate and calibrate_multistart.

Arguments

  • hessian_method: Method for computing the Hessian. Allowed options are:
    • nothing: The Hessian computed by the PEtabODEProblem is used.
    • :BB: Broyden's "bad" method.
    • :BFGS: Broyden-Fletcher-Goldfarb-Shanno update strategy.
    • :BG: Broyden's "good" method.
    • :Broyden: Broyden-class update scheme.
    • :SR1: Symmetric Rank 1 update.
    • :SSM: Structured Secant Method.
    • :TSSM: Totally Structured Secant Method.
  • verbose: Whether to print progress during the parameter estimation.

Description

Fides is a Newton-trust region optimizer for box-constrained optmization problems. More information on the algorithm can be found in [1]. Fides particularly excels when the full Hessian is too computationally expensive to compute, but a Gauss-Newton Hessian approximation can be computed (for more details see the documentation). In addition to supporting user Hessians via the PEtabODEProblem, it supports several Hessian approximation methods. Aa more extensive description than above see the Fides documentation.

  1. Fröhlich and Sorger, PLoS computational biology, pp e1010322 (2022)
source
PEtab.IpoptOptimizerType
IpoptOptimizer(LBFGS::Bool)

Setup the Ipopt Interior-point Newton method optmizer for parameter estimation.

Ipopt can be configured to use either the Hessian method from the PEtabODEProblem (LBFGS=false) or an LBFGS scheme (LBFGS=true). For setting other Ipopt options, see IpoptOptions.

See also calibrate and calibrate_multistart.

Description

Ipopt is an Interior-point Newton method for constrained non-linear optimization problems. More information on the algorithm can be found in [1].

  1. Wächter and Biegler, Mathematical programming, pp 25-57 (2006)
source
PEtab.IpoptOptionsType
IpoptOptions(; kwargs...)

Options for parameter estimation with IpoptOptimizer.

More details on the options can be found in the Ipopt documentation.

See also IpoptOptimizer, calibrate, and calibrate_multistart.

Keyword Arguments

  • print_level = 0: Output verbosity level (valid values are 0 ≤ print_level ≤ 12)
  • max_iter = 1000: Maximum number of iterations
  • tol = 1e-8: Relative convergence tolerance
  • acceptable_tol = 1e-6: Acceptable relative convergence tolerance
  • max_wall_time 1e20: Maximum wall time optimization is allowed to run
  • acceptable_obj_change_tol 1e20: Acceptance stopping criterion based on objective function change.
source

Parameter estimation results can be visualized using the plot-recipes detailed in this page, and with get_obs_comparison_plots:

PEtab.get_obs_comparison_plotsFunction
get_obs_comparison_plots(res, prob::PEtabODEProblem; kwargs...)::Dict

Plot the model fit against data for all simulation conditions and all observable ids for a PEtab.jl parameter estimation result (res).

Each entry in the returned Dict corresponds to plot(res, prob; obsids=[obsid], cid=cid, kwargs...) for all possible cid and all obsid.

source

As an alternative to the PEtab.jl interface to parameter estimation, a PEtabODEProblem can be converted to an OptimizationProblem to access the algorithms available via Optimization.jl:

PEtab.OptimizationProblemFunction
OptimizationProblem(prob::PEtabODEProblem; box_constraints::Bool = true)

Create an Optimization.jl OptimizationProblem from prob.

To use algorithms not compatible with box constraints (e.g., Optim.jl NewtonTrustRegion), set box_constraints = false. Note that with this option, optimizers may move outside the bounds, which can negatively impact performance. More information on how to use an OptimizationProblem can be found in the Optimization.jl documentation.

source

Bayesian Inference

PEtab.jl offers wrappers to perform Bayesian inference using state-of-the-art methods such as NUTS (the same sampler used in Turing.jl) or AdaptiveMCMC.jl. It should be noted that this part of PEtab.jl is planned to be moved to a separate package, so the syntax will change and be made more user-friendly in the future.

PEtab.PEtabLogDensityType

PEtabLogDensity(prob::PEtabODEProblem)

Create a LogDensityProblem using the posterior and gradient functions from prob.

This LogDensityProblem interface defines everything needed to perform Bayesian inference with packages such as AdvancedHMC.jl (which includes algorithms like NUTS, used by Turing.jl), and AdaptiveMCMC.jl.

source
PEtab.to_prior_scaleFunction
to_prior_scale(xpetab, target::PEtabLogDensity)

Transforms parameter x from the PEtab problem scale to the prior scale.

This conversion is needed for Bayesian inference, as in PEtab.jl Bayesian inference is performed on the prior scale.

Note

To use this function, the Bijectors, LogDensityProblems, and LogDensityProblemsAD packages must be loaded: using Bijectors, LogDensityProblems, LogDensityProblemsAD

source
PEtab.to_chainsFunction
to_chains(res, target::PEtabLogDensity; kwargs...)::MCMCChains

Converts Bayesian inference results obtained with PEtabLogDensity into an MCMCChains.

res can be the inference results from AdvancedHMC.jl or AdaptiveMCMC.jl. The returned chain has the parameters on the prior scale.

Keyword Arguments

  • start_time: Optional starting time for the inference, obtained with now().
  • end_time: Optional ending time for the inference, obtained with now().
Note

To use this function, the MCMCChains package must be loaded: using MCMCChains

source