如何将预测列表转换为R中的规范化表

Posted 2023-03-25

【中文标题】如何将预测列表转换为R中的规范化表

【英文标题】：How to convert a list of forecasts into a normalized table in R

【发布时间】：2021-11-21 16:56:06

【问题描述】：

我在 R 中使用 auto.arima() 和 lapply 进行一些时间序列预测，以生成大量这样的商店的预测：

my_data_set 是一个小标题列表，每个小标题都包含一个商店名称、日期（每月）和一个销售额

mod <- function(x) forecast(auto.arima(x$Value), h)
all_arima <- lapply(my_data_set, mod)

mod 是预测列表（每个商店一个）

我正在努力将其转化为更易消耗的输出，该输出将在每个存储/预测期间有一行，为上/下置信区间 (80、95) 提供一列，为平均预测提供一列。

如果有更好的方法从一开始就设置它，我会喜欢关于如何以不同方式处理它的建议。

【问题讨论】：

如果您包含一个简单的reproducible example，其中包含可用于测试和验证可能解决方案的示例输入和所需输出，则更容易为您提供帮助。

【参考方案1】：

您可以在这里使用两种方法：(1) 使用建议的预测包； (2) 使用专为这个问题设计的 fable 包。

首先，让我们创建一些示例合成数据。

library(tibble)
library(dplyr)
df <- tibble(
  Store = rep(c("A", "B"), c(200,200)),
  Month = rep(seq(as.Date("1995-01-01"), length=200, by="1 month"), 2),
  Value = rnorm(400)
)

对于预测包，我们会将数据拆分为一个小标题列表。我们可以使用as.data.frame() 函数来轻松将预测对象转换为数据框。

# Using forecast package
library(forecast)
my_data_set <- split(df, df$Store)
mod <- function(x) 
  x$Value %>%
    ts(frequency=12, start=lubridate::year(x$Month[1])) %>%
    auto.arima() %>%
    forecast() %>%
    as.data.frame()

lapply(my_data_set, mod)
#> $A
#>          Point Forecast     Lo 80    Hi 80  Lo 95 Hi 95
#> Sep 2011              0 -1.327999 1.327999 -2.031 2.031
#> Oct 2011              0 -1.327999 1.327999 -2.031 2.031
#> Nov 2011              0 -1.327999 1.327999 -2.031 2.031
#> Dec 2011              0 -1.327999 1.327999 -2.031 2.031
#> Jan 2012              0 -1.327999 1.327999 -2.031 2.031
#> Feb 2012              0 -1.327999 1.327999 -2.031 2.031
#> Mar 2012              0 -1.327999 1.327999 -2.031 2.031
#> Apr 2012              0 -1.327999 1.327999 -2.031 2.031
#> May 2012              0 -1.327999 1.327999 -2.031 2.031
#> Jun 2012              0 -1.327999 1.327999 -2.031 2.031
#> Jul 2012              0 -1.327999 1.327999 -2.031 2.031
#> Aug 2012              0 -1.327999 1.327999 -2.031 2.031
#> Sep 2012              0 -1.327999 1.327999 -2.031 2.031
#> Oct 2012              0 -1.327999 1.327999 -2.031 2.031
#> Nov 2012              0 -1.327999 1.327999 -2.031 2.031
#> Dec 2012              0 -1.327999 1.327999 -2.031 2.031
#> Jan 2013              0 -1.327999 1.327999 -2.031 2.031
#> Feb 2013              0 -1.327999 1.327999 -2.031 2.031
#> Mar 2013              0 -1.327999 1.327999 -2.031 2.031
#> Apr 2013              0 -1.327999 1.327999 -2.031 2.031
#> May 2013              0 -1.327999 1.327999 -2.031 2.031
#> Jun 2013              0 -1.327999 1.327999 -2.031 2.031
#> Jul 2013              0 -1.327999 1.327999 -2.031 2.031
#> Aug 2013              0 -1.327999 1.327999 -2.031 2.031
#> 
#> $B
#>          Point Forecast     Lo 80    Hi 80     Lo 95    Hi 95
#> Sep 2011              0 -1.274651 1.274651 -1.949411 1.949411
#> Oct 2011              0 -1.274651 1.274651 -1.949411 1.949411
#> Nov 2011              0 -1.274651 1.274651 -1.949411 1.949411
#> Dec 2011              0 -1.274651 1.274651 -1.949411 1.949411
#> Jan 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Feb 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Mar 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Apr 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> May 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Jun 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Jul 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Aug 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Sep 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Oct 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Nov 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Dec 2012              0 -1.274651 1.274651 -1.949411 1.949411
#> Jan 2013              0 -1.274651 1.274651 -1.949411 1.949411
#> Feb 2013              0 -1.274651 1.274651 -1.949411 1.949411
#> Mar 2013              0 -1.274651 1.274651 -1.949411 1.949411
#> Apr 2013              0 -1.274651 1.274651 -1.949411 1.949411
#> May 2013              0 -1.274651 1.274651 -1.949411 1.949411
#> Jun 2013              0 -1.274651 1.274651 -1.949411 1.949411
#> Jul 2013              0 -1.274651 1.274651 -1.949411 1.949411
#> Aug 2013              0 -1.274651 1.274651 -1.949411 1.949411

