使用 2 个数据框计算和绘制预测并在数据框中存储方法

Question

我想用不同的外生变量预测 v1 到 v6，exo1 到 exo3

    > ts.cnt
           v1   v2   v3   v4  v5   v6
Jan 2012  892 1091  175  615   9  139
Feb 2012  749 1057  148  701  20  123
Mar 2012 1021 1312  245  824  23  151
Apr 2012  878 1178  243  811  16  131
May 2012  894 1249  223  799   5  132
Jun 2012  925 1141  242  782  15  117

及其外生变量：

> head(ev)
 exo1    exo2    exo3
75155   77628   76113
82914   76113   77653
82719   77653   81304
79342   81304   81745
76493   81745   79163
76911   79163   77188

以下是一个 for 循环代码，通过遍历 v1-7 和 exo1-3 的每个组合来运行模型：

for (i in 1:6) {
  for (j in 1:ncol(ev)) {
  fit <- auto.arima(x = ts.cnt[,i], 
                    xreg = ev[,j])

  fcast <- forecast(fit, h = 12, xreg = ev[1:12,j])
  plot.forecast(fcast, 
                main = substitute(paste(b,' with ', a), 
                                  list(a = colnames(ev)[j], 
                                       b = colnames(ts.cnt[i])
                                  )
                )
  lines(fitted(fcast), col = 2)
  }
}

我还想将预测值 fcast$mean 保存到循环中的对象中，以便最终结果如下所示：

          exo   v1   v2   v3  v4   v5    v6
Jan 2015 exo1 1091  175  615   9  139  2924
Feb 2015 exo1 1057  148  701  20  123  2798
Mar 2015 exo1 1312  245  824  23  151  3576
....
Jan 2015 exo2 1178  243  811  16  131  3257
Feb 2015 exo2 1249  223  799   5  132  3304
Mar 2015 exo2 1141  242  782  15  117  3224
....
Jan 2015 exo3 1234  243  811  16  131  3257
Feb 2015 exo3 1249  223  799   5  132  3304
Mar 2015 exo3 1111  242  782  15  117  3224

提前致谢。

Answer 1

我或多或少使用了相同的循环逻辑，尽管：

我没有遍历列索引，而是遍历了列名称，这简化了很多事情
我在 level-1 循环中使用了 ev，在 level-2 中使用了 ts.cnt方便写入 output.df

我还在预测电话中将h=12更改为h=6，并将1:12替换为1:6，假设12和1:12 是对所提供样本数据的维度的疏忽。

由于我不熟悉这个字段，请确保生成的输出有意义。

output.df <- data.frame(per=rep(rownames(ts.cnt), 3), exo=rep(1:3, each=6), 
                        v1=NA, v2=NA, v3=NA, v4=NA, v5=NA, v6=NA)    

rows.offset <- 0 # defines row offset used when writing to df.output

for(ev.col in colnames(ev)) {
    for (ts.cnt.col in colnames(ts.cnt)) {
        fit <- auto.arima(x = ts.cnt[,ts.cnt.col], xreg = ev[,ev.col])
        fcast <- forecast(fit, h=6, xreg=ev[1:6, ev.col])
        plot.forecast(fcast, main = paste(ts.cnt.col, "with", ev.col))
        lines(fitted(fcast), col = 2)
        for(i in seq_along(fcast$mean)) {
            output.df[rows.offset + i, ts.cnt.col] <- fcast$mean[i]
        }
    }
    rows.offset <- rows.offset + 6
}

结果

        per exo       v1       v2       v3       v4       v5       v6
1  Jan 2012   1 848.9718 1114.304 202.0891 766.3091 14.14224 125.7534
2  Feb 2012   1 936.6196 1229.345 222.9528 835.6403 15.60229 138.7362
3  Mar 2012   1 934.4169 1226.453 222.4284 833.8978 15.56559 138.4099
4  Apr 2012   1 896.2693 1176.384 213.3478 803.7224 14.93013 132.7593
5  May 2012   1 864.0862 1134.142 205.6869 778.2649 14.39402 127.9922
6  Jun 2012   1 868.8081 1140.340 206.8109 782.0000 14.47268 128.6916
7  Jan 2012   2 878.2537 1152.167 209.6261 743.3959 14.32666 129.8641
8  Feb 2012   2 861.1136 1129.681 205.5350 728.8877 14.04706 127.3296
9  Mar 2012   2 878.5365 1152.538 209.6936 743.6353 14.33127 129.9059
10 Apr 2012   2 919.8426 1206.726 219.5528 778.5987 15.00509 136.0137
11 May 2012   2 924.8319 1213.272 220.7436 782.8219 15.08647 136.7514
12 Jun 2012   2 895.6201 1174.949 213.7712 758.0957 14.60995 132.4320
13 Jan 2012   3 862.2396 1131.347 205.7280 755.9879 14.20652 127.5946
14 Feb 2012   3 879.6854 1154.237 209.8905 793.2517 14.49396 130.1762
15 Mar 2012   3 921.0454 1208.506 219.7589 881.5961 15.17542 136.2967
16 Apr 2012   3 926.0412 1215.061 220.9509 892.2671 15.25773 137.0359
17 May 2012   3 896.7913 1176.682 213.9719 829.7897 14.77580 132.7075
18 Jun 2012   3 874.4177 1147.325 208.6336 782.0000 14.40717 129.3967

示例图

最后的笔记

如果需要，您可以使用以下代码将拟合和预测保存在列表中，例如：

fits <- list()
fcasts <- list()
for (ts.cnt.col in colnames(ts.cnt)) {
    for(ev.col in colnames(ev)) {
        fits[[paste(ts.cnt.col, ev.col, sep=".")]] <- 
            auto.arima(x = ts.cnt[,ts.cnt.col], xreg = ev[,ev.col])
        fcasts[[paste(ts.cnt.col, ev.col, sep=".")]] <- 
            forecast(fits[[paste(ts.cnt.col, ev.col, sep=".")]], h=6, xreg=ev[1:6, ev.col])
    }
}

这样您就可以使用 lapply 等函数对列表中的所有项目重复相同的操作。

数据

ts.cnt <- read.csv(text = "v1,v2,v3,v4,v5,v6
Jan 2012,892,1091,175,615,9,139
Feb 2012,749,1057,148,701,20,123
Mar 2012,1021,1312,245,824,23,151
Apr 2012,878,1178,243,811,16,131
May 2012,894,1249,223,799,5,132
Jun 2012,925,1141,242,782,15,117")

ev <- read.csv(text="exo1,exo2,exo3
75155,77628,76113
82914,76113,77653
82719,77653,81304
79342,81304,81745
76493,81745,79163
76911,79163,77188")

使用 2 个数据框计算和绘制预测并在数据框中存储方法

Calculating and plotting forecasts using 2 dataframes and storing means in a dataframe

r

loops

for-loop

time-series