0

我正在尝试使用 R 包来拟合一些时间序列,tsibble并且fable仍在建设中替代可信赖的 Rob Hyndmanforecast包。该系列全部组合成一个 tsibble,然后我将其与 ARIMA 相匹配,该函数可替换forecast::auto.arima.

我使用map_at,首先迭代除 之外的所有元素Date,然后再次使用 .从适合每个系列的模型中提取模型信息fablelite::components。(很多fable功能真的在fablelite)。

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

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

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 天,缺少周末

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

4

1 回答 1

3

这是一个例子,希望能做你想要的。

首先,您需要创建一个 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 天。

接下来,您可以通过一个简单的函数将 ARIMA 模型拟合到两个系列。

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(),但也许以下之一可以满足您的要求:

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()

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
于 2019-07-25T03:35:59.883 回答