在tidyverts/fable框架中为各时间序列配置不同变换参数

在tidyverts/fable预测框架中，当需要预测多个目标时间序列时，如何为每个序列配置不同的目标变换参数？具体来说，希望对每个时间序列执行Box-Cox变换，但为每个序列使用不同的lambda值（例如通过Guerrero方法为每个序列估算的lambda）。以下是尝试过程中遇到的错误，同时提供了一种变通方案，但不确定是否适用于层级时间序列，希望找到更优实现方式。

数据准备与lambda值估算 #

library(fpp3)

# 构造变换空间下的数据
z1 <- arima.sim(n=104,list(ar=0.9))
z2 <- arima.sim(n=104,list(ma=0.5))

# 逆变换得到原始空间数据
y1 <- fabletools::inv_box_cox(z1, lambda=0.25)
y2 <- fabletools::inv_box_cox(z2, lambda=0.75)

# 创建用于建模的tsibble
tibble(idx=1:104, y1=y1, y2=y2) %>% 
  pivot_longer(cols=c(y1,y2), names_to='series', values_to='value') %>%
  tsibble(index=idx, key=series) ->
  dat

# 用Guerrero方法为每个序列估算最优Box-Cox变换lambda
dat %>% 
  fabletools::features(value, features='guerrero') ->
  lambdas
# # A tibble: 2 × 2
# series lambda_guerrero
# <chr>            <dbl>
# 1 y1              0.0991
# 2 y2              0.751 

错误尝试：直接引入lambda参数 #

尝试1：将lambda作为外生变量

# 将最优lambda作为外生回归变量加入数据
dat %>% inner_join(lambdas, by=join_by(series)) -> dat.xrg

dat.xrg %>%
  model(arima=ARIMA(box_cox(value,lambda=lambda))) ->
  fit
# 错误信息：
# Error in `.g()`:
#   ! Response variable transformation has incompatible lengths, all arguments must be the length of the data 104 or 1.
# Run `rlang::last_trace()` to see where the error occurred.

尝试2：构造等长lambda向量

# 构造与总数据长度匹配的lambda向量
lambdas %>% pull(lambda_guerrero) %>% rep(each=104) -> lambda
length(lambda)
# [1] 208

dat %>%
  model(arima=ARIMA(box_cox(value, lambda=lambda))) ->
  fit
# 错误信息：
# Error in `.g()`:
#   ! Response variable transformation has incompatible lengths, all arguments must be the length of the data 208 or 1.
# Run `rlang::last_trace()` to see where the error occurred.

可行方案：嵌套分序列建模 #

通过嵌套-映射-解嵌套的方式，为每个序列单独配置lambda值并建模，最终仍能得到标准的mable对象，兼容fable的所有工具：

# 分序列嵌套建模，为每个序列指定对应lambda
dat %>% 
  nest(.by=series) %>% 
  inner_join(lambdas, by = "series") %>% 
  mutate(
    fit=map2(
      data, 
      lambda_guerrero, 
      \(.dat,.lambda) 
      model(
        .dat, 
        arima=ARIMA(box_cox(value, lambda=.lambda))
        )
      )
    ) %>% 
  unnest(cols=fit) %>% 
  select(series, arima) %>% 
  as_mable(key='series', model='arima') ->
  fit

# 查看拟合后的mable对象
fit
# # A mable: 2 x 2
# # Key:     series [2]
# series          arima
# <chr>         <model>
#   1 y1     <ARIMA(1,0,2)>
#   2 y2     <ARIMA(0,0,1)>

模型效果验证 #

# 计算模型精度指标
fit %>% accuracy()
# # A tibble: 2 × 11
# series .model .type       ME  RMSE   MAE   MPE  MAPE  MASE RMSSE     ACF1
# <chr>  <chr>  <chr>    <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>    <dbl>
# 1 y1     arima  Training 0.311  3.97 1.91  -450.  512. 0.956 0.961 -0.123  
# 2 y2     arima  Training 0.269  1.16 0.918 -884. 1003. 0.879 0.840  0.00342

# 绘制拟合结果可视化图
fit %>% 
  augment() %>%
  ggplot(aes(x=idx, y=value)) +
  geom_point() +
  geom_line(aes(y=.fitted),color='blue') +
  facet_grid(rows=vars(series))

内容的提问来源于stack exchange，提问作者lowndrul