对于模型时间中的递归面板数据，提取置信区间和残差失败

Posted 2023-02-14

【中文标题】对于模型时间中的递归面板数据，提取置信区间和残差失败

【英文标题】：Extracting confidence intervals and residuals is failing for recursive panel data in modeltime

【发布时间】：2022-01-09 17:11:41

【问题描述】：

在我的前一个问题中，我确实在模型时间中提取了递归集成的置信区间。但是，如果我将数据输入更改为时间序列面板，则预测（事前）对象中生成的置信区间列仅显示 NA。

在工作流程过程中肯定有一些错误，我猜到目前为止我没有发现。也许你们中的某个人可以提供帮助。

# Time Series ML
library(tidymodels)
library(modeltime)
library(modeltime.ensemble)

# Core
library(tidyverse)
library(timetk)


FORECAST_HORIZON <- 24

data = m750    # -----------> Single time series input works fine!

# changing the input to panel results in NA for residuals conf_hi and conf_lo...
data <- m4_monthly %>%
  group_by(id) %>%
  future_frame(
    .length_out = FORECAST_HORIZON,
    .bind_data  = TRUE
  ) %>%
  ungroup()

# data def
lag_transformer_grouped <- function(data)
  data %>%
    group_by(id) %>%
    tk_augment_lags(value, .lags = 1:FORECAST_HORIZON) %>%
    ungroup()


data_lags <- data %>%
  lag_transformer_grouped()

test_data <- data_lags %>%
  group_by(id) %>%
  slice_tail(n = 12) %>%
  ungroup()

train_data <- data_lags %>%
  drop_na()

future_data <- data_lags %>%
  filter(is.na(value))

# rec
recipe_spec <- recipe(value ~ date, train_data) %>%
  step_timeseries_signature(date) %>%
  step_rm(matches("(.iso$)|(.xts$)")) %>%
  step_normalize(matches("(index.num$)|(_year$)")) %>%
  step_dummy(all_nominal()) %>%
  step_fourier(date, K = 1, period = 12)

recipe_spec %>% prep() %>% juice()

# elnet 
model_fit_glmnet <- linear_reg(penalty = 1) %>%
  set_engine("glmnet") 

wflw_fit_glmnet <- workflow() %>%
  add_model(model_fit_glmnet) %>%
  add_recipe(recipe_spec %>% step_rm(date)) %>%
  fit(train_data)  %>%
  recursive(
    transform  = lag_transformer_grouped,
    train_tail = panel_tail(train_data, id, FORECAST_HORIZON),
    id         = "id"
  )

# xgboost    
model_fit_xgboost <- boost_tree("regression", learn_rate = 0.35) %>%
  set_engine("xgboost")

wflw_fit_xgboost <- workflow() %>%
  add_model(model_fit_xgboost) %>%
  add_recipe(recipe_spec %>% step_rm(date)) %>%
  fit(train_data) %>%
  recursive(
    transform  = lag_transformer_grouped,
    train_tail = panel_tail(train_data, id, FORECAST_HORIZON),
    id         = "id"
  )

# mtbl
data_models <- modeltime_table(
  wflw_fit_xgboost,
  wflw_fit_glmnet
)

# ensemble
ensemble_fit <- data_models %>% 
  ensemble_weighted(loadings = c(4, 6)) 

# doing calibration for panel with id, as indicated in the documentation 
calibration_tbl <- modeltime_table(
 ensemble_fit) %>%
 modeltime_calibrate(test_data, id = "id", quiet = FALSE) 

# but no residuals for me :(
calib_out <- calibration_tbl$.calibration_data[[1]] %>% 
  left_join(test_data %>% select(id, date, value))

# this function is explicitly for residual extraction. But without panel capability
residuals_method2 <- data_models %>%
  modeltime_calibrate(new_data = test_data) %>%
  modeltime_residuals() 

# Forecast ex ante
data_prepared_tbl <- bind_rows(train_data, test_data)

