# Multivariate Distributions & Parameter Inference

**URL:** https://discourse.julialang.org/t/multivariate-distributions-parameter-inference/57896
**Category:** New to Julia
**Created:** [March 24, 2021, 9:33pm UTC](https://discourse.julialang.org/t/multivariate-distributions-parameter-inference/57896 "2021-03-24T21:33:02Z")
**Posts on this page:** 2
**Page:** 1

<div class="post-metadata">

### Author: ![willsharpless](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/willsharpless/32/36397_2.png) [@willsharpless](https://discourse.julialang.org/u/willsharpless)
#### Post date: [March 24, 2021, 9:33pm UTC](https://discourse.julialang.org/t/multivariate-distributions-parameter-inference/57896/1 "2021-03-24T21:33:02Z")

</div>

Hi all,

I was trying to adapt the example from [this tutorial on parameter inference](https://turing.ml/dev/tutorials/10-bayesiandiffeq/#direct-handling-of-bayesian-estimation-with-turing)  
to a generalized Lotka-Volterra model with n equations when I ran into some issues. In order to write multivariate distributions into the Turing model (and avoid the much discussed ([Turing type error: expected Float64, got ForwardDiff.Dual](https://discourse.julialang.org/t/turing-type-error-expected-float64-got-forwarddiff-dual/47670))) I followed the instructions outlined at the [performance tips](https://turing.ml/dev/docs/using-turing/performancetips) note, however, in this workflow this pushed the error to the subsequent “predicted = …” line when evaluating with the sampled values because my ODE function is not fit to handle the ForwardDiff.Dual{…} type.

How should I restructure this to get around the type errors?

Here is the code:

```julia
##
using DifferentialEquations
using Plots
using Random
using Turing
# Turing.turnprogress(false)
## Generate Model and Synthetic Data
# Random.seed!(14);

function gLV!(dx,x,p,t)
    x[x .< 1e-8] .= 0
    mu, A, np = p
    dx .= x.*(mu + A*x)

    # no cc catch
    if any(x .> 100.0)
        dx .= 0
    end
end

function mult_noise!(dx,x,p,t)
    mu, A, np = p
    dx .= x.*np
end

mu = [0.4; 0.8; 0.2].*0.1
A = [-1.5 0.5 -0.2; -0.7 -1.3 0.1; 0.8 -0.2 -0.9]
np = [0.01; 0.01; 0.01] # species noise percentages

x0 = [0.1; 0.1; 0.1].*0.1
tspan = (0.0, 144.0)

# Definite Model with Read-Out Noise
np.=0

prob1 = ODEProblem(gLV!, x0, tspan, (mu, A, np))
sol = solve(prob1, Tsit5(), saveat = 24)
syndata = Array(sol) + 0.003*randn(size(Array(sol)))
plot(sol, alpha=0.3)
scatter!(sol.t, syndata')
title!("Definite")

## MCMC Fittin
Turing.setadbackend(:forwarddiff)

@model function fitglv(data, prob1, ::Type{T} = Float64) where {T}
    σ ~ InverseGamma(2, 3)

    n = size(prob1.u0,1)
    mu = Vector{T}(undef, n)
    A = Matrix{T}(undef, n, n)

    for i = 1:n
        mu[i] ~ truncated(Normal(0.4, 0.2), 0, 1.5)
    end

    for i = 1:n
        for j = 1:n
            if i == j
                A[i,j] ~ truncated(Normal(-0.8,1),-3,0)
            else
                A[i,j] ~ truncated(Normal(-0.3,1),-3,3)
            end
        end
    end

    np .= 0
    p = (mu, A, np)
    prob = remake(prob1, p=p)
    predicted = solve(prob, Tsit5(), saveat=24)

    for i = 1:length(predicted)
        data[:,i] ~ MvNormal(predicted[i], σ)
    end
end

model = fitglv(syndata, prob1)
chain = mapreduce(c -> sample(model, NUTS(.65), 1000), chainscat, 1:3)

```

and the full error I am gettting:

```julia
TypeError: in typeassert, expected Float64, got a value of type ForwardDiff.Dual{Nothing,Float64,7}
setindex!(::Array{Float64,1}, ::ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.VarInfo{NamedTuple{(:σ, :mu, :A),Tuple{DynamicPPL.Metadata{Dict{DynamicPPL.VarName{:σ,Tuple{}},Int64},Array{InverseGamma{Float64},1},Array{DynamicPPL.VarName{:σ,Tuple{}},1},Array{Float64,1},Array{Set{DynamicPPL.Selector},1}},DynamicPPL.Metadata{Dict{DynamicPPL.VarName{:mu,Tuple{Tuple{Int64}}},Int64},Array{Truncated{Normal{Float64},Continuous,Float64},1},Array{DynamicPPL.VarName{:mu,Tuple{Tuple{Int64}}},1},Array{Float64,1},Array{Set{DynamicPPL.Selector},1}},DynamicPPL.Metadata{Dict{DynamicPPL.VarName{:A,Tuple{Tuple{Int64,Int64}}},Int64},Array{Truncated{Normal{Float64},Continuous,Float64},1},Array{DynamicPPL.VarName{:A,Tuple{Tuple{Int64,Int64}}},1},Array{Float64,1},Array{Set{DynamicPPL.Selector},1}}}},Float64},DynamicPPL.Model{var"#103#104",(:data, :prob1, :T),(:T,),(),Tuple{Array{Float64,2},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Tuple{Array{Float64,1},Array{Float64,2},Array{Float64,1}},ODEFunction{true,typeof(gLV!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,typeof(SciMLBase.DEFAULT_OBSERVED),Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},SciMLBase.StandardODEProblem},Type{Float64}},Tuple{Type{Float64}}},DynamicPPL.Sampler{NUTS{Turing.Core.ForwardDiffAD{40},(),AdvancedHMC.DiagEuclideanMetric}},DynamicPPL.DefaultContext},Float64},Float64,7}, ::Int64) at array.jl:847
macro expansion at broadcast.jl:932 [inlined]
macro expansion at simdloop.jl:77 [inlined]
copyto! at broadcast.jl:931 [inlined]
copyto! at broadcast.jl:886 [inlined]
materialize! at broadcast.jl:848 [inlined]
materialize! at broadcast.jl:845 [inlined]
gLV!(::Array{Float64,1}, ::Array{Float64,1}, ::Tuple{Array{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.VarInfo{NamedTuple{(:σ, :mu, :A),Tuple{DynamicPPL.Metadata{Dict{DynamicPPL.VarName{:σ,Tuple{}},Int64},Array{InverseGamma{Float64},1},Array{DynamicPPL.VarName{:σ,Tuple{}},1},Array{Float64,1},Array{Set{DynamicPPL.Selector},1}},DynamicPPL.Metadata{Dict{DynamicPPL.VarName{:mu,Tuple{Tuple{Int64}}},Int64},Array{Truncated{Normal{Float64},Continuous,Float64},1},Array{DynamicPPL.VarName{:mu,Tuple{Tuple{Int64}}},1},Array{Float64,1},Array{Set{DynamicPPL.Selector},1}},DynamicPPL.Metadata{Dict{DynamicPPL.VarName{:A,Tuple{Tuple{Int64,Int64}}},Int64},Array{Truncated{Normal{Float64},Continuous,Float64},1},Array{DynamicPPL.VarName{:A,Tuple{Tuple{Int64,Int64}}},1},Array{Float64,1},Array{Set{DynamicPPL.Selector},1}}}},Float64},DynamicPPL.Model{var"#103#104",(:data, :prob1, :T),(:T,),(),Tuple{Array{Float64,2},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Tuple{Array{Float64,1},Array{Float64,2},Array{Float64,1}},ODEFunction{true,typeof(gLV!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,typeof(SciMLBase.DEFAULT_OBSERVED),Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},SciMLBase.StandardODEProblem},Type{Float64}},Tuple{Type{Float64}}},DynamicPPL.Sampler{NUTS{Turing.Core.ForwardDiffAD{40},(),AdvancedHMC.DiagEuclideanMetric}},DynamicPPL.DefaultContext},Float64},Float64,7},1},Array{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.VarInfo{NamedTuple{(:σ, :mu, :A),Tuple{DynamicPPL.Metadata{Dict{DynamicPPL.VarName{:σ,Tuple{}},Int64},Array{InverseGamma{Float64},1},Array{DynamicPPL.VarName{:σ,Tuple{}},1},Array{Float64,1},Array{Set{DynamicPPL.Selector},1}},DynamicPPL.Metadata{Dict{DynamicPPL.VarName{:mu,Tuple{Tuple{Int64}}},Int64},Array{Truncated{Normal{Float64},Continuous,Float64},1},Array{DynamicPPL.VarName{:mu,Tuple{Tuple{Int64}}},1},Array{Float64,1},Array{Set{DynamicPPL.Selector},1}},DynamicPPL.Metadata{Dict{DynamicPPL.VarName{:A,Tuple{Tuple{Int64,Int64}}},Int64},Array{Truncated{Normal{Float64},Continuous,Float64},1},Array{DynamicPPL.VarName{:A,Tuple{Tuple{Int64,Int64}}},1},Array{Float64,1},Array{Set{DynamicPPL.Selector},1}}}},Float64},DynamicPPL.Model{var"#103#104",(:data, :prob1, :T),(:T,),(),Tuple{Array{Float64,2},ODEProblem{Array{Float64,1},Tuple{Float64,Float64},true,Tuple{Array{Float64,1},Array{Float64,2},Array{Float64,1}},ODEFunction{true,typeof(gLV!),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,typeof(SciMLBase.DEFAULT_OBSERVED),Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},SciMLBase.StandardODEProblem},Type{Float64}},Tuple{Type{Float64}}},DynamicPPL.Sampler{NUTS{Turing.Core.ForwardDiffAD{40},(),AdvancedHMC.DiagEuclideanMetric}},DynamicPPL.DefaultContext},Float64},Float64,7},2},Array{Float64,1}}, ::Float64) at PostBin_Analysis_v1.1.jl:33
ODEFunction at scimlfunctions.jl:334 [inlined]
initialize!(::Ordinar...

```

 ![image](https://global.discourse-cdn.com/julialang/original/3X/0/7/075f6202f8992dd9cb2104931cee2699f3f96bdb.png)

---

<div class="post-metadata">

### Author: ![willsharpless](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/willsharpless/32/36397_2.png) [@willsharpless](https://discourse.julialang.org/u/willsharpless)
#### Post date: [March 24, 2021, 9:52pm UTC](https://discourse.julialang.org/t/multivariate-distributions-parameter-inference/57896/2 "2021-03-24T21:52:30Z")

</div>

Note, I tried using `filldist()` to see if it would work instead of the `{T}` fix, but my problem with this is that I need a different distribution for the diagonal of the matrix. My overwriting with CartesianIndex is not working:

 ![image](https://global.discourse-cdn.com/julialang/original/3X/a/6/a66b02feb31538a7f3266650cb107309683ef9f7.png)

Error:  
`MethodError: no method matching axes(::DistributionsAD.MatrixOfUnivariate{Continuous,Truncated{Normal{Float64},Continuous,Float64},FillArrays.Fill{Truncated{Normal{Float64},Continuous,Float64},2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}}}, ::Int64)`