R:使用寓言、tsibble 和地图预测多个时间序列 [英] R: Forecasting multiple time series with fable, tsibble and map

查看:113
本文介绍了R:使用寓言、tsibble 和地图预测多个时间序列的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我正在尝试使用 R 包 tsibblefable 来拟合一些时间序列,这是令人敬畏的 Rob Hyndman 的 forecast 的仍在建设中的替代品 包.该系列全部合并为一个 tsibble,然后我将其与 ARIMA 配合使用,该函数除其他外还替换了 forecast::auto.arima.

I am trying to fit some time series using the R packages tsibble and fable, the still-under-construction replacement for the redoubtable Rob Hyndman's forecast package. The series are all combined into one tsibble, which I then fit with ARIMA, a function which replaces, among other things, forecast::auto.arima.

我使用map_at,首先迭代除Date之外的所有元素,然后再次从已经拟合到每个系列的模型中提取模型信息使用 fablelite::components.(很多fable 函数确实在fablelite 中).

I use map_at, first to iterate over all the elements except the Date, and then again to extract the model information from the models that have been fit to each series using fablelite::components. (A lot of the fable functions are really in fablelite).

这失败了,显然是因为组件需要 mdl_df 类的对象,而我的模型对象有 mdl_defn

This fails, apparently because components expects an object of class mdl_df and my model objects have class mdl_defn

这是一个(几乎)重现错误的玩具示例:

Here is a toy example that (almost) reproduces the error:

library(tidyverse)
library(tsibble)
library(fable)
set.seed(1)
ar1  <-  arima.sim(model=list(ar=.6), n=10)
ma1 <- arima.sim(model=list(ma=0.4), n=10)
Date  <- c(ymd("2019-01-01"):ymd("2019-01-10"),  ymd("2019-01-01"):ymd("2019-01-10"))
tb <- tibble(Date, ar1, ma1)

# Fit the whole series
tb_all <- tb   %>% 
map_at(.at =  c("ar1", "ma1"), .f = ARIMA)
names(arima_all[2:3])<- c("ar1", "ma1")

# Extract model components
tb_components <- tb %>%  
  map_at(.at = c("ar1", "ma1"), 
         .f = fablelite::components)

请注意,在这个玩具中,就像我的真实数据一样,时间是以 5 天为单位的,缺少周末

Note that in this toy, like my real data, time is in 5-day weeks with missing weekends

在这个玩具示例中,错误消息说组件函数拒绝列表元素,理由是没有用于 ts 类的方法.在我的真实案例中,它使用更长的系列和更多的系列,但在我看来,其他方面相同,元素被拒绝,因为它们属于 mdl_defn 类.请注意,如果我使用 str( ) 检查 tb_all 的第二个和第三个元素,它们也会显示为 Classes 'mdl_defn', 'R6' 不确定错误信息中的 ts 来自哪里.

In this toy example, the the error message says the components function rejects list elements on grounds there is no method for class ts. In my real case, which uses longer series and more of them, but is to my eye otherwise identical, elements are rejected because they are of class mdl_defn. Note that if I examine the 2nd and third elements of tb_all with str( ), they also display as of Classes 'mdl_defn', 'R6' Not sure where the ts in the error message comes from.

推荐答案

这里有一个例子,希望能达到你想要的效果.

Here is an example that hopefully does something like what you want.

首先,您需要创建一个 tsibble:

First, you need to create a tsibble:

library(tidyverse)
library(tsibble)
library(fable)
library(lubridate)
set.seed(1)
ar1  <-  arima.sim(model=list(ar=.6), n=30)
ma1 <- arima.sim(model=list(ma=0.4), n=30)
Date  <- ymd(paste0("2019-01-",1:30))
tb <- bind_cols(Date=Date, ar1=ar1, ma1=ma1) %>%
  gather("Series", "value", -Date) %>%
  as_tsibble(index=Date, key=Series)
tb
#> # A tsibble: 60 x 3 [1D]
#> # Key:       Series [2]
#>    Date       Series   value
#>    <date>     <chr>    <dbl>
#>  1 2019-01-01 ar1    -2.07  
#>  2 2019-01-02 ar1    -0.118 
#>  3 2019-01-03 ar1    -0.116 
#>  4 2019-01-04 ar1    -0.0856
#>  5 2019-01-05 ar1     0.892 
#>  6 2019-01-06 ar1     1.36  
#>  7 2019-01-07 ar1     1.41  
#>  8 2019-01-08 ar1     1.76  
#>  9 2019-01-09 ar1     1.84  
#> 10 2019-01-10 ar1     1.18  
#> # … with 50 more rows

