Turing model fail due to extreme outliers sampled in Exponential distribution?

[amit1]

I am modelling some textbook models (minimal example given below)

using Turing

# Data ===================
contact = [1, 1, 1, 2, 2, 2, 2, 1, 2, 1]
total_tools = [ 13 ,22 ,24 ,43 ,33 ,19 ,40 ,28 ,55 ,71]
population = [1100, 1500, 3600, 4791, 7400, 8000, 9200, 13000, 17500, 275000]
population = (population .- mean(population))./std(population) # normalize
population = population .+ abs(minimum(population)) .+ 0.1    # shift to positive values


# model =============================
@model function m12_2(tools, population, cid)
	ϕ ~ Exponential(1.0)
	g ~ Exponential(1.0)
	a = Vector{Real}(undef,2)
	b = Vector{Real}(undef,2)
	a .~ Normal(1.0, 1.0)
	b .~ Exponential(1.0)
	for i = 1 : length(cid)
		λ = (exp(a[cid[i]]) * population[i]^ b[cid[i]])/g
		tools[i] ~ NegativeBinomial(λ/(ϕ + 0.000001), 1.0 / (1.0 + ϕ))
	end
end
describe(sample(m12_2(total_tools, population, contact), NUTS(),1000))

Problem is that above model fails randomly 2 out 3 times with error message:

