Simple example of ODE parameters inference

[denisshepelin]

Dear all, I'm trying to work with toy ODE that is predator prey.

using DifferentialEquations, Plots, ParameterizedFunctions, Random, Distributions, Turing, StatsPlots

predator_prey = @ode_def begin
    du = (alpha - beta*v) * u
    dv = (-gamma + delta * u) * v
    end alpha beta gamma delta

#  "Experimental" data generation
u0 = [10.0, 10.0]
tspan = (0.0, 40.0)
p = [1.1,0.4,0.4, 0.1]
prob = ODEProblem(predator_prey, u0, tspan, p)
sol = solve(prob, saveat=0.5)
# optional plot
# plot(sol, linestyle=:dot)

# modify data with adding gaussian noise
d1 = Normal(0, 0.5)
d2 = Normal(0, 0.5)
newsol = sol
newsol[1,:] = sol[1,:] + rand(d1,81)
newsol[2,:] = sol[2,:] + rand(d2,81)

# define statistical model
@model bayes_predator_prey(x) = begin
  alpha ~ Normal(1.0, 0.2)
  beta ~ Normal(0.5, 0.2)
  gamma ~ Normal(0.5, 0.2)
  delta ~ Normal(0.1, 0.1)
  # gather parameters and solve equation
  p = [alpha, beta, gamma ,delta]
  u0 = [10.0, 10.0]
  tspan = (0.0, 40.0)
  m = solve(ODEProblem(predator_prey, u0, tspan, p), saveat=0.5)
  # hope that it's likelihood
  N = length(x)
  for i in 1:N  
    x[i] ~ MvNormal(m[i], [0.5,0.5])
  end
end

# run the sampling 
chain = sample(bayes_predator_prey(newsol), HMC(1000, 0.1, 5))

I get the following error -

TypeError: in typeassert, expected Float64, got ForwardDiff.Dual{Nothing,Float64,4}

Am I defining the model incorrectly or there is some issue with ODE problem formulation that prevents me from sampling?

[cpfiffer]

I think this is because the solve function doesn't know what to do with Turing's automatic differentiation. I'm also not really the AD guy, so I'm tagging the autodiff people: @mohamed82008 @willtebbutt.

I don't know if it's applicable to your model, but Turing does have some interface with DifferentialEquations --- check out DiffEqBayes.jl and see if it might help.

[denisshepelin]

Thank you for your help @cpfiffer! I've looked at DiffEqBayes.jl implementation and Turing is supported there only with Importance Sampler (Turing.IS()) regime. I would really like to use HMC or NUTS for my work.

Regarding AD for Differential Equations solve - I've seen package DiffEqFlux which supposed to work with neural networks and then I assume should be differentiable with Flux.jl. Unfortunately I'm not advanced enough in Julia yet and I was stuck trying to get it work with Turing.

I would appreciate any comments where and how can I fill a better issue describing my problem.

[mohamed82008]

Sorry missed the ping. I will look into it.

[xukai92]

It should work with AD. At least I played with some models around 1 month ago. It's likely you will need to do something like u0 = convert.(typeof(alpha), (u0)) as that's what I see from a notebook I was helping. (I will ask who owns the notebook if he is happy to share the code.)

[jonasmac16]

Have a look at http://docs.juliadiffeq.org/latest/analysis/sensitivity.html#Local-Forward-Sensitivity-Analysis-1 where they describe how to make ODEs work with AD. The solver needs to use Dual numbers throughout for it to work with AD which is why you need to convert the type of the initial conditions. Also have a look here https://github.com/JuliaDiffEq/DiffEqBayes.jl/blob/master/src/turing_inference.jl which gives you a good template of a Turing model to start from as it handles unstable ODE solution for which the solver crashes or never finished gracefully.

I have adapted your code to implement this. I had to change the solver as it struggled to solve the problem within 2000 iter. You might also stick with Tsit5() and increase maxiter and play around with abstol and reltol.

using DifferentialEquations, Plots, ParameterizedFunctions, Random, Distributions, Turing, StatsPlots, RecursiveArrayTools

predator_prey = @ode_def begin
    du = (alpha - beta*v) * u
    dv = (-gamma + delta * u) * v
    end alpha beta gamma delta

