Write a more general forecast decomposition function

[nicholasjclark]

plot_mvgam_uncertainty was originally designed to handle models with correlated residuals. It doesn't behave well for State-Space models (where predictors may be used in the process model) and it certainly will fail with dimension-reduced models (i.e. shared latent states). A more general forecast uncertainty decomposition function would be very valuable, but this will take some thought to design and implement.

Something to also consider while doing this is a proportion of explained / prediction variance:

# Attempt at a proportion of prediction variance
# using term predictions
dat <- mgcv::gamSim(1, n = 30, scale = 2)

fit1 <- mvgam(y ~ x0 +
                s(x1, k = 9, bs = 'cr') +
                s(x2, k = 22, bs = 'cr') +
                ti(x1, x2, k = 4, bs = 'cr'), 
              data = dat, 
              family = gaussian(),
              chains = 2)

tpreds <- predict(fit1, type = 'terms', summary = FALSE)

do.call(rbind, 
        lapply(seq_len(NROW(tpreds$obs_effects[[1]])),
               function(x){
                 data.frame(
                   do.call(cbind, 
                           lapply(tpreds$obs_effects, `[`,x,))) %>%
                   dplyr::mutate_all(.funs = var) %>%
                   dplyr::mutate(totvar = rowSums(.)) %>%
                   dplyr::summarise_all(.funs = ~ . / totvar) %>%
                   dplyr::distinct() %>%
                   dplyr::select(-totvar)
               }
        )) -> predvars
colnames(predvars) <- names(tpreds$obs_effects)
predvars %>%
  dplyr::summarise_all(.funs = quantile, 
                       probs = c(0.025, 0.5, 0.975)) %>%
  round(., 3)