1

我是 R 的初学者，感谢任何帮助或提示开发函数/循环以自动执行以下预测过程：这是一个虚拟数据集

> class(stack_help) 
[1] "data.frame" 
> stack_help
    OO    GG      CC      DD
1   198.12 60.56   265.5  271.24
2   145.68 52.28   328.9  427.68
3   106.48 47.08  380.24  695.60
4    83.16 43.52  443.94  934.30
5    89.72 46.68   484.6 1084.26
6    86.48 35.46  415.56  924.68
7    93.68 24.40  376.42  798.14
8   101.70 22.68  260.42  427.72
9   115.88 22.00  228.26  245.72
10  137.24 21.60   212.7  140.64
11  129.82 18.78  230.02   46.04
12  145.00 17.62  220.74   47.16
13  135.38 18.84  245.52  143.28
14  146.38 20.68  322.18  490.20
15  154.08 19.48  374.16  621.48
16  149.34 22.68  484.28  999.50
17  152.74 28.90  533.26 1223.58
18  148.62 27.44  456.76  974.44
19  158.54 23.90  417.52  820.54
20  169.96 27.08  306.16  498.02
21  152.50 33.74   283.1  309.22
22  149.68 38.44  224.54  123.82
23  149.48 38.94  215.28   30.26
24  153.38 36.24  193.18   75.46
25  155.58 37.88  243.34  228.92
26  165.84 37.00  318.08  528.58
27  171.34 38.96   393.6  707.04
28  183.60 48.20  531.62 1169.40
29  192.58 44.46  507.96 1037.22
30  207.92 43.52     435  956.96
31  228.88 47.44  399.58  788.78
32  246.14 45.74  262.84  397.66
33  228.92 45.98   240.8  255.32
34  227.52 45.22  211.44   96.02
35  232.92 43.02  203.08   62.18
36  220.16 43.88  188.56   63.74
37  221.76 46.78  210.58  131.28
38  218.94 45.10  272.36  438.64
39  221.00 47.48  351.58  689.90
40  215.82 44.68  402.82  854.80
41  222.32 43.74  435.06 1013.92
42  239.40 52.26  474.24 1128.04
43  249.86 47.62  324.92  689.40
44  240.92 49.60  289.82  538.98
45  221.04 48.40  218.74  256.80
46  191.18 47.34  192.36  136.84
47  206.28 48.66  188.22   60.60
48  226.68 48.12  174.54   58.36
49  226.76 51.66  204.26  190.58
50  223.94 53.40  272.22  454.56
51  219.42 54.50  339.26  647.94
52  219.36 54.68 #VALUE! 1040.08
53  225.94 53.06  462.82 1066.12
54  233.04 52.64  425.32  916.22
55  218.48 64.22  438.06  961.36
56  205.76 56.44  292.24  534.28
57  206.06 53.42  225.32  272.24
58  206.22 52.50   190.2  117.16
59  215.44 52.14  182.12   32.56
60  221.92 51.10  175.82   47.50

感谢任何有关改进以下过程的建议，并且很高兴使用应用功能或循环功能来自动化它。

  1. OO 列是我用来创建预测模型的变量。

  2. 其他列是预测变量，我想测试预测是否与它们一起更好，或者仅与 OO 的过去数据一起使用。

  3. 我进行了 36 次观察，以拟合带有“forecast”包中函数 auto.arima 的 Arima 模型。

  4. 该函数提供了一些模型参数 p,d,q, , 比如说 0,1,0

现在我想以自动化方式测试模型并执行以下操作：

一个。预测未来的下一个时期，在上面的数据表中将相当于第 37 行。

湾。将预测结果与历史数据进行比较，第 37 行，OO 列。

C。从包“forecast”中调用准确度函数并与数据点行 37 进行比较。PLus ，将误差度量存储在向量中。

d。更新 'xdata' 参数，添加历史点 37 和 'xreg' 参数，再增加一个月作为预测变量，并为下一个周期调用另一个预测并重做此过程，直到我完成对 24 个预测的测试。

虽然我为模型安装了包“forecast”，但我发现使用包 astsa 中的函数“sarima.for”更容易。