这包含两个系列:ar1ma1 在相同的 30 天内.

This contains two series: ar1 and ma1 over the same 30 days.

接下来,您可以在一个简单的函数中将 ARIMA 模型安装到两个系列中.

Next you can fit ARIMA models to both series in one simple function.

tb_all <- tb %>% model(arima = ARIMA(value))
tb_all
#> # A mable: 2 x 2
#> # Key:     Series [2]
#>   Series arima                 
#>   <chr>  <model>               
#> 1 ar1    <ARIMA(0,0,2)>        
#> 2 ma1    <ARIMA(0,0,0) w/ mean>

最后,尚不清楚您要使用 components() 提取什么,但也许以下其中一项可以满足您的需求:

Finally, it is not clear what you are trying to extract using components(), but perhaps one of the following does what you want:

tidy(tb_all)
#> # A tibble: 3 x 7
#>   Series .model term     estimate std.error statistic  p.value
#>   <chr>  <chr>  <chr>       <dbl>     <dbl>     <dbl>    <dbl>
#> 1 ar1    arima  ma1         0.810     0.198      4.09 0.000332
#> 2 ar1    arima  ma2         0.340     0.181      1.88 0.0705  
#> 3 ma1    arima  constant    0.295     0.183      1.61 0.118
glance(tb_all)
#> # A tibble: 2 x 9
#>   Series .model sigma2 log_lik   AIC  AICc   BIC ar_roots  ma_roots 
#>   <chr>  <chr>   <dbl>   <dbl> <dbl> <dbl> <dbl> <list>    <list>   
#> 1 ar1    arima   0.695   -36.4  78.9  79.8  83.1 <cpl [0]> <cpl [2]>
#> 2 ma1    arima   1.04    -42.7  89.4  89.8  92.2 <cpl [0]> <cpl [0]>
augment(tb_all)
#> # A tsibble: 60 x 6 [1D]
#> # Key:       Series, .model [2]
#>    Series .model Date         value .fitted  .resid
#>    <chr>  <chr>  <date>       <dbl>   <dbl>   <dbl>
#>  1 ar1    arima  2019-01-01 -2.07    -0.515 -1.56  
#>  2 ar1    arima  2019-01-02 -0.118   -1.21   1.09  
#>  3 ar1    arima  2019-01-03 -0.116    0.511 -0.627 
#>  4 ar1    arima  2019-01-04 -0.0856  -0.155  0.0690
#>  5 ar1    arima  2019-01-05  0.892   -0.154  1.05  
#>  6 ar1    arima  2019-01-06  1.36     0.871  0.486 
#>  7 ar1    arima  2019-01-07  1.41     0.749  0.659 
#>  8 ar1    arima  2019-01-08  1.76     0.699  1.06  
#>  9 ar1    arima  2019-01-09  1.84     1.09   0.754 
#> 10 ar1    arima  2019-01-10  1.18     0.973  0.206 
#> # … with 50 more rows

要以传统方式查看模型输出,请使用 report():

To see model outputs in the traditional way, use report():

tb_all %>% filter(Series=='ar1') %>% report()
#> Series: value 
#> Model: ARIMA(0,0,2) 
#> 
#> Coefficients:
#>          ma1     ma2
#>       0.8102  0.3402
#> s.e.  0.1982  0.1809
#> 
#> sigma^2 estimated as 0.6952:  log likelihood=-36.43
#> AIC=78.86   AICc=79.78   BIC=83.06
tb_all %>% filter(Series=='ma1') %>% report()
#> Series: value 
#> Model: ARIMA(0,0,0) w/ mean 
#> 
#> Coefficients:
#>       constant
#>         0.2950
#> s.e.    0.1833
#> 
#> sigma^2 estimated as 1.042:  log likelihood=-42.68
#> AIC=89.36   AICc=89.81   BIC=92.17

这篇关于R:使用寓言、tsibble 和地图预测多个时间序列的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