#  "Experimental" data generation
u0 = [10.0, 10.0]
tspan = (0.0, 40.0)
p = [1.1,0.4,0.4, 0.1]
tp= collect(range(0,stop=40))
prob = ODEProblem(predator_prey, u0, tspan, p)
sol = solve(prob, saveat=tp)
# optional plot
# plot(sol, linestyle=:dot)
# modify data with adding gaussian noise
d1 = Normal(0, 0.5)
d2 = Normal(0, 0.5)
newsol = VectorOfArray(sol.u)[:,:]

newsol[1,:] = newsol[1,:] + rand(d1,length(tp))
newsol[2,:] = newsol[2,:] + rand(d2,length(tp))

function test_f(p, tp_save)
  _prob = remake(prob;u0=convert.(eltype(p),prob.u0),p=p)
  solve(_prob,Rodas4(),saveat=tp_save;abstol=1e-6, reltol= 1e-6)
end

plot(sol)
scatter!(tp,newsol[1,:])
scatter!(tp,newsol[2,:])

# define statistical model
@model bayes_predator_prey(x, tp_x) = begin
  alpha ~ Normal(1.0, 0.2)
  beta ~ Normal(0.5, 0.2)
  gamma ~ Normal(0.5, 0.2)
  delta ~ Normal(0.1, 0.1)
  # gather parameters and solve equation
  p = [alpha, beta, gamma ,delta]
  sol_tmp = test_f(p, tp_x)
  # hope that it's likelihood
  N = length(tp_x)

  fill_length = length(tp_x) - length(sol_tmp.u)

  for i in 1:fill_length
    if eltype(sol_tmp.u) <: Number
      push!(sol_tmp.u, Inf)
    else
      push!(sol_tmp.u, fill(Inf, size(sol_tmp[1])))
    end
  end

  for i in 1:N
    x[:,i] ~ MvNormal(sol_tmp.u[i], [0.5,0.5])
  end
end

# run the sampling
chain = sample(bayes_predator_prey(newsol, tp), HMC(1000, 0.1, 5))

[denisshepelin]

Thank you so much! Unfortunately the code by @jonasmac16 is not working fully for me and results in error which seems to be linked to DifferentialEquations.jl package:

ERROR: MethodError: no method matching OrdinaryDiffEq.Tsit5ConstantCache(::Float64, ::Float64, ::Float64, ::Float64, ::Float64, ::Float64, ::Rational{Int64}, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::Rational{Int64}, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::Rational{Int64}, ::BigFloat, ::BigFloat, ::BigFloat, ::Float64, ::BigFloat, ::BigFloat, ::Float64, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::BigFloat, ::Float64, ::Int64, ::Float64)
Closest candidates are:
  OrdinaryDiffEq.Tsit5ConstantCache(::T2, ::T2, ::T2, ::T2, ::T2, ::T2, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T, ::T) where {T, T2} at /Users/denshe/.julia/packages/OrdinaryDiffEq/tat1c/src/tableaus/low_order_rk_tableaus.jl:498
