Predicting with a TuringGLM model

[h-spiess]

Hey,

I can't figure out an easy way to use a Turing model produced using turing_model to predict on some hold-out data. In Turing one usually passes missing for the predicted values. However, this seems to not be possible in the turing_model step. What workflow do you recommend? I would like to avoid building a second equivalent model which uses the sampled parameter values as this somehow defeats the purpose of this package.

Thanks in advance,
Heiner

[penelopeysm]

hi @h-spiess I think this is not currently possible because turing_model returns something with y (the responses) as one of the model arguments, which means that it's always treated as an observation. Sorry! I think it could be fixed relatively easily though.

[h-spiess]

Hey. Thanks for your answer. I could help with that but I think the fix is not that trivial.

Sure, it "easy" to accept a Vector containing missings. However, it's not clear to me how the API should adapt to that. E.g. how should standardization computed from the training data be passed? Additionally, there has to be a way to extract the corresponding priors from the fitted model to initialize the prediction model correctly.

Have you thought about the necessary API changes before?

Best,
Heiner

[penelopeysm]

Apologies, my comment was very brief and unclear. I hardly intend to claim that it is trivial --- but I think at least a direct fix would hopefully not require super deep digging into the internals of Turing / DynamicPPL.

You are right that there are a few components. I think it depends on how deep you want a fix to be.

The general point is that a Turing model can't be modified after it's been created. However, you can associate with it a context which modifies how it is executed. For example doing model | (; a=2.0) would create a new context within which the value of a is 2.0 (or to be precise, a new context within which DynamicPPL.getconditioned(..., @varname(a)) would return 2.0).

Right now the models are defined like:

@model function(X, y, ...; std_y=std(y))
    params ~ priors
    µ = f(params, X)
    y ~ dist(µ)
end

the simplest way would be to change the model something like this:

@model function(X, y, ...; std_y=std(y))
    params ~ priors
    if DynamicPPL.hasconditioned(__context__, @varname(X_pred))
        µ_pred = f(params, DynamicPPL.getconditioned(__context__, @varname(X_pred)))
        y_pred ~ dist(µ_pred)
    else
        µ = f(params, X)
        y ~ dist(X)
    end
end

then I think one could do

m = turing_model(X, y, ...) # X and y are stored in the model, so std(y) strictly refers to training data
chain = sample(m, ...)
pred_m = m | (; X_pred=...)
predict(pred_m, chain) # should generate a new chain with y_pred in it

That's obviously a bit hacky and would need to be documented carefully, but on the grand scheme of things it's quite a small level of hackiness (compared to many other things in Turing).

A full solution (and actually, maybe this is better for Turing in general, so maybe this would be better to do) would involve making something like a PredictionContext that automated the above. It could either live in TuringGLM or it could live in DynamicPPL; I think keeping it in TuringGLM would be safer because here we know the variable names and the structure of the model, so there would be fewer edge cases to cover.