r - 使用 Fable 聚合预测

Question

问题： 使用 fable 我可以轻松地生成具有分组结构的时间序列的预测，甚至可以使用 Fable 的aggregate_key/reconcile语法来生成连贯的顶级预测。但是，我无法使用此方法轻松访问汇总预测，而我使用的替代方法涉及放弃寓言（预测表）结构。谁能告诉我是否有更简单/有意的方法来使用这个包？正如您在示例中看到的那样，我可以使用其他方法到达那里，但我想知道是否有更好的方法。任何帮助都感激不尽！

方法 1：aggregate_key我在不使用/ 的情况下总结预测的努力reconcile主要是使用 dplyrgroup_by和summarise，但是预测的预测区间被格式化为正态分布对象，似乎不支持使用这种方法求和。为了解决这个问题，我一直在使用hilo并unpack_hilo提取不同预测区间的界限，然后可以使用通常的方法对其进行求和。但是我真的很想保留寓言结构和分布对象，使用这种方法是不可能的。

方法 2： 替代方案，使用aggregate_key/ reconcileonly 似乎支持使用聚合min_trace。我知道这种方法是为了实现最佳对账，而我想要的是一个简单的自下而上的汇总预测。感觉应该有一种简单的方法可以使用这种语法获得自下而上的预测，但到目前为止我还没有找到。此外，即使使用min_trace我也不确定如何访问聚合预测本身，如您在示例中所见！

使用方法 1 的示例：

library(fable)
#> Loading required package: fabletools
library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union

lung_deaths_agg <- as_tsibble(cbind(mdeaths, fdeaths))
  
fc_1 <- lung_deaths_agg %>% 
  model(lm = TSLM(value ~ trend() + season())) %>% 
  forecast()

fc_1
#> # A fable: 48 x 5 [1M]
#> # Key:     key, .model [2]
#>    key     .model    index        value .mean
#>    <chr>   <chr>     <mth>       <dist> <dbl>
#>  1 fdeaths lm     1980 Jan N(794, 5940)  794.
#>  2 fdeaths lm     1980 Feb N(778, 5940)  778.
#>  3 fdeaths lm     1980 Mar N(737, 5940)  737.
#>  4 fdeaths lm     1980 Apr N(577, 5940)  577.
#>  5 fdeaths lm     1980 May N(456, 5940)  456.
#>  6 fdeaths lm     1980 Jun N(386, 5940)  386.
#>  7 fdeaths lm     1980 Jul N(379, 5940)  379.
#>  8 fdeaths lm     1980 Aug N(335, 5940)  335.
#>  9 fdeaths lm     1980 Sep N(340, 5940)  340.
#> 10 fdeaths lm     1980 Oct N(413, 5940)  413.
#> # ... with 38 more rows

fc_1 %>%
  hilo() %>% 
  unpack_hilo(c(`80%`, `95%`)) %>% 
  as_tibble() %>% 
  group_by(index) %>% 
  summarise(across(c(.mean, ends_with("upper"), ends_with("lower")), sum))
#> `summarise()` ungrouping output (override with `.groups` argument)
#> # A tibble: 24 x 6
#>       index .mean `80%_upper` `95%_upper` `80%_lower` `95%_lower`
#>       <mth> <dbl>       <dbl>       <dbl>       <dbl>       <dbl>
#>  1 1980 Jan 2751.       3089.       3267.       2414.       2236.
#>  2 1980 Feb 2687.       3024.       3202.       2350.       2171.
#>  3 1980 Mar 2535.       2872.       3051.       2198.       2020.
#>  4 1980 Apr 2062.       2399.       2577.       1725.       1546.
#>  5 1980 May 1597.       1934.       2113.       1260.       1082.
#>  6 1980 Jun 1401.       1738.       1916.       1064.        885.
#>  7 1980 Jul 1343.       1680.       1858.       1006.        827.
#>  8 1980 Aug 1200.       1538.       1716.        863.        685.
#>  9 1980 Sep 1189.       1527.       1705.        852.        674.
#> 10 1980 Oct 1482.       1819.       1998.       1145.        967.
#> # ... with 14 more rows

使用方法 2 的示例：

fc_2 <- lung_deaths_agg %>%
  aggregate_key(key, value = sum(value)) %>% 
  model(lm = TSLM(value ~ trend() + season())) %>%
  reconcile(lm = min_trace(lm)) %>% 
  forecast()

fc_2
#> # A fable: 72 x 5 [1M]
#> # Key:     key, .model [3]
#>    key     .model    index        value .mean
#>    <chr>   <chr>     <mth>       <dist> <dbl>
#>  1 fdeaths lm     1980 Jan N(794, 5606)  794.
#>  2 fdeaths lm     1980 Feb N(778, 5606)  778.
#>  3 fdeaths lm     1980 Mar N(737, 5606)  737.
#>  4 fdeaths lm     1980 Apr N(577, 5606)  577.
#>  5 fdeaths lm     1980 May N(456, 5606)  456.
#>  6 fdeaths lm     1980 Jun N(386, 5606)  386.
#>  7 fdeaths lm     1980 Jul N(379, 5606)  379.
#>  8 fdeaths lm     1980 Aug N(335, 5606)  335.
#>  9 fdeaths lm     1980 Sep N(340, 5606)  340.
#> 10 fdeaths lm     1980 Oct N(413, 5606)  413.
#> # ... with 62 more rows

fc_2 %>% as_tibble() %>% select(key) %>% slice(50:55)
#> # A tibble: 6 x 1
#>   key         
#>   <chr>       
#> 1 <aggregated>
#> 2 <aggregated>
#> 3 <aggregated>
#> 4 <aggregated>
#> 5 <aggregated>
#> 6 <aggregated>

