# Problem with Turing.jl

**URL:** <https://discourse.julialang.org/t/problem-with-turing-jl/117413>\
**Category:** General Usage\
**Tags:** bug\
**Created:** [July 24, 2024, 11:38am UTC](https://discourse.julialang.org/t/problem-with-turing-jl/117413 "2024-07-24T11:38:16Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![Houda\_Yaq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/houda_yaq/32/50747_2.png) [@Houda\_Yaq](https://discourse.julialang.org/u/Houda_Yaq)\
**Post date:** [July 24, 2024, 11:38am UTC](https://discourse.julialang.org/t/problem-with-turing-jl/117413/1 "2024-07-24T11:38:17Z")

</div>

Hello everyone,  
I am trying to estimate parameters of a SDE system using Turing.jl package, but I keep running into this issue and I cannot figure out how to solve it.  
This is my code:

```julia
using Turing
using Plots
using StatsPlots
using RecursiveArrayTools
#Prepare the data
t_end = 100.0
n_points = length(observed_data)
dt = t_end / (n_points - 1)
t_eval = 0:dt:t_end
#Set up the SDE problem with initial epsilon values
initial_ε = [0.1, 0.5, 0.4, 0.7, 0.2]  

noise_dim = 2*(K+L)

noise_prototype = zeros(noise_dim, noise_dim)

prob = SDEProblem(drift!, diffusion!, u0, (0,100), initial_ε, 
noise_rate_prototype=noise_prototype)

@model function fit_sde(data, prob)
    #Priors for epsilon parameters
    alpha ~ Uniform(0.0, 1.0)
    beta ~ filldist(Uniform(0.0, 1.0), K) # Assuming K epsilon values
    #Construct p 
    p = vcat(alpha, beta)
    
    #Sample from the prior
    _prob = remake(prob, p=p)

    #Solve the SDE ensemble
    
    predicted = solve(_prob, ISSEulerHeun(), tstops=t; p=p, saveat=0.1)
    #Likelihood (assuming Gaussian noise on observations)
    for j in eachindex(t)
        data[:, j] ~ MvNormal(predicted[j], Diagonal(fill(0.01^2, 12)))
    end
end
#Create the model instance
model = fit_sde(observed_data, prob)
#Sample from the posterior distribution (using NUTS sampler)
chain = sample(model, NUTS(0.65), 1000)

```

and this is the error I keep getting:

```julia
MethodError: no method matching Float64(::ForwardDiff.Dual{ForwardDiff.Tag{DynamicPPL.DynamicPPLTag, Float64}, Float64, 5})

Closest candidates are:
  (::Type{T})(::Real, ::RoundingMode) where T<:AbstractFloat
   @ Base rounding.jl:207
  (::Type{T})(::T) where T<:Number
   @ Core boot.jl:792
  Float64(::IrrationalConstants.Sqrt2π)
   @ IrrationalConstants ~/.julia/packages/IrrationalConstants/vp5v4/src/macro.jl:112
  ...

Stacktrace:
  [1] convert(::Type{Float64}, x::ForwardDiff.Dual{ForwardDiff.Tag{DynamicPPL.DynamicPPLTag, Float64}, Float64, 5})
    @ Base ./number.jl:7
  [2] setindex!
    @ ./array.jl:1024 [inlined]
  [3] macro expansion
    @ ./multidimensional.jl:960 [inlined]
  [4] macro expansion
    @ ./cartesian.jl:64 [inlined]
  [5] _unsafe_setindex!(::IndexLinear, ::Matrix{Float64}, ::Diagonal{ForwardDiff.Dual{ForwardDiff.Tag{DynamicPPL.DynamicPPLTag, Float64}, Float64, 5}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{DynamicPPL.DynamicPPLTag, Float64}, Float64, 5}}}, ::UnitRange{Int64}, ::UnitRange{Int64})
    @ Base ./multidimensional.jl:955
  [6] _setindex!
    @ ./multidimensional.jl:944 [inlined]
  [7] setindex!(::Matrix{Float64}, ::Diagonal{ForwardDiff.Dual{ForwardDiff.Tag{DynamicPPL.DynamicPPLTag, Float64}, Float64, 5}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{DynamicPPL.DynamicPPLTag, Float64}, Float64, 5}}}, ::UnitRange{Int64}, ::UnitRange{Int64})
    @ Base ./abstractarray.jl:1396
  [8] diffusion!(dW::Matrix{Float64}, u::Vector{ForwardDiff.Dual{ForwardDiff.Tag{DynamicPPL.DynamicPPLTag, Float64}, Float64, 5}}, p::Vector{ForwardDiff.Dual{ForwardDiff.Tag{DynamicPPL.DynamicPPLTag, Float64}, Float64, 5}}, t::Float64)
    @ Main ./In[195]:485

```

Does anyone know how to fix this? I would really appreciate it.  
Thank you!

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [July 24, 2024, 11:45am UTC](https://discourse.julialang.org/t/problem-with-turing-jl/117413/2 "2024-07-24T11:45:53Z")

</div>

How is your drift and diffusion defined?

---

<div class="post-metadata">

**Author:** ![Houda\_Yaq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/houda_yaq/32/50747_2.png) [@Houda\_Yaq](https://discourse.julialang.org/u/Houda_Yaq)\
**Post date:** [July 24, 2024, 11:49am UTC](https://discourse.julialang.org/t/problem-with-turing-jl/117413/3 "2024-07-24T11:49:03Z")

</div>

> [@Houda\_Yaq](#):
>
> `MethodError: no method matching Float64(::ForwardDiff.Dual`

The drift and the diffusion terms are a bit complxe and defined as follows:

```julia
function drift!(du, u, p, t) 
    alpha = p[1]
    beta = p[2:(K+1)]
    
    du[1:K] = mu_S(K, ch, N, u[(K+1):(2K)], u[1:K], u[(2K+L+1):(2K+2L)], alpha, NH, epsilon, a, omega) 
    du[(K+1):(2K)] = mu_I(K, ch, N, u[(K+1):(2K)], u[1:K], u[(2K+L+1):(2K+2L)], alpha, beta, NH, epsilon, a, omega)   
    du[(2K+1):(2K+L)] = mu_S_H(u[(2K+1):(2K+L)],u[(2K+L+1):(2K+2L)], N, u[1:K], u[(K+1):(2K)], alpha, H, L, A, epsilon, a, omega)
    du[(2K+L+1):(2K+2L)] = mu_I_H(u[(2K+L+1):(2K+2L)], u[(2K+1):(2K+L)], N, u[1:K], u[(K+1):(2K)], alpha, H, L, A, nu, epsilon, a, omega)

end

function diffusion!(dW, u, p, t)
    
    dW[1:K, 1:K] = sigma_SS(K, ch, N, u[(K+1):(2K)], u[1:K], u[(2K+L+1):(2K+2L)], alpha, NH, epsilon, a, omega)
    dW[1:K, (K+1):(2K)] = sigma_SI(K, ch, N, u[(K+1):(2K)], u[1:K], u[(2K+L+1):(2K+2L)], H, alpha, NH, epsilon, a, omega)
    dW[1:K, (2K+1):(2K+L)] = zeros(K,L)
    dW[1:K, (2K+L+1):(2K+2L)] = zeros(K,L)

    dW[(K+1):(2K), 1:K] = sigma_II(K, u[(K+1):(2K)], beta, N, NH)
    dW[(K+1):(2K), (K+1):(2K)] = zeros(K, K)
    dW[(K+1):(2K), (2K+1):(2K+L)] = zeros(K,L)
    dW[(K+1):(2K), (2K+L+1):(2K+2L)] = zeros(K,L)

    
    dW[(2K+1):(2K+L), 1:K] = zeros(L,K)
    dW[(2K+1):(2K+L), (K+1):(2K)] = zeros(L,K)
    dW[(2K+1):(2K+L), (2K+1):(2K+L)] = sigma_SS_HH(u[(2K+1):(2K+L)],u[(2K+L+1):(2K+2L)], N, u[1:K], u[(K+1):(2K)], alpha, H, L, A, epsilon, a, omega)
    dW[(2K+1):(2K+L), (2K+L+1):(2K+2L)] = sigma_SI_HH(u[(2K+1):(2K+L)],u[(2K+L+1):(2K+2L)], N, u[1:K], u[(K+1):(2K)], alpha, H, L, A, epsilon, a, omega)

    dW[(2K+L+1):(2K+2L), 1:K] = zeros(L,K)
    dW[(2K+L+1):(2K+2L), (K+1):(2K)] = zeros(L,K)
    dW[(2K+L+1):(2K+2L), (2K+1):(2K+L)] = sigma_II_HH(u[(2K+L+1):(2K+2L)], nu, H)
    dW[(2K+L+1):(2K+2L), (2K+L+1):(2K+2L)] = zeros(L, L)

end
```

---

<div class="post-metadata">

**Author:** ![Houda\_Yaq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/houda_yaq/32/50747_2.png) [@Houda\_Yaq](https://discourse.julialang.org/u/Houda_Yaq)\
**Post date:** [July 24, 2024, 11:55am UTC](https://discourse.julialang.org/t/problem-with-turing-jl/117413/4 "2024-07-24T11:55:35Z")

</div>

@ChrisRackauckas I think the problem is in the diffusion terms, all of the terms are the sqre root of a certain function, this is. an example

```julia
function sigma_II_tilde(j, i, beta, N, NH)
   
 # Compute the inside of the square root
    inside_sqrt = i[j] * beta[j]

  # Return the square root
  return(sqrt(inside_sqrt))
end

```

---

<div class="post-metadata">

**Author:** ![nilshg](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nilshg/32/2283_2.png) [@nilshg](https://discourse.julialang.org/u/nilshg)\
**Post date:** [July 24, 2024, 12:08pm UTC](https://discourse.julialang.org/t/problem-with-turing-jl/117413/5 "2024-07-24T12:08:04Z")

</div>

I haven’t looked at your code in detail, and as far as I can see it’s calling lots of functions which aren’t defined in what you’ve posted, but the error combined with your use of untyped `zeros` makes me think that you are using `Float64` arrays somewhere where you need to be more flexible on types, consider:

```julia
julia> using ForwardDiff

julia> f1(x) = x^2
f1 (generic function with 1 method)

julia> function f2(x)
           out = zeros(1)
           out[1] = x
           return only(x)^2
       end
f2 (generic function with 1 method)

julia> function f3(x)
           out = zeros(typeof(x), 1)
           out[1] = x
           return only(x)^2
       end
f3 (generic function with 1 method)

julia> ForwardDiff.derivative(f1, 5)
10

julia> ForwardDiff.derivative(f2, 5)
ERROR: MethodError: no method matching Float64(::ForwardDiff.Dual{ForwardDiff.Tag{typeof(f2), Int64}, Int64, 1})

julia> ForwardDiff.derivative(f3, 5)
10

```

---

<div class="post-metadata">

**Author:** ![fafafrens](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fafafrens/32/210932_2.png) [@fafafrens](https://discourse.julialang.org/u/fafafrens)\
**Post date:** [July 27, 2024, 9:11am UTC](https://discourse.julialang.org/t/problem-with-turing-jl/117413/6 "2024-07-27T09:11:30Z")

</div>

Hi as the previously answerd  
In diffusion! You use zeros(K,L) that is a float I would replace with something like

```julia
....
dW[1:K, (2K+1):(2K+L)] = zeros(eltype(u),K,L)
...

```