Difficulty using Submodels in Turing

Hi, I am trying to use Submodels as described in Submodels – Turing.jl. However I cannot figure out how to get it to work. I posted a deliberately simplified reproducable example below, and the error I was getting below it. Can someone explain how to use submodels? Thanks!

using Turing, Distributions

@model function my_submodel(n)
    sigma_rw ~ truncated(Cauchy(0, 1); lower=0)
    alpha = zeros(n)
    alpha[1] ~ Normal(0, sigma_rw)
    for i in 2:n
        alpha[i] ~ Normal(alpha[i-1], sigma_rw)
    end
    return alpha
end

@model function my_model(y)
    n = length(y)
    sigma_err ~ truncated(Cauchy(0, 1); lower=0)
    rw ~ to_submodel(my_submodel(n))
    for i in eachindex(y)
        y[i] ~ Normal(rw.alpha[i], sigma_err)
    end
end

n = 50
sub = my_submodel(50)
rw = sub()

model = my_model(rw)
chain = sample(model, NUTS(), 1_000)

Error when calling sampling:

julia> chain = sample(model, NUTS(), 1_000)
ERROR: type Array has no field alpha
Stacktrace:
  [1] getproperty
    @ ./Base.jl:49 [inlined]
  [2] macro expansion
    @ ~/.julia/packages/DynamicPPL/mOOQl/src/compiler.jl:566 [inlined]
  [3] my_model
    @ ~/Desktop/git/dssi-decsci-clinical-supply/dev/splines/tmp.jl:17 [inlined]
  [4] _evaluate!!
    @ ~/.julia/packages/DynamicPPL/mOOQl/src/model.jl:974 [inlined]
  [5] evaluate_threadunsafe!!(model::DynamicPPL.Model{…}, varinfo::DynamicPPL.VarInfo{…})
    @ DynamicPPL ~/.julia/packages/DynamicPPL/mOOQl/src/model.jl:940
  [6] check_model_and_trace(model::DynamicPPL.Model{…}, varinfo::DynamicPPL.VarInfo{…}; error_on_failure::Bool)
    @ DynamicPPL.DebugUtils ~/.julia/packages/DynamicPPL/mOOQl/src/debug_utils.jl:428
  [7] check_model_and_trace
    @ ~/.julia/packages/DynamicPPL/mOOQl/src/debug_utils.jl:418 [inlined]
  [8] check_model
    @ ~/.julia/packages/DynamicPPL/mOOQl/src/debug_utils.jl:451 [inlined]
  [9] _check_model
    @ ~/.julia/packages/Turing/ObWSF/src/mcmc/abstractmcmc.jl:5 [inlined]
 [10] _check_model
    @ ~/.julia/packages/Turing/ObWSF/src/mcmc/abstractmcmc.jl:8 [inlined]
 [11] sample(rng::Random.TaskLocalRNG, model::DynamicPPL.Model{…}, sampler::NUTS{…}, N::Int64; check_model::Bool, chain_type::Type, initial_params::InitFromUniform{…}, initial_state::Nothing, progress::Bool, nadapts::Int64, discard_adapt::Bool, discard_initial::Int64, kwargs::@Kwargs{})
    @ Turing.Inference ~/.julia/packages/Turing/ObWSF/src/mcmc/hmc.jl:101
 [12] sample
    @ ~/.julia/packages/Turing/ObWSF/src/mcmc/hmc.jl:86 [inlined]
 [13] #sample#1
    @ ~/.julia/packages/Turing/ObWSF/src/mcmc/abstractmcmc.jl:71 [inlined]
 [14] sample(model::DynamicPPL.Model{…}, spl::NUTS{…}, N::Int64)
    @ Turing.Inference ~/.julia/packages/Turing/ObWSF/src/mcmc/abstractmcmc.jl:68
 [15] top-level scope
    @ ~/Desktop/git/dssi-decsci-clinical-supply/dev/splines/tmp.jl:27
Some type information was truncated. Use `show(err)` to see complete types.

I don’t know much about Turing but it seems your submodel function returns alpha directly, not an object with a field named alpha. Does it work when you replace rw.alpha with rw?

1 Like

gdalle is right. try accesing rw and not rw.alpha, because you do return alpha. That should work:

@model function my_model(y)
    n = length(y)
    sigma_err ~ truncated(Cauchy(0, 1); lower=0)
    rw ~ to_submodel(my_submodel(n))
    for i in eachindex(y)
        y[i] ~ Normal(rw[i], sigma_err)
    end
end

Also have a look at the examples here. If you want to access model parameters with dot notation you should return a named tuple, i.e., (; alpha=alpha)

1 Like

For now the rw in rw ~ to_submodel(...) refers to the return value only. Maybe one point of confusion is that, by default, the names of the submodel latent variables are like rw.alpha[3] within the outer model scope (prefixed with the name given to the return value). Using this name, you could condition on the submodel latent variables like my_model() | (@varname(rw.alpha[3]) => 1), and so on.

1 Like

Thanks all! Those suggestions work. My specific point of confusion is from the Submodel docs, in the first code block it says

@model function inner()
    a ~ Normal()
    return a + 100
end

@model function outer()
    # This line adds the variable `x.a` to the chain.
    # The inner variable `a` is prefixed with the
    # left-hand side of the `~` operator, i.e. `x`.
    x ~ to_submodel(inner())
    # Here, the value of x will be `a + 100` because
    # that is the return value of the submodel.
    b ~ Normal(x)
end

It sounded to me like there was some behind the scenes magic going on to make a “namespace” x with all the random variables defined using ~ in the submodel on the right hand side of the ~ which I could access with dot notation, regardless of what I actually returned from the submodel. Later on in the docs it says

The reason for this is because it is entirely coincidental that the return value of the submodel is equal to a. In general, a return value can be anything, and conditioning on it is in general not a meaningful operation.

Which I understood to mean that one should avoid relying on the specific return values of a submodel and use this “namespace” dot method. But I see that I somehow still have a misunderstanding of what the docs mean.