ArgumentError: NegativeBinomial: the condition r > zero(r) is not satisfied.

    macro expansion@utils.jl:6[inlined]
    #NegativeBinomial#46@negativebinomial.jl:42[inlined]
    NegativeBinomial@negativebinomial.jl:41[inlined]
    #1@Other: 10[inlined]
    (::Main.workspace#160.var"#1#2")(::DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, ::DynamicPPL.ThreadSafeVarInfo{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#3"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, ...
...
    Model@model.jl:91[inlined]
    Model@model.jl:104[inlined]
    f@ad.jl:111[inlined]
    vector_mode_dual_eval@apiutils.jl:37[inlined]
    vector_mode_gradient!(::Vector{Float64}, ::Turing.Core.var"#f#3"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.DefaultContext}, ::Vector{Float64}, ::ForwardDiff.GradientConfig{ForwardDiff.Tag{Turing.Core.var"#f#3"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 6, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#3"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 6}}})@gradient.jl:113
    gradient!@gradient.jl:37[inlined]
    gradient!(::Vector{Float64}, ::Turing.Core.var"#f#3"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.DefaultContext}, ::Vector{Float64}, ::ForwardDiff.GradientConfig{ForwardDiff.Tag{Turing.Core.var"#f#3"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 6, Vector{ForwardDiff.Dual{ForwardDiff.Tag{Turing.Core.var"#f#3"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}, Float64, 6}}})@gradient.jl:35
    gradient_logp(::Turing.Core.ForwardDiffAD{40}, ::Vector{Float64}, ::DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, ::DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, ::DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, ::DynamicPPL.DefaultContext)@ad.jl:121
    gradient_logp@ad.jl:83[inlined]
    ∂logπ∂θ@hmc.jl:433[inlined]
    ∂H∂θ@hamiltonian.jl:31[inlined]
    macro expansion@UnPack.jl:100[inlined]
    var"#step#9"(::Bool, ::Val{false}, ::typeof(AdvancedHMC.step), ::AdvancedHMC.Leapfrog{Float64}, ::AdvancedHMC.Hamiltonian{AdvancedHMC.DiagEuclideanMetric{Float64, Vector{Float64}}, Turing.Inference.var"#logπ#52"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}}, Turing.Inference.var"#∂logπ∂θ#51"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}}}, ::AdvancedHMC.PhasePoint{Vector{Float64}, AdvancedHMC.DualValue{Float64, Vector{Float64}}}, ::Int64)@integrator.jl:88
    step@integrator.jl:66[inlined]
    build_tree(::Random._GLOBAL_RNG, ::AdvancedHMC.Trajectory{AdvancedHMC.MultinomialTS, AdvancedHMC.Leapfrog{Float64}, AdvancedHMC.GeneralisedNoUTurn{Float64}}, ::AdvancedHMC.Hamiltonian{AdvancedHMC.DiagEuclideanMetric{Float64, Vector{Float64}}, Turing.Inference.var"#logπ#52"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}}, Turing.Inference.var"#∂logπ∂θ#51"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}}}, ::AdvancedHMC.PhasePoint{Vector{Float64}, AdvancedHMC.DualValue{Float64, Vector{Float64}}}, ::AdvancedHMC.MultinomialTS{Float64}, ::Int64, ::Int64, ::Float64)@trajectory.jl:595
    transition(::Random._GLOBAL_RNG, ::AdvancedHMC.Trajectory{AdvancedHMC.MultinomialTS, AdvancedHMC.Leapfrog{Float64}, AdvancedHMC.GeneralisedNoUTurn{Float64}}, ::AdvancedHMC.Hamiltonian{AdvancedHMC.DiagEuclideanMetric{Float64, Vector{Float64}}, Turing.Inference.var"#logπ#52"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}}, Turing.Inference.var"#∂logπ∂θ#51"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}}}, ::AdvancedHMC.PhasePoint{Vector{Float64}, AdvancedHMC.DualValue{Float64, Vector{Float64}}})@trajectory.jl:641
    transition@sampler.jl:57[inlined]
    var"#step#42"(::Int64, ::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, ::typeof(AbstractMCMC.step), ::Random._GLOBAL_RNG, ::DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, ::DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, ::Turing.Inference.HMCState{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, AdvancedHMC.HMCKernel{AdvancedHMC.FullMomentumRefreshment, AdvancedHMC.Trajectory{AdvancedHMC.MultinomialTS, AdvancedHMC.Leapfrog{Float64}, AdvancedHMC.GeneralisedNoUTurn{Float64}}}, AdvancedHMC.Hamiltonian{AdvancedHMC.DiagEuclideanMetric{Float64, Vector{Float64}}, Turing.Inference.var"#logπ#52"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}}, Turing.Inference.var"#∂logπ∂θ#51"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}}}, AdvancedHMC.PhasePoint{Vector{Float64}, AdvancedHMC.DualValue{Float64, Vector{Float64}}}, AdvancedHMC.Adaptation.StanHMCAdaptor{AdvancedHMC.Adaptation.WelfordVar{Float64, Vector{Float64}}, AdvancedHMC.Adaptation.NesterovDualAveraging{Float64}}})@hmc.jl:237
    macro expansion@sample.jl:134[inlined]
    macro expansion@ProgressLogging.jl:328[inlined]
    macro expansion@logging.jl:8[inlined]
    var"#mcmcsample#20"(::Bool, ::String, ::Nothing, ::Int64, ::Int64, ::Type, ::Base.Iterators.Pairs{Symbol, Int64, Tuple{Symbol}, NamedTuple{(:nadapts,), Tuple{Int64}}}, ::typeof(AbstractMCMC.mcmcsample), ::Random._GLOBAL_RNG, ::DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, ::DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, ::Int64)@sample.jl:114
    var"#sample#40"(::Type, ::Nothing, ::Bool, ::Int64, ::Bool, ::Int64, ::Base.Iterators.Pairs{Union{}, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, ::typeof(StatsBase.sample), ::Random._GLOBAL_RNG, ::DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, ::DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, ::Int64)@hmc.jl:133
    sample@hmc.jl:116[inlined]
    #sample#2@Inference.jl:142[inlined]
    sample@Inference.jl:142[inlined]
    #sample#1@Inference.jl:132[inlined]
    sample(::DynamicPPL.Model{Main.workspace#160.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, ::Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}, ::Int64)@Inference.jl:132
    top-level scope@Local: 1[inlined]

To figure out the bug I used try ... catch blocks inside my model as

@model function m12_2(tools, population, cid)
	ϕ ~ Exponential(1.0)
	g ~ Exponential(1.0)
	a = Vector{Real}(undef,2)
	b = Vector{Real}(undef,2)
	a .~ Normal(1.0, 1.0)
	b .~ Exponential(1.0)
	for i = 1 : length(cid)
		λ = (exp(a[cid[i]]) * population[i]^ b[cid[i]])/g
		try
			tools[i] ~ NegativeBinomial(λ/(ϕ + 0.000001), 1.0 / (1.0 + ϕ))
		catch
			println(i)
			println(population[i], " ",b[cid[i]])
			return -1
		end
	end
end

It looks like the random variable b, drawn from Exponential(1.0) sample extremely large values for some reason (like 2.0947637423773094e13 ?), resulting in population[i]^b[cid] becoming numerically zero.
Given below is the output of above model.

      From worker 3:	1
      From worker 3:	0.1 Dual{ForwardDiff.Tag{Turing.Core.var"#f#3"{DynamicPPL.TypedVarInfo{NamedTuple{(:ϕ, :g, :a, :b), Tuple{DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:ϕ, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:ϕ, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:g, Tuple{}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:g, Tuple{}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Normal{Float64}}, Vector{AbstractPPL.VarName{:a, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}, DynamicPPL.Metadata{Dict{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}, Int64}, Vector{Distributions.Exponential{Float64}}, Vector{AbstractPPL.VarName{:b, Tuple{Tuple{Int64}}}}, Vector{Float64}, Vector{Set{DynamicPPL.Selector}}}}}, Float64}, DynamicPPL.Model{Main.workspace#173.var"#1#2", (:tools, :population, :cid), (), (), Tuple{Vector{Int64}, Vector{Float64}, Vector{Int64}}, Tuple{}}, DynamicPPL.Sampler{Turing.Inference.NUTS{Turing.Core.ForwardDiffAD{40}, (), AdvancedHMC.DiagEuclideanMetric}}, DynamicPPL.DefaultContext}, Float64}}(2.0947637423773094e13,0.0,0.0,0.0,0.0,2.0947637423773094e13,0.0)

1. Is it some bug or something wrong in my model?
2. Is there some sort of Turing debugger to assist in such matter?

[dlakelan]

Although there are plans to remove it eventually for now I think you should use filldist or arraydist to declare a and b. See the Turing docs

[amit1]

Is expected outcome of filldist different from .~ ? I tried it but the problem persist. Currently my only workaround is to use truncated(Exponential(1.0), 1, 10) .

[dlakelan]

Which sampling algorithm are you using?

[amit1]

NUTS. Though I tried SMC, and HMC as well. But the same error persisted.

[dlakelan]

So this is sampling as a parameter, if it’s heading off to infinity it might be because some other parameter collapses and you get divergent transitions. try NUTS(1000,0.99) to make it do a very careful solution to the HMC differential equations, this might help?

[dlakelan]

@amit1

λ/(ϕ + 0.000001)

you might also try making this fudge factor larger, like phi + 0.001

[amit1]

I made fudge factor = 0.001 and sampled using NUTS(1000,0.99) but the problem still persists.

[dlakelan]

weird. does this occur during the NUTS warmup or after?

[amit1]

I would assume its after warmup, as it only gives the error if i draw samples > 500 or so. At lower samples it rarely crashes if at all.

[dlakelan]

I suspect if phi goes to zero then lambda needs to go to zero to maintain a reasonable ratio. Try making the fudge factor something like 0.2 just to see if it fixes the issue.

[amit1]

No. Still the same with fudge factor of 0.2.
Only truncated(Exponential()) works as of now.

[dlakelan]

Does it do the same thing with RWMH algorithm?

Ok, so just to make sure there isn’t a basic problem with sampling I ran this:

using Turing, StatsPlots, MCMCChains

@model function f()
    x ~ Exponential(1.0)
end

ch = sample(f(),NUTS(1000,.8),2000)

And it sampled normally with 97.5%tile of 3.8621 under latest version of Turing

Then I copy pasta your original code and get the error

ERROR: LoadError: ArgumentError: NegativeBinomial: the condition r > zero(r) is not satisfied.
Stacktrace:

The fact that this is an intermittent bug suggests to me that the model itself is somehow unstable.

If I modify the code as follows:

		λ = (exp(a[cid[i]]) * (.00001 + population[i]^ b[cid[i]]))/g

The model samples just fine and doesn’t give strange values for b. I think this is a modeling issue rather than a sampling issue

[amit1]

When you would modify it as (.00001 + population[i]^ b[cid[i]])) issue will be mitigated as \lambda will no longer go to zero. Problem is that population[i]^ b[cid[i]] goes to zero if b goes higher. With addition of 0.0001, values of b will no longer reduce \lambda to zero.

To replicate it in simpler model i tried the following:

@model function f()
           x ~ Exponential(1.0)
           y = 0.1^x
           z ~ NegativeBinomial(1/y,0.1)
           end

But doing so does not reproduce the original error. I do get warnings though.

┌ Warning: The current proposal will be rejected due to numerical error(s).
│   isfinite.((θ, r, ℓπ, ℓκ)) = (true, true, false, true)
└ @ AdvancedHMC ~/.julia/packages/AdvancedHMC/HQHnm/src/hamiltonian.jl:47

So I think you are right this might be model issue. Will try looking in model more carefully and update in future if I notice something more.
Thank you

[dlakelan]

@amit1

With addition of 0.0001, values of b will no longer reduce \lambda λ\lambda to zero.

Yes but also in the samples that I did values of B never went particularly high the whole phenomenon disappeared. I think you’re right to look more carefully at the model and you will find the answer somehow in modifications to the model