Stacktrace:
 [1] OrdinaryDiffEq.Tsit5ConstantCache(::Type, ::Type) at /Users/denshe/.julia/packages/OrdinaryDiffEq/tat1c/src/tableaus/low_order_rk_tableaus.jl:647
 [2] alg_cache(::Tsit5, ::Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1}, ::Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1}, ::Type, ::Type, ::Type, ::Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1}, ::Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1}, ::getfield(Main, Symbol("##363")){getfield(Main, Symbol("##3#7")),getfield(Main, Symbol("##4#8")),getfield(Main, Symbol("##5#9")),Nothing,Nothing,getfield(Main, Symbol("##6#10")),Expr,Expr}, ::Float64, ::Float64, ::Float64, ::Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1}, ::Bool, ::Type{Val{true}}) at /Users/denshe/.julia/packages/OrdinaryDiffEq/tat1c/src/caches/low_order_rk_caches.jl:348
 [3] #__init#297(::Array{Int64,1}, ::Array{Float64,1}, ::Array{Float64,1}, ::Nothing, ::Bool, ::Nothing, ::Bool, ::Bool, ::Bool, ::Nothing, ::Bool, ::Bool, ::Float64, ::Bool, ::Rational{Int64}, ::Float64, ::Float64, ::Int64, ::Rational{Int64}, ::Int64, ::Int64, ::Rational{Int64}, ::Bool, ::Int64, ::Nothing, ::Nothing, ::Int64, ::Float64, ::Float64, ::typeof(DiffEqBase.ODE_DEFAULT_NORM), ::typeof(LinearAlgebra.opnorm), ::typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), ::typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Bool, ::Int64, ::String, ::typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), ::Nothing, ::Bool, ::Bool, ::Bool, ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(DiffEqBase.__init), ::ODEProblem{Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},getfield(Main, Symbol("##363")){getfield(Main, Symbol("##3#7")),getfield(Main, Symbol("##4#8")),getfield(Main, Symbol("##5#9")),Nothing,Nothing,getfield(Main, Symbol("##6#10")),Expr,Expr},Nothing,DiffEqBase.StandardODEProblem}, ::Tsit5, ::Array{Any,1}, ::Array{Any,1}, ::Array{Any,1}, ::Type{Val{true}}) at /Users/denshe/.julia/packages/OrdinaryDiffEq/tat1c/src/solve.jl:242
 [4] (::getfield(DiffEqBase, Symbol("#kw##__init")))(::NamedTuple{(:saveat, :abstol, :reltol),Tuple{Array{Int64,1},Float64,Float64}}, ::typeof(DiffEqBase.__init), ::ODEProblem{Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},getfield(Main, Symbol("##363")){getfield(Main, Symbol("##3#7")),getfield(Main, Symbol("##4#8")),getfield(Main, Symbol("##5#9")),Nothing,Nothing,getfield(Main, Symbol("##6#10")),Expr,Expr},Nothing,DiffEqBase.StandardODEProblem}, ::Tsit5, ::Array{Any,1}, ::Array{Any,1}, ::Array{Any,1}, ::Type{Val{true}}) at ./none:0
 [5] #__solve#296(::Base.Iterators.Pairs{Symbol,Any,Tuple{Symbol,Symbol,Symbol},NamedTuple{(:saveat, :abstol, :reltol),Tuple{Array{Int64,1},Float64,Float64}}}, ::Function, ::ODEProblem{Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},getfield(Main, Symbol("##363")){getfield(Main, Symbol("##3#7")),getfield(Main, Symbol("##4#8")),getfield(Main, Symbol("##5#9")),Nothing,Nothing,getfield(Main, Symbol("##6#10")),Expr,Expr},Nothing,DiffEqBase.StandardODEProblem}, ::Tsit5, ::Array{Any,1}, ::Array{Any,1}, ::Array{Any,1}, ::Type{Val{true}}) at /Users/denshe/.julia/packages/OrdinaryDiffEq/tat1c/src/solve.jl:6
 [6] (::getfield(DiffEqBase, Symbol("#kw##__solve")))(::NamedTuple{(:saveat, :abstol, :reltol),Tuple{Array{Int64,1},Float64,Float64}}, ::typeof(DiffEqBase.__solve), ::ODEProblem{Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},getfield(Main, Symbol("##363")){getfield(Main, Symbol("##3#7")),getfield(Main, Symbol("##4#8")),getfield(Main, Symbol("##5#9")),Nothing,Nothing,getfield(Main, Symbol("##6#10")),Expr,Expr},Nothing,DiffEqBase.StandardODEProblem}, ::Tsit5, ::Array{Any,1}, ::Array{Any,1}, ::Array{Any,1}, ::Type{Val{true}}) at ./none:0 (repeats 5 times)
 [7] #solve#456(::Base.Iterators.Pairs{Symbol,Any,Tuple{Symbol,Symbol,Symbol},NamedTuple{(:saveat, :abstol, :reltol),Tuple{Array{Int64,1},Float64,Float64}}}, ::Function, ::ODEProblem{Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},getfield(Main, Symbol("##363")){getfield(Main, Symbol("##3#7")),getfield(Main, Symbol("##4#8")),getfield(Main, Symbol("##5#9")),Nothing,Nothing,getfield(Main, Symbol("##6#10")),Expr,Expr},Nothing,DiffEqBase.StandardODEProblem}, ::Tsit5) at /Users/denshe/.julia/packages/DiffEqBase/6ewdP/src/solve.jl:39
 [8] (::getfield(DiffEqBase, Symbol("#kw##solve")))(::NamedTuple{(:saveat, :abstol, :reltol),Tuple{Array{Int64,1},Float64,Float64}}, ::typeof(solve), ::ODEProblem{Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},Tuple{Float64,Float64},true,Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1},getfield(Main, Symbol("##363")){getfield(Main, Symbol("##3#7")),getfield(Main, Symbol("##4#8")),getfield(Main, Symbol("##5#9")),Nothing,Nothing,getfield(Main, Symbol("##6#10")),Expr,Expr},Nothing,DiffEqBase.StandardODEProblem}, ::Tsit5) at ./none:0
 [9] test_f(::Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1}, ::Array{Int64,1}) at ./REPL[11]:3
 [10] macro expansion at ./REPL[15]:8 [inlined]
 [11] (::getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}})(::Turing.Core.RandomVariables.UntypedVarInfo, ::Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}, ::Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}}) at /Users/denshe/.julia/packages/Turing/RZOZ8/src/core/compiler.jl:388
 [12] #call#3 at /Users/denshe/.julia/packages/Turing/RZOZ8/src/Turing.jl:62 [inlined]
 [13] Model at /Users/denshe/.julia/packages/Turing/RZOZ8/src/Turing.jl:62 [inlined]
 [14] runmodel! at /Users/denshe/.julia/packages/Turing/RZOZ8/src/core/RandomVariables.jl:82 [inlined]
 [15] f at /Users/denshe/.julia/packages/Turing/RZOZ8/src/core/ad.jl:109 [inlined]
 [16] vector_mode_dual_eval(::getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}}, ::Array{Real,1}, ::ForwardDiff.GradientConfig{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4,Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1}}) at /Users/denshe/.julia/packages/ForwardDiff/N0wMF/src/apiutils.jl:37
 [17] vector_mode_gradient!(::Array{Real,1}, ::getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}}, ::Array{Real,1}, ::ForwardDiff.GradientConfig{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4,Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1}}) at /Users/denshe/.julia/packages/ForwardDiff/N0wMF/src/gradient.jl:103
 [18] gradient! at /Users/denshe/.julia/packages/ForwardDiff/N0wMF/src/gradient.jl:35 [inlined]
 [19] gradient!(::Array{Real,1}, ::getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}}, ::Array{Real,1}, ::ForwardDiff.GradientConfig{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4,Array{ForwardDiff.Dual{ForwardDiff.Tag{getfield(Turing.Core, Symbol("#f#26")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}},Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},Real},Real,4},1}}) at /Users/denshe/.julia/packages/ForwardDiff/N0wMF/src/gradient.jl:33
 [20] gradient_logp_forward(::Array{Real,1}, ::Turing.Core.RandomVariables.UntypedVarInfo, ::Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}}, ::Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}) at /Users/denshe/.julia/packages/Turing/RZOZ8/src/core/ad.jl:116
 [21] gradient_logp at /Users/denshe/.julia/packages/Turing/RZOZ8/src/core/ad.jl:80 [inlined]
 [22] #2 at /Users/denshe/.julia/packages/Turing/RZOZ8/src/inference/support/hmc_core.jl:12 [inlined]
 [23] #_leapfrog#26(::getfield(Turing.Inference, Symbol("##6#7")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}}, ::getfield(Turing.Inference, Symbol("##8#9")){Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}}, ::Function, ::Array{Real,1}, ::Array{Float64,1}, ::Int64, ::Float64, ::getfield(Turing.Inference, Symbol("##2#3")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}},Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}}}) at /Users/denshe/.julia/packages/Turing/RZOZ8/src/inference/support/hmc_core.jl:147
 [24] (::getfield(Turing.Inference, Symbol("#kw##_leapfrog")))(::NamedTuple{(:rev_func, :log_func),Tuple{getfield(Turing.Inference, Symbol("##6#7")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},getfield(Turing.Inference, Symbol("##8#9")){Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}}}}, ::typeof(Turing.Inference._leapfrog), ::Array{Real,1}, ::Array{Float64,1}, ::Int64, ::Float64, ::Function) at ./none:0
 [25] #_hmc_step#27(::Function, ::Function, ::Function, ::Array{Real,1}, ::Float64, ::getfield(Turing.Inference, Symbol("##4#5")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}},Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}}}, ::Function, ::getfield(Turing.Inference, Symbol("##16#17")), ::Int64, ::Float64, ::getfield(Turing.Inference, Symbol("##14#15")){Int64}) at /Users/denshe/.julia/packages/Turing/RZOZ8/src/inference/support/hmc_core.jl:195
 [26] (::getfield(Turing.Inference, Symbol("#kw##_hmc_step")))(::NamedTuple{(:rev_func, :log_func),Tuple{getfield(Turing.Inference, Symbol("##6#7")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},getfield(Turing.Inference, Symbol("##8#9")){Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}}}}, ::typeof(Turing.Inference._hmc_step), ::Array{Real,1}, ::Float64, ::Function, ::Function, ::Function, ::Int64, ::Float64, ::getfield(Turing.Inference, Symbol("##14#15")){Int64}) at ./none:0
 [27] (::getfield(Turing.Inference, Symbol("#kw##hmc_step")))(::NamedTuple{(:rev_func, :log_func),Tuple{getfield(Turing.Inference, Symbol("##6#7")){Turing.Core.RandomVariables.UntypedVarInfo,Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}},getfield(Turing.Inference, Symbol("##8#9")){Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}}}}, ::typeof(Turing.Inference.hmc_step), ::Array{Real,1}, ::Float64, ::Function, ::Function, ::Function, ::Float64, ::HMC{Turing.Core.ForwardDiffAD{40},Any}, ::getfield(Turing.Inference, Symbol("##14#15")){Int64}) at /Users/denshe/.julia/packages/Turing/RZOZ8/src/inference/hmc.jl:67
 [28] step(::Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}}, ::Turing.Sampler{HMC{Turing.Core.ForwardDiffAD{40},Any}}, ::Turing.Core.RandomVariables.UntypedVarInfo, ::Val{false}) at /Users/denshe/.julia/packages/Turing/RZOZ8/src/inference/hmc.jl:233
 [29] macro expansion at ./util.jl:213 [inlined]
 [30] #sample#37(::Bool, ::Nothing, ::Int64, ::Nothing, ::Function, ::Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}}, ::HMC{Turing.Core.ForwardDiffAD{40},Any}) at /Users/denshe/.julia/packages/Turing/RZOZ8/src/inference/hmc.jl:150
 [31] sample(::Turing.Model{Tuple{:alpha,:beta,:gamma,:delta},Tuple{:x},getfield(Main, Symbol("###inner_function#371#18")){Array{Int64,1}},NamedTuple{(:x,),Tuple{Array{Float64,2}}},NamedTuple{(:x,),Tuple{Symbol}}}, ::HMC{Turing.Core.ForwardDiffAD{40},Any}) at /Users/denshe/.julia/packages/Turing/RZOZ8/src/inference/hmc.jl:102
 [32] top-level scope at none:0

