Search code examples
predictrstanrstanarm

Predicting from the full posterior distribution using stan_glmer


Could I ask for some help please?

I have fit a binomial model using stan_glmer and have picked the model which I think best fits the data. I have used the posterior predict command to compare my observed data to data simulated by the model and it seems very similar.

I now want to predict the probability of an event for different levels of the predictors. I would usually use the predict command in glmer but I know I should use the posterior_predict command for stan_glmer to take into account the full uncertainty in the model. If x1 and x2 are continuous predictors for a binary event and I want a random intercept on group, the model formula would be:

      model <- stan_glmer(binary event ~ x1 + x2 +(1 | group), family="binomial"

My question is: I want to vary the predictors (x1 and x2) to see how the model predicts the observed data (and the variability in those predictions), maybe as a plot but I’m not sure how. Any help or guidance would be greatly appreciated.


Solution

  • In short, posterior_predict has a newdata argument that expects a data.frame with values of x1, x2, and group. This argument is similar to that in many other prediction functions and there is an example of using that can be executed via example(posterior_predict, package = "rstanarm").

    In your case, it might be something like nd <- with(original_data, expand.grid(x1 = seq(from = min(x1), to = max(x1), length.out = 20), x2 = seq(from = min(x2), to = max(x2), length.out = 20), group = levels(group))) PPD <- posterior_predict(model, newdata = nd) but you could choose the values of x1 and x2 in various other ways.