Turing: help with multinomial likelihood

Specific Domains Probabilistic programming
1

Hi
I’m running through a number of JAGS models, converting them to Turing. Mostly this is an easy process, but I’ve got stuck on what is probably something very simple.

Original JAGS model:

model{
  alpha ~ dbeta(1,1) 
  beta ~ dbeta(1,1)  
  gamma ~ dbeta(1,1)        
  pi[1] <- alpha*beta
  pi[2] <- alpha*(1-beta)
  pi[3] <- (1-alpha)*(1-gamma)
  pi[4] <- (1-alpha)*gamma   
  y[1:4] ~ dmulti(pi[],n)
}

where y = [14 4 5 210] and n = sum(y).

Attempt with Turing:

using Turing, StatsPlots

@model function model(y)
    n = sum(y)
    # priors
    α ~ Beta(1,1)
    β ~ Beta(1,1)
    γ ~ Beta(1,1)
    # rates
    p = Vector(undef, 4)
    p[1] = α * β
    p[2] = α * (1 - β)
    p[3] = (1 - α) * (1 - γ)
    p[4] = (1 - α) * γ
    # likelihood
    for i in 1:4
        y[i] ~ Multinomial(n, p[i])
    end
end

y = [14 4 5 210]
chain = sample(model(y), MH(), 1000)

I get the error ERROR: MethodError: no method matching Multinomial(::Int64, ::Float64)

I tried y ~ Multinomial(n, p) without the loop first, but no luck with that.

2

Try

@model function model(y, ::Type{T} = Float64) where {T}
    n = sum(y)
    # priors
    α ~ Beta(1,1)
    β ~ Beta(1,1)
    γ ~ Beta(1,1)
    # rates
    p = Vector{T}(undef, 4)
    p[1] = α * β
    p[2] = α * (1 - β)
    p[3] = (1 - α) * (1 - γ)
    p[4] = (1 - α) * γ
    # likelihood
    y ~ Multinomial(n, p)
end

y = [14, 4, 5, 210]
chain = sample(model(y), MH(), 1000)

You can do this as well:

@model function model(y)
    n = sum(y)
    # priors
    α ~ Beta(1,1)
    β ~ Beta(1,1)
    γ ~ Beta(1,1)
    # rates
    p = [α * β
         α * (1 - β)
         (1 - α) * (1 - γ)
         (1 - α) * γ]
    # likelihood
    y ~ Multinomial(n, p)
end

y = [14, 4, 5, 210]
chain = sample(model(y), MH(), 1000)
2 Likes
3

Perfect. That second model has a much cleaner look. Both worked, although the chains were bad with MH(). But works just great with HMC(0.05, 10) or NUTS(). Thanks very much @jangevaare

4

No problem. I don’t know off the top of my head if there’s performance implications of one specification over the other in the current version of Turing. Providing the type explicitly may be better. See the following link for more info:

https://turing.ml/dev/docs/using-turing/performancetips#ensure-that-types-in-your-model-can-be-inferred

1 Like