要使用 fable 包，我们可以只取包含所有商店的原始数据框，并将其转换为 tsibble 对象，然后将其通过管道传递给模型和预测，如下所示。

# Using fable
library(tsibble)
library(fable)
df %>%
  mutate(Month = yearmonth(Month)) %>%
  as_tsibble(index=Month, key=Store) %>%
  model(ARIMA(Value)) %>%
  forecast() %>%
  mutate(
    pi80 = hilo(Value, 80),
    pi95 = hilo(Value, 95)
  ) %>%
  unpack_hilo(cols = c(pi80, pi95))
#> # A fable: 48 x 9 [1M]
#> # Key:     Store, .model [2]
#>    Store .model          Month     Value .mean pi80_lower pi80_upper pi95_lower
#>    <chr> <chr>           <mth>    <dist> <dbl>      <dbl>      <dbl>      <dbl>
#>  1 A     ARIMA(Value) 2011 Sep N(0, 1.1)     0      -1.33       1.33      -2.03
#>  2 A     ARIMA(Value) 2011 Oct N(0, 1.1)     0      -1.33       1.33      -2.03
#>  3 A     ARIMA(Value) 2011 Nov N(0, 1.1)     0      -1.33       1.33      -2.03
#>  4 A     ARIMA(Value) 2011 Dec N(0, 1.1)     0      -1.33       1.33      -2.03
#>  5 A     ARIMA(Value) 2012 Jan N(0, 1.1)     0      -1.33       1.33      -2.03
#>  6 A     ARIMA(Value) 2012 Feb N(0, 1.1)     0      -1.33       1.33      -2.03
#>  7 A     ARIMA(Value) 2012 Mar N(0, 1.1)     0      -1.33       1.33      -2.03
#>  8 A     ARIMA(Value) 2012 Apr N(0, 1.1)     0      -1.33       1.33      -2.03
#>  9 A     ARIMA(Value) 2012 May N(0, 1.1)     0      -1.33       1.33      -2.03
#> 10 A     ARIMA(Value) 2012 Jun N(0, 1.1)     0      -1.33       1.33      -2.03
#> # … with 38 more rows, and 1 more variable: pi95_upper <dbl>

^{由reprex package 创建于 2021-09-30 (v2.0.1)}

这种方法也更加灵活，因为可能有多个分组变量（例如，商店和产品）。寓言方法记录在开放获取教科书https://OTexts.com/fpp3。

【讨论】：

感谢 Rob - 这太棒了！我不熟悉这些其他软件包。快速提问 - 这个实现模型（ARIMA(Value)）是否允许我指定预测期的数量？

预测周期数在预测函数中指定，而不是模型函数。

【参考方案2】：

forecast 的输出是list。我们可以提取组件并转换为data.frame

library(forecast)
mod <- function(x) 
       frcst <- forecast(auto.arima(x$Value), h)
       data.frame(Mean = as.numeric(frcst$mean), 
                  lower = as.numeric(frcst$lower[, "95%"]),
                  upper = as.numeric(frcst$upper[, "95%]))
   

然后应用函数

lapply(my_data_set, mod)