[SOLVED] Error when predicting partial effects using new data for gamlss model

Error when predicting partial effects using new data for gamlss model

I'm here re-raising the issue of predicting CI's for gamlss models using the newdata argument. A further complication is that I'm interested in partial effects as well.

A closely related issue (without partial effects) was un-resolved in 2018: Error when predicting new fitted values from R gamlss object.

I'm wondering if there has been updates that also extend to partial effects. The example below reproduces the error (notice the `type = "terms" specifying I'm interested in the effects of each model term)".

library(gamlss)
library(tidyverse)

#example data 
test_df <- tibble(x = rnorm(1e4),
                  x2 = rnorm(n = 1e4),
                  y = x2^2 + rnorm(1e4, sd = 0.5)) 

#fitting gamlss model
gam_test = gamlss(formula = y ~ pb(x2) + x,
                  sigma.fo= y ~ pb(x2) + x,
                  data = test_df)
#data I want predictions for
pred_df <- tibble(x = seq(-0.5, 0.5, length.out = 300),
                  x2 = seq(-0.5, 0.5, length.out = 300)) 

#returns error when se.fit = TRRUE
pred <- predictAll(object = gam_test,
                   type = "terms",
                   se.fit = TRUE, #works if se.fit = FALSE
                   newdata = pred_df)

Many thanks in advance!

Solution

I talked to the main developer of the gamlss software (who is responsible for this function). He says that the option se.fit=TRUE with type="terms" has not yet been implemented, and unfortunately he is too busy at present.

One idea is to bootstrap the original data, and predict terms for each bootstrap sample, and then use the results to obtain CI's.