Parameter Estimation of DelayDiffEq (and other Non-ODE) Problems

Estimation
1

I am currently trying to model and simulate a system of delay differential equations(DDEs) in Pumas. While I can currently simulate, when I try to do parameter estimation I get the error
“Incompatible problem+solver pairing”, which I presume means that Pumas does not currently support non-ode parameter estimation.

The closest advice I’ve seen for non-ode parameter estimation portion came from @ChrisRackauckas in a Pumas discussion thread from July 2020 found here. Are there any updates on if/when Pumas will support non-ode parameter estimation? Also, the aforementioned thread said there is a workaround for parameter estimation of SDE problems using another Julia Library. Does such a workaround exist for the parameter estimation of a system of DDEs? Thank you in advance.

2

This error means that an ODE solver was specified on a non-ODE. When you did the fit call, did you pass a DDE solver algorithm?

3

Apologies, yeah that was the issue, thanks!

Now I’m onto a new issue with the fit, specifically I am able to simulate it correctly with the following settings:

diffeq_options = (alg = MethodOfSteps(Rodas4()),
                    reltol = 1e-15, abstol = 1e-14, dtmax = 0.001,
                    saveat = 0:600, maxiters = 1e7, tstops = discontinuities)

However, applying the same settings to the fit (FOCE) results in an error that states the dde is too stiff to solve.

Increasing maxiters to 1e9 doesn’t fix it (haven’t tested past that since 1e9 already takes a while). Dropping both of the tolerances to 1e-3 gets a different error:

┌ Info: Checking the initial parameter values.
└ @ Pumas /build/_work/PumasSystemImages/PumasSystemImages/julia_depot/packages/Pumas/9AGx1/src/estimation/likelihoods.jl:3995

MethodError: no method matching _valid_foce_distribution(::Matrix{ForwardDiff.Dual{ForwardDiff.Tag{Pumas.Tag, Float64}, Float64, 1}})

I get a similar error when using LaplaceI.

Any idea what I could be setting up incorrectly? At the moment, everything in the model aspect is exactly the same as the simulation. Below is my current fit function. Let me know if you need to see other aspects of the code. Thank you again.

ada_fit = fit(p_model,ada_fig5a_prefit,init_params(p_model), Pumas.FOCE();
                optimize_fn=Pumas.DefaultOptimizeFN(show_trace=true), 
                diffeq_options=diffeq_options,
                ensemblealg = EnsembleDistributed())
4

It could look like you are missing a dot or @. in the derived block. Notice that it is checking a Matrix in _valid_foce_distribution(::Matrix...). E.g.

julia> Central = [1.0, 2.0, 3.0];

julia> Vc = [1.1, 1.2, 1.1];

julia> Central / Vc
3×3 Matrix{Float64}:
 0.284974  0.310881  0.284974
 0.569948  0.621762  0.569948
 0.854922  0.932642  0.854922

vs

julia> @. Central / Vc
3-element Vector{Float64}:
 0.9090909090909091
 1.6666666666666667
 2.727272727272727