Error for vglm triple interaction

[iago-pssjd]

I was just trying to reproduce the example in https://cran.r-project.org/web/packages/sjPlot/vignettes/plot_interactions.html, but with VGAM::vglm, which seems to be allowed for the ggpredict function (thanks to the existence of get_predictions_vglm), but I get the next error:

library(VGAM)
library(sjPlot)
library(sjmisc)
data(efc)
# make categorical
efc$c161sex <- to_factor(efc$c161sex)
# fit model with 3-way-interaction
fit <- lm(neg_c_7 ~ c12hour * barthtot * c161sex, data = efc)
# it is the same as
fit <- vglm(neg_c_7 ~ c12hour * barthtot * c161sex, family = uninormal, data = efc)
# select only levels 30, 50 and 70 from continuous variable Barthel-Index
plot_model(fit, type = "pred", terms = c("c12hour", "barthtot [30,50,70]", "c161sex"))
Error in 1:ncoly : argument of length 0
 traceback()
10: eta2theta(eta[, M1 * (1:ncoly) - 1], "identitylink", earg = structure(list(
        theta = , inverse = FALSE, deriv = 0, short = TRUE, tag = FALSE), function.name = "identitylink"))
9: linv(prdat$fitted.values)
8: withCallingHandlers(expr, warning = function(w) if (inherits(w, 
       classes)) tryInvokeRestart("muffleWarning"))
7: suppressWarnings(linv(prdat$fitted.values))
6: get_predictions_vglm(model, data_grid, ci.lvl, linv, ...)
5: select_prediction_method(model_class = model_class, model = model, 
       data_grid = data_grid, ci.lvl = ci.lvl, type = type, model_info = model_info, 
       ppd = ppd, terms = original_terms, value_adjustment = typical, 
       vcov.fun = vcov.fun, vcov.type = vcov.type, vcov.args = vcov.args, 
       condition = condition, interval = interval, ...)
4: ggpredict_helper(model = model, terms = terms, ci.lvl = ci.lvl, 
       type = type, typical = typical, ppd = ppd, condition = condition, 
       back.transform = back.transform, vcov.fun = vcov.fun, vcov.type = vcov.type, 
       vcov.args = vcov.args, interval = interval, ...)
3: ggeffects::ggpredict(model = model, terms = terms, ci.lvl = ci.lvl, 
       type = pred.type, ...)
2: plot_type_eff(type = type, model = model, terms = terms, ci.lvl = ci.lvl, 
       pred.type = pred.type, facets = grid, show.data = show.data, 
       jitter = jitter, geom.colors = colors, axis.title = axis.title, 
       title = title, legend.title = legend.title, axis.lim = axis.lim, 
       case = case, show.legend = show.legend, dot.size = dot.size, 
       line.size = line.size, ...)
1: plot_model(fit, type = "pred", terms = c("c12hour", "barthtot [30,50,70]", 
       "c161sex"))

[strengejacke]

The family you specified has a link-inverse function that requires an "extra" argument. I'm not sure which values/objects this argument expects. I'll need some time to find out more...

library(VGAM)
#> Loading required package: stats4
#> Loading required package: splines
library(sjPlot)
library(sjmisc)
data(efc)
# make categorical
efc$c161sex <- to_factor(efc$c161sex)
# fit model with 3-way-interaction
fit <- vglm(neg_c_7 ~ c12hour * barthtot * c161sex, family = uninormal, data = efc)

link_inv <- fit@family@linkinv
link_inv
#> function (eta, extra = NULL) 
#> {
#>     M1 <- extra$M1
#>     ncoly <- extra$ncoly
#>     if ("identitylink" == "explink") {
#>         if (any(eta[, M1 * (1:ncoly) - 1] <= 0)) {
#>             warning("turning some columns of 'eta' positive in @linkinv")
#>             for (ii in 1:ncoly) eta[, M1 * ii - 1] <- pmax(1e-05, 
#>                 eta[, M1 * ii - 1])
#>         }
#>     }
#>     eta2theta(eta[, M1 * (1:ncoly) - 1], "identitylink", earg = list(
#>         theta = , inverse = FALSE, deriv = 0, short = TRUE, tag = FALSE))
#> }
#> <bytecode: 0x000001e171294bd0>
#> <environment: 0x000001e170acaa00>

^{Created on 2022-07-31 by the reprex package (v2.0.1)}