future_tbl <- data_prepared_tbl %>%
  group_by(id) %>%
  future_frame(.length_out = "2 years") %>%
  ungroup()

# documentation indicates, that calculation of conf_lo and conf_hi should work
ex_ante_obj <- 
  calibration_tbl %>%
  modeltime_forecast(
    new_data    = future_tbl,
    actual_data = data,
    keep_data = TRUE,
    conf_interval = 0.95,
    conf_by_id = TRUE
  )

【参考方案1】：

幸运的是，我通过 github 获得了一些帮助。再次感谢马特丹乔！

我将在下面演示正确的工作流程。这只是我使用的面板数据的错误定义，但你自己看看：

# Time Series ML
library(tidymodels)
library(modeltime)
library(modeltime.ensemble)

# Core
library(tidyverse)
library(timetk)


FORECAST_HORIZON <- 24

# changing the input to panel results in NA for conf_hi and conf_lo...
data <- m4_monthly %>%
  group_by(id) %>%
  ungroup()

# data def
lag_transformer_grouped <- function(data)
  data %>%
    group_by(id) %>%
    tk_augment_lags(value, .lags = 1:FORECAST_HORIZON) %>%
    ungroup()


data_lags <- data %>%
  lag_transformer_grouped()

test_data <- data_lags %>%
  group_by(id) %>%
  slice_tail(n = FORECAST_HORIZON) %>%
  ungroup()

train_data <- data_lags %>%
  drop_na()


# rec
recipe_spec <- recipe(value ~ date, train_data) %>%
  step_timeseries_signature(date) %>%
  step_rm(matches("(.iso$)|(.xts$)")) %>%
  step_normalize(matches("(index.num$)|(_year$)")) %>%
  step_dummy(all_nominal()) %>%
  step_fourier(date, K = 1, period = 12)

recipe_spec %>% prep() %>% juice()

# elnet 
model_fit_glmnet <- linear_reg(penalty = 1) %>%
  set_engine("glmnet") 

wflw_fit_glmnet <- workflow() %>%
  add_model(model_fit_glmnet) %>%
  add_recipe(recipe_spec %>% step_rm(date)) %>%
  fit(train_data)  %>%
  recursive(
    transform  = lag_transformer_grouped,
    train_tail = panel_tail(train_data, id, FORECAST_HORIZON),
    id         = "id"
  )

# xgboost    
model_fit_xgboost <- boost_tree("regression", learn_rate = 0.35) %>%
  set_engine("xgboost")

wflw_fit_xgboost <- workflow() %>%
  add_model(model_fit_xgboost) %>%
  add_recipe(recipe_spec %>% step_rm(date)) %>%
  fit(train_data) %>%
  recursive(
    transform  = lag_transformer_grouped,
    train_tail = panel_tail(train_data, id, FORECAST_HORIZON),
    id         = "id"
  )

# mtbl
data_models <- modeltime_table(
  wflw_fit_xgboost,
  wflw_fit_glmnet
)

# ensemble
ensemble_fit <- data_models %>% 
  ensemble_weighted(loadings = c(4, 6)) 

# doing calibration for panel with id, as indicated in the documentation 
calibration_tbl <- modeltime_table(
  ensemble_fit) %>%
  modeltime_calibrate(test_data, id = "id", quiet = FALSE) 

# but no residuals for me :(
calib_out <- calibration_tbl$.calibration_data[[1]] %>% 
  left_join(test_data %>% select(id, date, value))


# Forecast ex ante
data_prepared_tbl <- bind_rows(train_data, test_data)

future_tbl <- data_prepared_tbl %>%
  group_by(id) %>%
  future_frame(.length_out = "2 years") %>%
  ungroup()

# documentation indicates, that calculation of conf_lo and conf_hi should work
ex_ante_obj <- 
  calibration_tbl %>%
  modeltime_forecast(
    new_data    = future_tbl,
    actual_data = data,
    keep_data = TRUE,
    conf_interval = 0.95,
    conf_by_id = TRUE
  )