Julia version is 1.1.1

julia> Pkg.status()
    Status `~/.julia/environments/v1.1/Project.toml`
  [ebbdde9d] DiffEqBayes v1.1.0
  [aae7a2af] DiffEqFlux v0.5.0
  [41bf760c] DiffEqSensitivity v3.2.3
  [0c46a032] DifferentialEquations v6.4.0
  [31c24e10] Distributions v0.20.0
  [bbc10e6e] DynamicHMC v1.0.5
  [587475ba] Flux v0.8.3
  [f6369f11] ForwardDiff v0.10.3
  [7073ff75] IJulia v1.18.1
  [6fdf6af0] LogDensityProblems v0.8.2
  [65888b18] ParameterizedFunctions v4.1.1
  [d96e819e] Parameters v0.10.3
  [91a5bcdd] Plots v0.25.1
  [731186ca] RecursiveArrayTools v0.20.0
  [f3b207a7] StatsPlots v0.11.0
  [84d833dd] TransformVariables v0.3.3
  [fce5fe82] Turing v0.6.17

Please let me know if I can provide additional information.

[mohamed82008]

Perhaps an issue on a DiffEq package will get more attention from DiffEq people.

[xukai92]

Please also check DiffEqBayes.jl if you want to use DiffEq and Turing in a blackbox manner.. Here is a recent notebook https://github.com/JuliaDiffEq/DiffEqTutorials.jl/blob/master/notebook/models/pendulum_bayesian_inference.ipynb

[denisshepelin]

Please also check DiffEqBayes.jl if you want to use DiffEq and Turing in a blackbox manner.. Here is a recent notebook https://github.com/JuliaDiffEq/DiffEqTutorials.jl/blob/master/notebook/models/pendulum_bayesian_inference.ipynb

Looks like DIffEqBayes was updated to better fit current Turing architecture. I've tried their turing_inference function and it worked quite fine!