1

I'm trying to create a simple AR(2) process Y where:

Y[t] = ϵ[t] + b1 * y[t - 1] + b2 * y[t - 2]
ϵ[t] ~ Normal(0, σ)

and b1 and b2 are parameters drawn from some prior distributions.

The code is as follows:

using Statistics
using Turing, Distributions

# This is supposed to be AR(2) model
# Y[t] = ε[t] + b1 * Y[t-1] + b2 * Y[t-2]
@model function AR2(y, b1_mean, b2_mean, ::Type{T}=Float64; σ=.3, σ_b=.01) where T <: Real
    N = length(y)
    
    b1 ~ Normal(b1_mean, σ_b)
    b2 ~ Normal(b2_mean, σ_b)
    ϵ = Vector{T}(undef, N)
    
    y[1] ~ Normal(0, σ)
    y[2] ~ Normal(0, σ)
    
    for k ∈ 3:N
        ϵ[k] ~ Normal(0, σ)
        # Maybe I shouldn't be assigning here?
        y[k] = ϵ[k] + b1 * y[k - 1] + b2 * y[k - 2] # this is LINE 19
    end
    
    y
end

# These are the parameters
b1, b2 = .3, .6

println("Constructing random process with KNOWN parameters: b1=$b1, b2=$b2")
AR2_process = AR2(fill(missing, 500), b1, b2)()
AR2_process_orig = copy(AR2_process)

@show typeof(AR2_process)

println("Estimating parameters...")
model_ar = AR2(AR2_process, b1, b2)
chain_ar = sample(model_ar, HMC(0.001, 10), 200)

@show mean(chain_ar[:b1]), mean(chain_ar[:b2])
@show all(AR2_process .≈ AR2_process_orig)

println("Now try to estimate Float64")
model_ar = AR2(Float64.(AR2_process), b1, b2)
chain_ar = sample(model_ar, HMC(0.001, 10), 200)

@show mean(chain_ar[:b1]), mean(chain_ar[:b2])

The output is like this:

