24

我在时间序列上使用 lm，实际上效果很好，而且超级超级快。

假设我的模型是：

> formula <- y ~ x

我在训练集上训练它：

> train <- data.frame( x = seq(1,3), y = c(2,1,4) )
> model <- lm( formula, train )

...我可以对新数据进行预测：

> test <- data.frame( x = seq(4,6) )
> test$y <- predict( model, newdata = test )
> test
  x        y
1 4 4.333333
2 5 5.333333
3 6 6.333333

这非常好用，而且速度非常快。

我想将滞后变量添加到模型中。现在，我可以通过增加我的原始训练集来做到这一点：

> train$y_1 <- c(0,train$y[1:nrow(train)-1])
> train
  x y y_1
1 1 2   0
2 2 1   2
3 3 4   1

更新公式：

formula <- y ~ x * y_1

...并且培训会很好地工作：

> model <- lm( formula, train )
> # no errors here

但是，问题在于无法使用“预测”，因为无法以批处理方式在测试集中填充 y_1。

现在，对于很多其他回归的东西，有很方便的方式在公式中表达它们，比如poly(x,2)等等，这些直接使用未经修改的训练和测试数据。

所以，我想知道是否有某种方法可以在公式中表达滞后变量，以便predict可以使用？理想情况下：

formula <- y ~ x * lag(y,-1)
model <- lm( formula, train )
test$y <- predict( model, newdata = test )

...无需增加（不确定这是否是正确的词）训练和测试数据集，并且能够predict直接使用？

4

4 回答 4

17

看看例如为您提供滞后运算符的dynlm包。更一般地说，计量经济学和时间序列的任务视图将有更多内容供您查看。

以下是其示例的开头——滞后 1 个月和 12 个月：

R>      data("UKDriverDeaths", package = "datasets")
R>      uk <- log10(UKDriverDeaths)
R>      dfm <- dynlm(uk ~ L(uk, 1) + L(uk, 12))
R>      dfm

Time series regression with "ts" data:
Start = 1970(1), End = 1984(12)

Call:
dynlm(formula = uk ~ L(uk, 1) + L(uk, 12))

Coefficients:
(Intercept)     L(uk, 1)    L(uk, 12)  
      0.183        0.431        0.511  

R>
于 2012-10-27T02:53:26.893 回答
6

按照 Dirk 的建议dynlm，我不太清楚如何预测，但搜索它导致我通过https://stats.stackexchange.com/questions/6758/1-step-ahead-predictions-with-dynlmdyn打包-r-包

然后经过几个小时的实验，我想出了以下函数来处理预测。路上有很多“问题”，例如，您似乎无法进行时间序列，并且 predict 的结果被一大堆类似的东西rbind所抵消，所以我觉得这个答案与仅命名 a 相比显着增加start包，虽然我赞成德克的回答。

因此，一个可行的解决方案是：

  • 使用dyn包
  • 使用以下方法进行预测

predictDyn 方法：

# pass in training data, test data,
# it will step through one by one
# need to give dependent var name, so that it can make this into a timeseries
predictDyn <- function( model, train, test, dependentvarname ) {
    Ntrain <- nrow(train)
    Ntest <- nrow(test)
    # can't rbind ts's apparently, so convert to numeric first
    train[,dependentvarname] <- as.numeric(train[,dependentvarname])
    test[,dependentvarname] <- as.numeric(test[,dependentvarname])
    testtraindata <- rbind( train, test )
    testtraindata[,dependentvarname] <- ts( as.numeric( testtraindata[,dependentvarname] ) )
    for( i in 1:Ntest ) {
       result <- predict(model,newdata=testtraindata,subset=1:(Ntrain+i-1))
       testtraindata[Ntrain+i,dependentvarname] <- result[Ntrain + i + 1 - start(result)][1]
    }
    return( testtraindata[(Ntrain+1):(Ntrain + Ntest),] )
}

示例用法：

library("dyn")

# size of training and test data
N <- 6
predictN <- 10

# create training data, which we can get exact fit on, so we can check the results easily
traindata <- c(1,2)
for( i in 3:N ) { traindata[i] <- 0.5 + 1.3 * traindata[i-2] + 1.7 * traindata[i-1] }
train <- data.frame( y = ts( traindata ), foo = 1)

# create testing data, bunch of NAs
test <- data.frame( y = ts( rep(NA,predictN) ), foo = 1)

# fit a model
model <- dyn$lm( y ~ lag(y,-1) + lag(y,-2), train )
# look at the model, it's a perfect fit. Nice!
print(model)

test <- predictDyn( model, train, test, "y" )
print(test)

# nice plot
plot(test$y, type='l')

输出：

> model

Call:
lm(formula = dyn(y ~ lag(y, -1) + lag(y, -2)), data = train)

Coefficients:
(Intercept)   lag(y, -1)   lag(y, -2)  
        0.5          1.7          1.3  

> test
             y foo
7     143.2054   1
8     325.6810   1
9     740.3247   1
10   1682.4373   1
11   3823.0656   1
12   8686.8801   1
13  19738.1816   1
14  44848.3528   1
15 101902.3358   1
16 231537.3296   1

编辑：嗯，这虽然超级慢。即使我将数据限制在数据集中subset的固定几行，每次预测也需要大约 24 毫秒，或者，对于我的任务，0.024*7*24*8*20*10/60/60= 1.792 hours:-O

于 2012-10-27T09:28:03.090 回答
3

试试 ARIMA 函数。AR 参数用于自回归，这意味着滞后 y。xreg = 允许您添加其他 X 变量。您可以使用 predict.ARIMA 获得预测。

于 2018-08-30T17:35:40.113 回答
1

这里有一个想法：

你为什么不创建一个新的数据框？用您需要的回归量填充数据框。对于您想要的任何变量的所有滞后，您可以使用 L1、L2、...、Lp 之类的列，然后，您可以像使用横截面类型的回归一样使用您的函数。

因为您不必每次调用拟合和预测函数时都对数据进行操作，而是将数据转换一次，因此速度会快得多。我知道 Eviews 和 Stata 提供了滞后的运营商。确实有一些方便。但是，如果您不需要像“lm”计算这样的所有功能，它也是低效的。如果您有数十万次迭代要执行，而您只需要预测，或预测和信息标准（如 BIC 或 AIC）的值，您可以通过避免进行您不会进行的计算来在速度上击败“lm”使用 -- 只需在函数中编写一个 OLS 估计器，就可以了。

于 2017-12-03T14:35:11.210 回答