fc_2 %>% as_tibble() %>% select(key) %>% filter(key == "<aggregated>")
#> # A tibble: 0 x 1
#> # ... with 1 variable: key <chr>

Answer 1

方法一：

在将事物加在一起时，使用分布需要更多的注意（而不是数字）。更具体地说，可以毫无问题地添加正态分布的平均值：

library(distributional)
mean(dist_normal(2,3) + dist_normal(4,1))
#> [1] 6
mean(dist_normal(2,3)) + mean(dist_normal(4,1))
#> [1] 6

^{由reprex 包（v0.3.0）于 2020-07-03 创建}

然而，分位数（用于产生 80% 和 95% 的区间）不能：

library(distributional)
quantile(dist_normal(2,3) + dist_normal(4,1), 0.9)
#> [1] 10.05262
quantile(dist_normal(2,3), 0.9) + quantile(dist_normal(4,1), 0.9)
#> [1] 11.12621

^{由reprex 包（v0.3.0）于 2020-07-03 创建}

如果要聚合分布，则需要计算分布本身的总和：

library(fable)
library(dplyr)
lung_deaths_agg <- as_tsibble(cbind(mdeaths, fdeaths))

fc_1 <- lung_deaths_agg %>% 
  model(lm = fable::TSLM(value ~ trend() + season())) %>% 
  forecast()
fc_1 %>% 
  summarise(value = sum(value), .mean = mean(value))
#> # A fable: 24 x 3 [1M]
#>       index          value .mean
#>       <mth>         <dist> <dbl>
#>  1 1980 Jan N(2751, 40520) 2751.
#>  2 1980 Feb N(2687, 40520) 2687.
#>  3 1980 Mar N(2535, 40520) 2535.
#>  4 1980 Apr N(2062, 40520) 2062.
#>  5 1980 May N(1597, 40520) 1597.
#>  6 1980 Jun N(1401, 40520) 1401.
#>  7 1980 Jul N(1343, 40520) 1343.
#>  8 1980 Aug N(1200, 40520) 1200.
#>  9 1980 Sep N(1189, 40520) 1189.
#> 10 1980 Oct N(1482, 40520) 1482.
#> # … with 14 more rows

^{由reprex 包（v0.3.0）于 2020-07-03 创建}

请注意，这将需要 fabletools (>=0.2.0.9000) 和分布式 (>=0.1.0.9000) 的开发版本，因为我添加了新功能以使此示例正常工作。

方法二：

自下而上协调的实验性支持可使用fabletools:::bottom_up(). 这目前是一个内部功能，因为我仍在研究如何在 fabletools 中更普遍地进行协调的一些细节。

匹配聚合值应使用is_aggregated().

fc_2 <- lung_deaths_agg %>%
  aggregate_key(key, value = sum(value)) %>% 
  model(lm = TSLM(value ~ trend() + season())) %>%
  reconcile(lm = min_trace(lm)) %>% 
  forecast()

fc_2 %>% 
  filter(is_aggregated(key))
#> # A fable: 24 x 5 [1M]
#> # Key:     key, .model [1]
#>    key          .model    index          value .mean
#>    <chr>        <chr>     <mth>         <dist> <dbl>
#>  1 <aggregated> lm     1980 Jan N(2751, 24989) 2751.
#>  2 <aggregated> lm     1980 Feb N(2687, 24989) 2687.
#>  3 <aggregated> lm     1980 Mar N(2535, 24989) 2535.
#>  4 <aggregated> lm     1980 Apr N(2062, 24989) 2062.
#>  5 <aggregated> lm     1980 May N(1597, 24989) 1597.
#>  6 <aggregated> lm     1980 Jun N(1401, 24989) 1401.
#>  7 <aggregated> lm     1980 Jul N(1343, 24989) 1343.
#>  8 <aggregated> lm     1980 Aug N(1200, 24989) 1200.
#>  9 <aggregated> lm     1980 Sep N(1189, 24989) 1189.
#> 10 <aggregated> lm     1980 Oct N(1482, 24989) 1482.
#> # … with 14 more rows

^{由reprex 包（v0.3.0）于 2020-07-03 创建}

将聚合向量与 进行比较"<aggregated>"是不明确的，因为您的键的字符值可能"<aggregated>"没有值是<aggregated>. 我现在更新了 fabletools 以匹配"<aggregated>"带有警告和提示的聚合值，所以这段代码现在给出：

fc_2 %>% 
  filter(key == "<aggregated>")
#> Warning: <aggregated> character values have been converted to aggregated values.
#> Hint: If you're trying to compare aggregated values, use `is_aggregated()`.
#> # A fable: 24 x 5 [1M]
#> # Key:     key, .model [1]
#>    key          .model    index          value .mean
#>    <chr>        <chr>     <mth>         <dist> <dbl>
#>  1 <aggregated> lm     1980 Jan N(2751, 24989) 2751.
#>  2 <aggregated> lm     1980 Feb N(2687, 24989) 2687.
#>  3 <aggregated> lm     1980 Mar N(2535, 24989) 2535.
#>  4 <aggregated> lm     1980 Apr N(2062, 24989) 2062.
#>  5 <aggregated> lm     1980 May N(1597, 24989) 1597.
#>  6 <aggregated> lm     1980 Jun N(1401, 24989) 1401.
#>  7 <aggregated> lm     1980 Jul N(1343, 24989) 1343.
#>  8 <aggregated> lm     1980 Aug N(1200, 24989) 1200.
#>  9 <aggregated> lm     1980 Sep N(1189, 24989) 1189.
#> 10 <aggregated> lm     1980 Oct N(1482, 24989) 1482.
#> # … with 14 more rows

^{由reprex 包（v0.3.0）于 2020-07-03 创建}