r - 交叉验证季节性线性模型

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 处开始下降，此时地平线足够大，可以使用下一个点。