1

我正在尝试基于此处的示例构建滚动回归函数，但除了返回预测值之外，我还想返回一些滚动模型诊断（即系数、t 值和 mabye R^2）。我希望根据结果类型以离散对象的形式返回结果。上面链接中提供的示例成功地创建了滚动预测，但我需要一些帮助打包和写出滚动模型诊断：

最后，我希望该函数返回三 (3) 个对象：

  1. 预测
  2. 系数
  3. T 值
  4. R^2

下面是代码：

require(zoo)
require(dynlm)

## Create Some Dummy Data
set.seed(12345)
x <- rnorm(mean=3,sd=2,100)
y <- rep(NA,100)
y[1] <- x[1]
for(i in 2:100) y[i]=1+x[i-1]+0.5*y[i-1]+rnorm(1,0,0.5)
int <- 1:100
dummydata <- data.frame(int=int,x=x,y=y)
zoodata <- as.zoo(dummydata)


rolling.regression <- function(series) {
  mod <- dynlm(formula = y ~ L(y) + L(x), data = as.zoo(series)) # get model

  nextOb <- max(series[,'int'])+1 # To get the first row that follows the window
  if (nextOb<=nrow(zoodata)) {   # You won't predict the last one

    # 1) Make Predictions
    predicted <- predict(mod,newdata=data.frame(x=zoodata[nextOb,'x'],y=zoodata[nextOb,'y']))
    attributes(predicted) <- NULL
    c(predicted=predicted,square.res <-(predicted-zoodata[nextOb,'y'])^2)

    # 2) Extract coefficients
    #coefficients <- coef(mod)

    # 3) Extract rolling coefficient t values
    #tvalues <- ????(mod)

    # 4) Extract rolling R^2
    #rsq <-


  }
}    

rolling.window <- 20
results.z <-  rollapply(zoodata, width=rolling.window, FUN=rolling.regression, by.column=F, align='right')

因此，在弄清楚如何从模型（即 mod）中提取 t 值之后，我需要做什么才能使函数返回三 (3) 个单独的对象（即预测、系数和 T 值）？

我对 R 相当陌生，对功能真的很陌生，对动物园也很陌生，我被困住了。

任何帮助将不胜感激。

4

1 回答 1

2

我希望我能正确理解您，但这里是您功能的一个小编辑：

rolling.regression <- function(series) {
  mod <- dynlm(formula = y ~ L(y) + L(x), data = as.zoo(series)) # get model

  nextOb <- max(series[,'int'])+1 # To get the first row that follows the window
  if (nextOb<=nrow(zoodata)) {   # You won't predict the last one
    # 1) Make Predictions
    predicted=predict(mod,newdata=data.frame(x=zoodata[nextOb,'x'],y=zoodata[nextOb,'y']))
    attributes(predicted)<-NULL
    #Solution 1; Quicker to write
    #     c(predicted=predicted, 
    #       square.res=(predicted-zoodata[nextOb,'y'])^2,
    #       summary(mod)$coef[, 1],
    #       summary(mod)$coef[, 3],
    #       AdjR = summary(mod)$adj.r.squared)

    #Solution 2; Get column names right
    c(predicted=predicted, 
      square.res=(predicted-zoodata[nextOb,'y'])^2,
      coef_intercept = summary(mod)$coef[1, 1],
      coef_Ly = summary(mod)$coef[2, 1],
      coef_Lx = summary(mod)$coef[3, 1],
      tValue_intercept = summary(mod)$coef[1, 3],
      tValue_Ly = summary(mod)$coef[2, 3],
      tValue_Lx = summary(mod)$coef[3, 3],
      AdjR = summary(mod)$adj.r.squared)
  }
}



rolling.window <- 20
results.z <-  rollapply(zoodata, width=rolling.window, FUN=rolling.regression, by.column=F, align='right')

    head(results.z)
   predicted square.res coef_intercept   coef_Ly  coef_Lx tValue_intercept tValue_Ly tValue_Lx      AdjR
20 10.849344   0.721452     0.26596465 0.5798046 1.049594       0.38309211  7.977627  13.59831 0.9140886
21 12.978791   2.713053     0.26262820 0.5796883 1.039882       0.37741499  7.993014  13.80632 0.9190757
22  9.814676  11.719999     0.08050796 0.5964808 1.073941       0.12523824  8.888657  15.01353 0.9340732
23  5.616781  15.013297     0.05084124 0.5984748 1.077133       0.08964998  9.881614  16.48967 0.9509550
24  3.763645   6.976454     0.26466039 0.5788949 1.068493       0.51810115 11.558724  17.22875 0.9542983
25  9.433157  31.772658     0.38577698 0.5812665 1.034862       0.70969330 10.728395  16.88175 0.9511061

要了解它是如何工作的，请用回归做一个小例子：

x <- rnorm(1000); y <- 2*x + rnorm(1000)
reg <- lm(y ~ x)
summary(reg)$coef
              Estimate Std. Error    t value Pr(>|t|)
(Intercept) 0.02694322 0.03035502  0.8876033 0.374968
x           1.97572544 0.03177346 62.1816310 0.000000

如您所见，summary首先调用然后获取它的系数（coef(summary(reg))也可以）为您提供一个包含估计值、标准误差和 t 值的表格。因此估计值保存在该表的第 1 列中，t 值保存在第 3 列中。这就是我在更新rolling.regression函数中获取它们的方式。

编辑

我更新了我的解决方案；现在它还包含调整后的 R2。如果您只想要普通的 R2，请摆脱.adj.

编辑 2

快速而肮脏的 hack 如何命名列：

rolling.regression <- function(series) {
  mod <- dynlm(formula = y ~ L(y) + L(x), data = as.zoo(series)) # get model

  nextOb <- max(series[,'int'])+1 # To get the first row that follows the window
  if (nextOb<=nrow(zoodata)) {   # You won't predict the last one
    # 1) Make Predictions
    predicted=predict(mod,newdata=data.frame(x=zoodata[nextOb,'x'],y=zoodata[nextOb,'y']))
    attributes(predicted)<-NULL
    #Get variable names
    strVar <- c("Intercept", paste0("L", 1:(nrow(summary(mod)$coef)-1)))
    vec <- c(predicted=predicted, 
             square.res=(predicted-zoodata[nextOb,'y'])^2,
             AdjR = summary(mod)$adj.r.squared,
             summary(mod)$coef[, 1],
             summary(mod)$coef[, 3])
    names(vec)[4:length(vec)] <- c(paste0("Coef_", strVar), paste0("tValue_", strVar))

    vec
  }
}
于 2013-02-06T17:04:50.490 回答