将不同的预测方法传递给R中的分层时间序列预测？答案

【问题标题】：Passing different forecasting method to hierarchical time series forecast in R?将不同的预测方法传递给R中的分层时间序列预测？
【发布时间】：2015-04-24 13:31:35
【问题描述】：

我有一个分层时间序列，其底层序列都表现出间歇性需求。使用 Hyndman 的 HTS 包在层次结构中进行最佳组合似乎是有利的。使用 Kourentzes 的 MAPA 包对间歇性需求进行多重聚合预测似乎也很有利。本质上，我想做这样的事情：

forecast(my_hts, method='comb', fmethod='MAPA')

但是，我不清楚是否/如何将两者结合起来，因为forecast.gts() 只接受fmethod=c("ets", "arima", "rw").

有没有一种巧妙的方法可以将不同的预测方法传递给forecast.gts()，而无需撕毁代码？

举例说明我的意思：

library(hts)
library(MAPA)
set.seed(1)

#note intermittent demand of bottom level time series
x <- ts(rpois(365, lambda=0.05), frequency=365, start=2014)
y <- ts(rpois(365, lambda=0.07), frequency=365, start=2014)

#it's easy to make a MAPA forecast for the top-level time series
#but this isn't an optimal hierarchical forecast
mapasimple(x+y)

#it's also easy to make this a HTS and make an optimal hierarchical forecast
#but now I cannot use MAPA
z <- hts(data.frame(x,y)))
z_arima <- forecast(z, fmethod="arima")
z_rw <- forecast(z, fmethod="rw")
z_ets <- forecast(z, fmethod="ets")

#z_MAPA <- ?

【问题讨论】：

感谢您添加可重现的示例。我们现在可以尝试将其迁移到Stack Overflow，在那里编程问题属于并且很容易回答。
我显然遵从您的要求，但由于上述软件包的两位作者（据我所知）都活跃在 CrossValidated 中，我认为这是一个比一般 SO 更好的发布位置。
这是一个合理的猜测，@user1569317，这是一个棘手且有争议的问题。事实上，SO 上有一大群活跃的 R 用户（我相信比这里还多）。我的标准是 OP 需要解释什么；如果是统计概念，Q就属于这里，如果是代码如何工作，Q就属于SO。我将您的 Q 读为“是否有一种巧妙的方法可以将不同的预测方法传递给 predict.gts() 而无需破坏代码？”，而不是“预测如何工作”。因此我认为SO更合适。（我也相信你会在那里获得更好/更快的 A。）

标签： r time-series forecasting hierarchical

【解决方案1】：

我在过去使用 MAPA 预测和 combinef 的小型实验中成功运行。斯蒂芬关于负面预测的观点仍然是一个问题。您可能有兴趣查看 tsintermittent 软件包中的 iMAPA，该软件包专门用于生成间歇性需求时间序列的预测。与实现 Croston 方法的 MAPA 相比，SBA 和指数平滑可以在各种时间聚合级别中进行选择。

我很想进一步了解您的结果。也许 E. Spiliotis 等人在去年的 ISF 上提出了一项相关研究。 Examining the effect of temporal aggregation on forecasting accuracy for hierarchical energy consumption time series。

【讨论】：

感谢您的提醒 - iMAPA 和演示文稿看起来都非常有希望。如果我想出一些有趣的事情，我会通过电子邮件发送结果。

【解决方案2】：

我发帖是因为在仔细查看 hts 文档（在此处插入当之无愧的 RTFM）后，我想我找到了使用 hts 中的 combinef() 函数的解决方法，该函数可用于优化组合预测在forecast.gts() 环境之外。在接受作为答案之前，我会先搁置一段时间，以便其他人可以告诉我我是否错了。

fh <- 8

library(hts)
library(MAPA)
set.seed(1)

x <- ts(rpois(365, lambda=0.05), frequency=365, start=2014)
y <- ts(rpois(365, lambda=0.07), frequency=365, start=2014)

my_hts <- hts(data.frame(x,y))

ally <- aggts(my_hts)

allf <- matrix(NA, nrow = fh, ncol = ncol(ally))

for(i in 1:ncol(ally)){
    allf[,i] <- mapafor(ally[,i], 
                        mapaest(ally[,i], outplot=0), 
                        fh = fh, 
                        outplot=0)$outfor
}

allf <- ts(allf)

y.f <- combinef(allf, my_hts$nodes, weights=NULL, keep="bottom")

#here's what the non-reconciled, bottom-level MAPA forecasts look like
print(allf[1,1:2])

 Series 1   Series 2
1 0.1343304 0.06032574

#here's what the reconciled MAPA bottom-level forecasts look like
#notice that they're different
print(y.f[1,])

[1] 0.06030926 0.07402938

【讨论】：

这对我来说完全正确。（我冒昧地纠正了一些错别字。）它应该在您的设置中工作，并具有时间不变的预测。如果您的预测随时间而变化（例如，趋势、季节性或因果因素），那么您可能会遇到我上面提到的问题：您的“最佳”预测可能会变成负数。否则很好，+1！

【解决方案3】：

“有没有一种巧妙的方法可以将不同的预测方法传递给forecast.gts()，而不必撕毁代码？” - 几乎可以肯定不是。 forecast.gts() 不允许您按照 family 参数到 glm() 和类似的行插入预测方法。

运行 MAPA 预测然后自己重新实现最佳组合可能是最简单的方法。这真的不是那么难，在hts 包出现之前我自己做了几次。查看hts 包中的引用；其中大部分都可以在 Rob Hyndman 的网站上以预印本的形式获得。

将最优组合方法与间歇性需求预测相结合的一个关键问题是最优组合很可能产生负预测。在此处使用的 GLS 意义上，这些可能是“最佳的”，但它们对需求没有意义。我的解决方案是使用 mgcv 包中的 pcls() 进行实际对账，以将解决方案（即底层对账预测）限制为非负数。

【讨论】：

感谢您的回答。介意看看下面我提出的解决方案吗？如果无效，我会接受你的回答。