Constructing random process with KNOWN parameters: b1=0.3, b2=0.6
typeof(AR2_process) = Vector{Union{Missing, Float64}}
Estimating parameters...
Sampling 100%|█████████████████| Time: 0:00:47
(mean(chain_ar[:b1]), mean(chain_ar[:b2])) = (0.29404359023752596, 0.5912394358365198)
all(AR2_process .≈ AR2_process_orig) = true
Now try to estimate Float64
Sampling 100%|█████████████████| Time: 0:00:02
ERROR: LoadError: TypeError: in typeassert, expected Float64, got a value of type ForwardDiff.Dual{Nothing, Float64, 10}
Stacktrace:
  [1] setindex!(A::Vector{Float64}, x::ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}, i1::Int64)
    @ Base ./array.jl:839
  [2] #2
    @ ~/test/turing_question.jl:19 [inlined]
  [3] (::var"#2#3")(__rng__::Random._GLOBAL_RNG, __model__::DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, __varinfo__::DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}, Vector{Set{DynamicPPL.Selector}}}}}, ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}, __sampler__::DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, __context__::DynamicPPL.DefaultContext, y::Vector{Float64}, b1_mean::Float64, b2_mean::Float64, #unused#::Type{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}, σ::Float64, σ_b::Float64)
    @ Main ./none:0
  [4] macro expansion
    @ ~/.julia/packages/DynamicPPL/SgzCy/src/model.jl:0 [inlined]
  [5] _evaluate
    @ ~/.julia/packages/DynamicPPL/SgzCy/src/model.jl:154 [inlined]
  [6] evaluate_threadunsafe
    @ ~/.julia/packages/DynamicPPL/SgzCy/src/model.jl:127 [inlined]
  [7] (::DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}})(rng::Random._GLOBAL_RNG, varinfo::DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}, Vector{Set{DynamicPPL.Selector}}}}}, ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}, sampler::DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, context::DynamicPPL.DefaultContext)
    @ DynamicPPL ~/.julia/packages/DynamicPPL/SgzCy/src/model.jl:92
  [8] Model
    @ ~/.julia/packages/DynamicPPL/SgzCy/src/model.jl:98 [inlined]
  [9] f
    @ ~/.julia/packages/Turing/28kgo/src/core/ad.jl:111 [inlined]
 [10] chunk_mode_gradient!(result::Vector{Float64}, f::Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, x::Vector{Float64}, cfg::ForwardDiff.GradientConfig{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}})
    @ ForwardDiff ~/.julia/packages/ForwardDiff/QOqCN/src/gradient.jl:150
 [11] gradient!
    @ ~/.julia/packages/ForwardDiff/QOqCN/src/gradient.jl:39 [inlined]
 [12] gradient!(result::Vector{Float64}, f::Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, x::Vector{Float64}, cfg::ForwardDiff.GradientConfig{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#1"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 10}}})
    @ ForwardDiff ~/.julia/packages/ForwardDiff/QOqCN/src/gradient.jl:35
 [13] gradient_logp(::Turing.Core.ForwardDiffAD{40}, θ::Vector{Float64}, vi::DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, model::DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, sampler::DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, ctx::DynamicPPL.DefaultContext)
    @ Turing.Core ~/.julia/packages/Turing/28kgo/src/core/ad.jl:121
 [14] gradient_logp (repeats 2 times)
    @ ~/.julia/packages/Turing/28kgo/src/core/ad.jl:83 [inlined]
 [15] ∂logπ∂θ
    @ ~/.julia/packages/Turing/28kgo/src/inference/hmc.jl:433 [inlined]
 [16] ∂H∂θ
    @ ~/.julia/packages/AdvancedHMC/MIxdK/src/hamiltonian.jl:31 [inlined]
 [17] phasepoint(h::AdvancedHMC.Hamiltonian{AdvancedHMC.UnitEuclideanMetric{Float64, Tuple{Int64}}, Turing.Inference.var"#logπ#52"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}}, Turing.Inference.var"#∂logπ∂θ#51"{DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}}}, θ::Vector{Float64}, r::Vector{Float64})
    @ AdvancedHMC ~/.julia/packages/AdvancedHMC/MIxdK/src/hamiltonian.jl:69
 [18] phasepoint
    @ ~/.julia/packages/AdvancedHMC/MIxdK/src/hamiltonian.jl:139 [inlined]
 [19] initialstep(rng::Random._GLOBAL_RNG, model::DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, spl::DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, vi::DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}; init_params::Nothing, nadapts::Int64, kwargs::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ Turing.Inference ~/.julia/packages/Turing/28kgo/src/inference/hmc.jl:167
 [20] initialstep(rng::Random._GLOBAL_RNG, model::DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, spl::DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}, vi::DynamicPPL.TypedVarInfo{NamedTuple{(:b1, :b2, :ϵ), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b1, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b1, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b2, Tuple{}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:b2, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}, Int64}, Vector{Normal{Float64}}, Vector{AbstractPPL.VarName{:ϵ, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64})
    @ Turing.Inference ~/.julia/packages/Turing/28kgo/src/inference/hmc.jl:153
 [21] step(rng::Random._GLOBAL_RNG, model::DynamicPPL.Model{var"#2#3", (:y, :b1_mean, :b2_mean, :T, :σ, :σ_b), (:T, :σ, :σ_b), (), Tuple{Vector{Float64}, Float64, Float64, Type{Float64}, Float64, Float64}, Tuple{Type{Float64}, Float64, Float64}}, spl::DynamicPPL.Sampler{HMC{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.UnitEuclideanMetric}}; resume_from::Nothing, kwargs::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ DynamicPPL ~/.julia/packages/DynamicPPL/SgzCy/src/sampler.jl:91
 [22] step
    @ ~/.julia/packages/DynamicPPL/SgzCy/src/sampler.jl:66 [inlined]
 [23] macro expansion
    @ ~/.julia/packages/AbstractMCMC/ByHEr/src/sample.jl:123 [inlined]
...

(I had to cut the error message because the post exceeded the 30000 characters limit)

So it's saying that line 19:

y[k] = ϵ[k] + b1 * y[k - 1] + b2 * y[k - 2]

...attempted to insert "a value of type ForwardDiff.Dual{Nothing, Float64, 10}" into y, but y is of type Vector{Float64}, since I explicitly cast it to Float64 with Float64.(AR2_process), so the types didn't match which caused the error. This kinda makes sense.

Question

Why did the previous run work, though? Why did sampling from AR2(AR2_process, b1, b2), where AR2_process is a Vector{Union{Missing, Float64}}, not raise any errors? Wouldn't it attempt to assign ForwardDiff.Dual to y[k] here as well?

Potential solution

Perhaps I shouldn't be assigning stuff to y[k] because observations are supposed to be drawn from some distribution, so it should probably be y[k] ~ Normal(...) since I'm adding values drawn from normal distributions, but I'm not particularly sure what the parameters to that normal distribution should be...

ForceBru
  • 43,482
  • 10
  • 63
  • 98

1 Answers1

2

I don't know why your code fails the way it does, but I can explain and fix your real problem.

Yes, you must always use tildes for observations; but the crucial is that the right hand side of your observation statement is not a distribution, but deterministic. If you don't actually need to sample/estimate the epsilons, you can use the following model instead:

@model function AR2(y, b1_mean, b2_mean; σ=.3, σ_b=.01) where T <: Real
    N = length(y)
           
    b1 ~ Normal(b1_mean, σ_b)
    b2 ~ Normal(b2_mean, σ_b)
           
    y[1] ~ Normal(0, σ)
    y[2] ~ Normal(0, σ)
           
    for k ∈ 3:N
        y[k] ~ Normal(b1 * y[k - 1] + b2 * y[k - 2], σ) # this is LINE 19
    end
           
    y
end

This only works because we know the properties of Normal, though. The question of tracking determinstic functions of random variables is still an open issue.

Or, you try out the new hot stuff and write a submodel based implementation.

phipsgabler
  • 20,535
  • 4
  • 40
  • 60
  • Yeah, I ended up doing it this way too. However, I don't really like relying on properties of the normal distribution. What if `ϵ` had such a distribution that sums like `b1 * y[k - 1] + ... + ϵ` had a really complex non-standard distribution? Guess I'll wait for `@submodel` to land (it doesn't look like it got into the current v0.15.24, did it?) – ForceBru May 29 '21 at 15:33