在代码之前，还有更多信息：

  • Train.OO 将是上述数据表的前 36 个观察值的时间序列对象
  • n.ahead = 预测范围的参数：在这种情况下为 1 个周期
  • 0,1,0 将是 ARIMA 模型 (p,d,q)
  • Train.GG 将是预测变量，即 GG 列的前 36 个观察值
  • newxreg 只是数据表中 TS 对象中的一个数据点的一部分，该数据点将成为预测的预测变量。

现在代码

fc.1 <- sarima.for( 
xdata = Train.00,    
n.ahead = 1, 0, 1, 0 , 
xreg = Train.GG, 
newxreg = window(ts(slack_help$GG, start = c(2009,1), 
frequency = 12), start = c(2012,1) , end = c(2012,11)))
fc.1                      
fc.1.acc <- accuracy(fc.1$pred, 
                 window(ts((slack_help$OO), start = c(2009,1),frequency = 
12),    start = c(2012,1), end = c(2012,1), frequency =12)

现在第二个命令：

fc.2 <- sarima.for( 
xdata = window(ts((slack_help$OO), start = c(2009,1),frequency = 12), 
start = c(2009,1), end = c(2012,1), frequency =12),  
n.ahead = 1, 0, 1, 0 , 
xreg = window(ts((slack_help$GG), start = c(2009,1),frequency = 12), 
start = c(2009,1), end = c(2012,1), frequency =12),
newxreg = window(ts((slack_help$GG), start = c(2009,1),frequency = 12), 
start = c(2009,2), end = c(2012,2), frequency =12),

fc.2
fc.2.acc <- accuracy(fc.2$pred, 
                 window(ts((slack_help$OO), start = c(2009,1),frequency = 
12),  start = c(2012,2), end = c(2012,2), frequency =12)

fc.2.acc

我为以下预测做了这个。基本相同的代码，只是更新了窗口函数的日期，以削减预测中要考虑的正确时间序列。

总共 24 个电话。

我知道这是低效的“蛮力”。但是，我对如何开始开发函数/循环有点迷茫。感谢有关如何自动化上述步骤的任何评论或提示。提前致谢！

4

1 回答 1

0

1/12迭代每月时间序列的最简单方法之一是使用一个月这一事实。举个例子，如果数据开始于 2009 年 1 月，那么我们可以使它等价于2009.000. 这可以通过timeProp <-tsp(stack_help)[1]在使您的数据成为ts-object ( stack_help <- ts(stack_help, start=c(2009,1), freq=12)) 之后使用来提取。那么 2009 年 2 月是timeProp + 1/12= 2009.083，2009 年 3 月是timeProp+2/12= 2009.167。2011 年 1 月是timeProp+24/12=2011.000等等。让我们应用这个：

library(forecast)

#first define the ts-object so we don't have to repeat it every time we use it
stack_help <- ts(stack_help, start=c(2009,1), freq=12)

#extract the start date
timeProp <- tsp(stack_help)[1]

#set vector for storing accuracy measures
accuracyMeasure <- 0
#set counter for above vector
k <- 1

#star the loop with the minimum length of data you want to use to estimate model.
#In this case start with a length of 36 (Jan 2009 - Des 2011)
for(i in 36:(nrow(stack_help)-1)){

  #create new train and test set for each iteration. This way makes your code
  #clearer and more transparent and easier to maintain in the future
  train <- window(stack_help, end=timeProp+((i-1)/12))
  test <- window(stack_help, start=timeProp+i/12, end=(timeProp+i/12))

  #Estimate the model. This can be changed to e.g. auto.arima(). Haven't tried
  #with sarima.for, but should be straight forward to use that as well.
  arrdarr <- Arima(train[,1], order=c(0,1,0), xreg=train[,2])

  #Forecast h=1 with the new xreg (h=1 is automatic since nrow(test)=1)
  foreArima <- forecast(arrdarr, xreg=test[,2])

  #Extract the test MAPE accuracy. 2 selects the test accuracy. "MAPE" can be changed to extract
  #others, e.g. #ME", "MASE". Remember that with h=1, you can scratch the "mean" part of the measures.
  accuracyMeasure[k] <- accuracy(foreArima, test[,1])[2,"MAPE"]

  k <- k+1
}
于 2018-03-08T08:36:08.483 回答