TuringLang / Turing.jl

Bayesian inference with probabilistic programming.
https://turinglang.org
MIT License
2.01k stars 218 forks source link

@model treats an observed variable as a random variable if it's provided by unpacking #1978

Open GBarnsley opened 1 year ago

GBarnsley commented 1 year ago

When I provide values to a Turing model via unpacking from an array or tuple, they don't seem to be recognised and are treated as random variables. See the following code snippet:

using Turing, Optim
x = 3;
N = 10;

data = (x, N)

@model function multinom_model_tuple(data)
    x, N = data
    p ~ Uniform(0, 1)

    x ~ Binomial(N, p)
    #data[1] ~ Binomial(data[2], p) #this fixes the issue
end

@model function multinom_model(x, N)
    p ~ Uniform(0, 1)

    x ~ Binomial(N, p)
end

ppl_model_tuple = multinom_model_tuple(data)
ppl_model = multinom_model(x, N)

#Demonstration with MLE/MAP
mle_estimate_tuple = optimize(ppl_model_tuple, MLE());
mle_estimate = optimize(ppl_model, MLE());

mle_estimate.values
mle_estimate_tuple.values

#demonstration with sample
mh_estimate_tuple = sample(ppl_model_tuple, MH(), 1000);
mh_estimate = sample(ppl_model, MH(), 1000);

mh_estimate
mh_estimate_tuple

Is this just a limitation of the language?

It's no major issue since unpacking just helps to keep more complex models neat but it might be worth noting this in the documentation somewhere (if its not already).

torfjelde commented 1 year ago

That is indeed a restriction with the modelling language. We basically inspect the arguments of the model to determine if the thing on the left-hand-side of ~ should be interpreted as a observations or random variable.

One thing you can do is the following:

@model function multinom_model_tuple(data)
    x, N = data
    p ~ Uniform(0, 1)
    x ~ Binomial(N, p)

    return x
end

ppl_model_tuple = multinom_model_tuple() | (x = data[1],)

# Even make a constructor which automatically does this:
function make_model(::typeof(multinom_model_tuple), data)
    return multinom_model_tuple(data) | (x = data[1],)
end

ppl_model_tuple = make_model(multinom_model_tuple, data)

This is using the alternative way of conditioning (you can do it either by passing it as an argument, or you can do model | (x=1, y=2, ...)) as a way circumventing this restriction of the observations-by-passing-as-argument approach.

devmotion commented 1 year ago

Alternatively, you can also define the model with the unpacked arguments and define convenience functions that unpack for you. Something like

@model function multinom_model(x, N)
    p ~ Uniform(0, 1)
    x ~ Binomial(N, p)
    return x
end

multinom_model(data) = multinom_model(data.x, data.N)