交叉验证季节性线性模型

Question

我正在尝试对我的线性模型执行 CV，该模型具有季节性虚拟变量，因此我无法随机抽样。

y = rnorm(120,0,3) + 20*sin(2*pi*(1:120)/12) 
x = months(ISOdate(2012,1:12,1))
reg.data = data.frame(y, x)
model = lm(y ~ x, data = reg.data)

我的简历功能是：

cross.valid = function(model, min.fit = as.integer(nrow(model$model)*0.7), h = 1)
{
  dados = model$model
  n.rows = nrow(dados)

  results = data.frame(pred = numeric(), actual = numeric())

  for (i in seq(1, n.rows - min.fit - h + 1, by = h))
  {
   dados.train = dados[1:(i + min.fit - 1), ]
   model <- update(model, data = dados.train)

   dados.pred = dados[(i + min.fit):(i + min.fit + h - 1), -1, drop = FALSE]

   predic = predict(model, newdata = dados.pred, interval = 'prediction')
   actual = dados[(i + min.fit):(i + min.fit + h - 1), 1]
   results = rbind(results, data.frame(pred = predic[1:h, 'fit'], actual = actual))
  }

  results
}

例子：

cv1 = cross.valid(model, h = 1)
mae = with(cv1, mean(abs(actual - pred )))
print(mae)

不同视野 (h) 的 MAE 值太接近。代码本身是否有效？有没有更好的解决方案/包来做这件事？

谢谢！

Answer 1

我认为您的功能没有任何不妥之处。调查forecast 包；我怀疑它会提供很多你需要的功能。

我已经简洁地重写了你的函数：

set.seed(1)
y = rnorm(120,0,3) + 20*sin(2*pi*(1:120)/12) 
x = months(ISOdate(2012,1:12,1))
reg.data = data.frame(y, x)

pred.set<-function(i,h) {
  train<-reg.data[1:(i + min.fit - 1),]
  test<-reg.data[(i + min.fit):(i + min.fit + h - 1),]
  pred<-predict(lm(y~x, data=train), newdata=test)
  abs(test$y - pred)
}

pred.by.horiz<-function(h) 
               mean(sapply(seq(1, nrows - min.fit - h + 1, by = h),pred.set,h=h))

pred.by.horiz 与您的函数（和后处理）的输出完全匹配。

正如你提到的，地平线似乎不会影响 MAE：

mae.by.h<-sapply(seq(nrows-min.fit),pred.by.horiz)
plot(mae.by.h,type='l',col='red',lwd=2,xlab='Horizon',ylab='Mean absolute error')

也许您预计平均误差会随着预测范围的增加而增加。对于许多时间序列模型而言，这是正确的，但在您的线性模型中，添加更多数据并不能帮助您预测序列中的下一个点（除非您添加 12 个月或更长时间）。

例如，考虑当h 为 1 时会发生什么。您从 84 个月的数据开始，每个月有 7 个数据点。现在，您添加一个数据点，即明年 1 月，并尝试预测 2 月的结果。但是你的额外数据点只会帮助你预测下一个一月，这就是你的线性函数的工作方式。看模型总结：

lm(y ~ x, data = reg.data)
Coefficients:
(Intercept)      xAugust    xDecember    xFebruary     xJanuary  
   17.11380    -32.74962    -17.81076     -0.03235     -6.63998  
      xJuly        xJune       xMarch         xMay    xNovember  
  -26.69203    -17.41170      2.96735     -7.11166    -25.43532  
   xOctober   xSeptember  
  -33.56517    -36.93474 

每个预测仅基于两个变量：截距和预测月份。所以预测提前一分并不比预测提前五分容易。这就是为什么 MAE 没有随着视野的增加而上升，问题在于您对数据建模的方式，而不是 MAE 函数。

关于您的函数，我没有完全理解的一件事是，为什么您决定在每次迭代时将训练集的大小增加 h。看看当你尝试增加 1 时会发生什么，这很有启发性：

# Code to increment by 1
pred.by.horiz2<-
  function(h) mean(sapply(seq(1, nrows - min.fit - h + 1, by = 1),pred.set,h=h))
mae.by.h2<-sapply(seq(nrows-min.fit),pred.by.horiz2)
plot(mae.by.h2,type='l',col='red',lwd=2,xlab='Horizon',ylab='Mean absolute error')

这里的模式很复杂，但您会注意到 MAE 从 12 开始下降，此时地平线足够大，可以使用下一个点。