Predict for communityPGLMM

[vmikk]

Dear PGLMM gurus,

Unfortunately, currently there is no predict(object, newdata = ...) method for communityPGLMM class, and it is only possible to extract the fitted values with fitted.communityPGLMM.

In particular, I'm very interested in how to predict for the zero-inflated binomial model fitted with INLA:

modd <- phyr::pglmm(
    cbind(PositiveSamples, TotalSamples - PositiveSamples) ~ X * F + (1|sp__) + (1|site) + (1|sp__@site),
    family = "zeroinflated.binomial",
    data = dat,
    tree = tree,
    bayes = TRUE
    )

where X is a continuous predictor and F is a factor with 2 levels. In order to inspect the model behavior, I would like to plot the effect of X on the proportion of PositiveSamples for each level of F. This is very similar to the effect displays from effects package or marginal_effects from brms:

N <- 200
dummydata <- rbind(
  data.frame(
    PositiveSamples = sample(x = 0:5, size = N/2, replace = T),
    TotalSamples = 12,
    X = rnorm(n = N/2, mean = 10, sd = 1),
    F = "F1"
  ),
  data.frame(
    PositiveSamples = sample(x = 6:10, size = N/2, replace = T),
    TotalSamples = 12,
    X = rnorm(n = N/2, mean = 10, sd = 1),
    F = "F2"
  ))

mod <- glm(
  cbind(PositiveSamples, TotalSamples - PositiveSamples) ~ X * F,
  data = dummydata,
  family = binomial)

plot(effects::allEffects(mod))

There is also no built-in predict method in INLA which would allow predicting for new data, e.g., in the same way as for lm or glm with predict(mod, newdata = expand.grid(X = 1:3, F = c("F1", "F2")). Therefore, it is not easy to construct such plots.

For prediction, I also tried to append additional data cases with missing (NA) response to the model input data, similar to the INLA tutorial, and call communityPGLMM.bayes directly:

communityPGLMM.bayes(
  formula = modd$formula,
  data = dat,            # main data + data for prediction
  family = "zeroinflated.binomial",
  sp = modd$sp,
  site = modd$site,
  random.effects = modd$random.effects,
  s2.init = modd$s2.init,
  B.init = modd$B.init,
  marginal.summ = modd$marginal.summ,
  verbose = TRUE,
  calc.DIC = FALSE,
  prior = "inla.default", prior_alpha = 1, prior_mu = 0.1
  )

but I've got an error "variable lengths differ (found for 'inla_effects[[1]]') ". Running phyr::pglmm on the data with NAs or inlabru:::predict.inla(modd$inla.model) also doesn't work.

Maybe you know how to predict using the result of communityPGLMM class?

I would be grateful for any advice you could provide. With best regards, Vladimir Mikryukov