This question derived from this question Join residual by group in summary Forecast table in R where forecast performed using ets function, for each group(1 and 2) The only one and serious difference is that it works with one variable. if I have a lot of variables, i must perform forecast for all of them at once. Let's take example
df=structure(list(Variable = structure(c(1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L), .Label = c("x", "y"), class = "factor"), group = c(1L,
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L,
1L, 2L, 2L, 2L, 2L, 2L, 2L), year = c(1973L, 1974L, 1975L, 1976L,
1977L, 1978L, 1973L, 1974L, 1975L, 1976L, 1977L, 1978L, 1973L,
1974L, 1975L, 1976L, 1977L, 1978L, 1973L, 1974L, 1975L, 1976L,
1977L, 1978L), Jan = c(9007L, 7750L, 8162L, 7717L, 7792L, 7836L,
9007L, 7750L, 8162L, 7717L, 7792L, 7836L, 9007L, 7750L, 8162L,
7717L, 7792L, 7836L, 9007L, 7750L, 8162L, 7717L, 7792L, 7836L
), Feb = c(8106L, 6981L, 7306L, 7461L, 6957L, 6892L, 8106L, 6981L,
7306L, 7461L, 6957L, 6892L, 8106L, 6981L, 7306L, 7461L, 6957L,
6892L, 8106L, 6981L, 7306L, 7461L, 6957L, 6892L), Mar = c(8928L,
8038L, 8124L, 7767L, 7726L, 7791L, 8928L, 8038L, 8124L, 7767L,
7726L, 7791L, 8928L, 8038L, 8124L, 7767L, 7726L, 7791L, 8928L,
8038L, 8124L, 7767L, 7726L, 7791L), Apr = c(9137L, 8422L, 7870L,
7925L, 8106L, 8192L, 9137L, 8422L, 7870L, 7925L, 8106L, 8192L,
9137L, 8422L, 7870L, 7925L, 8106L, 8192L, 9137L, 8422L, 7870L,
7925L, 8106L, 8192L), May = c(10017L, 8714L, 9387L, 8623L, 8890L,
9115L, 10017L, 8714L, 9387L, 8623L, 8890L, 9115L, 10017L, 8714L,
9387L, 8623L, 8890L, 9115L, 10017L, 8714L, 9387L, 8623L, 8890L,
9115L), Jun = c(10826L, 9512L, 9556L, 8945L, 9299L, 9434L, 10826L,
9512L, 9556L, 8945L, 9299L, 9434L, 10826L, 9512L, 9556L, 8945L,
9299L, 9434L, 10826L, 9512L, 9556L, 8945L, 9299L, 9434L), Jul = c(11317L,
10120L, 10093L, 10078L, 10625L, 10484L, 11317L, 10120L, 10093L,
10078L, 10625L, 10484L, 11317L, 10120L, 10093L, 10078L, 10625L,
10484L, 11317L, 10120L, 10093L, 10078L, 10625L, 10484L), Aug = c(10744L,
9823L, 9620L, 9179L, 9302L, 9827L, 10744L, 9823L, 9620L, 9179L,
9302L, 9827L, 10744L, 9823L, 9620L, 9179L, 9302L, 9827L, 10744L,
9823L, 9620L, 9179L, 9302L, 9827L), Sep = c(9713L, 8743L, 8285L,
8037L, 8314L, 9110L, 9713L, 8743L, 8285L, 8037L, 8314L, 9110L,
9713L, 8743L, 8285L, 8037L, 8314L, 9110L, 9713L, 8743L, 8285L,
8037L, 8314L, 9110L), Oct = c(9938L, 9129L, 8466L, 8488L, 8850L,
9070L, 9938L, 9129L, 8466L, 8488L, 8850L, 9070L, 9938L, 9129L,
8466L, 8488L, 8850L, 9070L, 9938L, 9129L, 8466L, 8488L, 8850L,
9070L), Nov = c(9161L, 8710L, 8160L, 7874L, 8265L, 8633L, 9161L,
8710L, 8160L, 7874L, 8265L, 8633L, 9161L, 8710L, 8160L, 7874L,
8265L, 8633L, 9161L, 8710L, 8160L, 7874L, 8265L, 8633L), Dec = c(8927L,
8680L, 8034L, 8647L, 8796L, 9240L, 8927L, 8680L, 8034L, 8647L,
8796L, 9240L, 8927L, 8680L, 8034L, 8647L, 8796L, 9240L, 8927L,
8680L, 8034L, 8647L, 8796L, 9240L)), .Names = c("Variable", "group",
"year", "Jan", "Feb", "Mar", "Apr", "May", "Jun", "Jul", "Aug",
"Sep", "Oct", "Nov", "Dec"), class = "data.frame", row.names = c(NA,
-24L))
Variable columns has variable x and y, and each variable has groups 1 and 2 to perfrom forecast
load libraries
load_pkgs <- c("forecast", "zoo", "timetk", "tidyverse")
sapply(load_pkgs, function(x) suppressPackageStartupMessages(library(x, character.only = T)))
Step 1: Pre-processing
# perform split by group
ld <- split(df[, -1], df$group)
# Tidy-up the splits
library(lubridate)
ld <- lapply(ld, function(x) {
x %>%
gather(key, value, -year) %>%
unite(date, year, key, sep = "-") %>%
mutate(date = paste0(date, "-01")) %>%
mutate(date =ymd(date))
})
Then
# Transform time series to ts objects
ld <- lapply(ld, function(x) {
yr <- lubridate::year(min(x$date))
mth <- lubridate::month(min(x$date))
timetk::tk_ts(data = x, select = value, frequency = 12,
start = c(yr, mth))
})
Step 2: Train and forecast with ets
# helping function
make_df <- function(ts_obj) {
ts_df <- timetk::tk_tbl(preserve_index = TRUE, ts_obj) %>%
mutate(index = zoo::as.Date(x = .$index, frac = 0)) %>%
dplyr::rename(date = index)
return(ts_df)
}
The following function trains ets and forecasts the next 12 months; then, it prepares tables with the fitted and forecasting values:
lts <- lapply(ld, function(ts_obj) {
# train ets model and get fitted results
res_model <- ets(ts_obj, model = "ZZZ")
res_fit <- ts(as.numeric(res_model$fitted), start = start(ts_obj), frequency = 12)
# add extra metrics you may be interested in
model <- res_model[["method"]]
mse <- res_model[["mse"]]
# get forecasts for the next 12 months
res_fct <- forecast(res_model, h = 12)
res_fcst <- ts(res_fct$mean, start = end(ts_obj) + 1/12, frequency = 12)
# transform results to tbl
# for fitted output we keep the residuals and the 95% CI
res_fit_tbl <- make_df(res_fit) %>%
mutate(residuals = as.numeric(res_model[["residuals"]])) %>%
mutate(CI95_upper = value + 1.96*sqrt(res_model$sigma2),
CI95_lower = value - 1.96*sqrt(res_model$sigma2))
# the forecast output does not have residuals
res_fcst_tbl <- make_df(res_fcst)
return(list(res_fit_tbl = res_fit_tbl, res_fcst_tbl = res_fcst_tbl, model = model, mse = mse)) # don't forget to pass the extra metrics as output
})
Step 3: Bring together the fitted and forecasting outputs across different groups
# add groups back + other metrics of interest
lts_all <- lapply(names(lts), function(x, lts) {
output_fit <- lts[[x]][["res_fit_tbl"]] %>%
mutate(group = x,
model = lts[[x]][["model"]],
mse = lts[[x]][["mse"]])
output_fcst <- lts[[x]][["res_fcst_tbl"]] %>%
mutate(group = x)
return(list(output_fit=output_fit, output_fcst=output_fcst))
}, lts)
Then
# bring together the fitted respectively forecasting results
output_fit_all <- lapply(lts_all, function(x) x[[1]])
output_fit_all <- bind_rows(output_fit_all)
output_fcst_all <- lapply(lts_all, function(x) x[[2]])
output_fcst_all <- bind_rows(output_fcst_all)
How to do, that this code perform forecast for all variable, like reproducible example