代码之家 › 专栏 › 技术社区 › Coolio2654

自定义函数[closed]上的“下标越界”错误

lm function r

Coolio2654 · 技术社区 · 7 年前

我决定尝试手工编写一个简单的自定义函数,在一些 lm()

这是我第一次创建自己的代码,但我已经用R工作了几个月,我认为我对它有相当的理解,所以我不明白为什么我在运行它时总是得到一个“下标超出界限”。

我的代码:

custom_test<-function(data,coeff,alt,alternative=c("two.sided","greater","less"),clevel=.95){
  dof<-data$df.residual
  top<-data$coefficients["coeff"]-alt
  bottom=coef(summary(data))["coeff","Std. Error"]
  stat<-abs(top/bottom)
  if (alternative=="two.sided") {
    tstat<-qt(clevel/2,dof)
    pstat<-2*pt(tstat,dof)
    return(pstat)
  } else if (alternative=="greater") {
      tstat<-qt(clevel/2,dof)
      pstat<-pt(tstat,dof)
      return(pstat)
  } else if (alternative=="less") {
      tstat<-qt(clevel/2,dof)
      pstat<-pt(tstat,dof)
      return(pstat)
  } else {
      return("Error")
  }

}

lm() hrsemp 作为变量,并获取错误:

custom_test(fit9,hrsemp,0,alternative="less")
Error in coef(summary(data))["coeff", "Std. Error"] : 
  subscript out of bounds

但每次我自己手动运行有问题的代码时,我都会得到一个答案:

> coef(fit9)
(Intercept)      hrsemp  log(sales) log(employ) 
12.45837237 -0.02926893 -0.96202698  0.76147045 
> coef(summary(fit9))["hrsemp", "Std. Error"]
[1] 0.02280484

关于这个错误的其他堆栈交换问题似乎都有细微的不同,到目前为止,我还无法将它们的经验总结到我的代码中。

1 回复 | 直到 5 年前

duckmayr 7 年前

Frank 是正确的;您在base遇到这个错误的原因与其他人相同:您试图访问一个不存在的对象的元素。更具体地说,在您的示例中,您试图访问 "coeff" 行和 "Std. Error" 第列,共列 coef(summary(data)) . 您要执行以下操作:

custom_test<-function(data,coeff,alt,alternative=c("two.sided","greater","less"),clevel=.95){
    dof<-data$df.residual
    top<-data$coefficients[coeff]-alt
    bottom=coef(summary(data))[coeff,"Std. Error"]
    stat<-abs(top/bottom)
    if (alternative=="two.sided") {
        tstat<-qt(clevel/2,dof)
        pstat<-2*pt(tstat,dof)
        return(pstat)
    } else if (alternative=="greater") {
        tstat<-qt(clevel/2,dof)
        pstat<-pt(tstat,dof)
        return(pstat)
    } else if (alternative=="less") {
        tstat<-qt(clevel/2,dof)
        pstat<-pt(tstat,dof)
        return(pstat)
    } else {
        return("Error")
    }

}

并将变量名作为字符串传递:

set.seed(42)
hrsemp <- rnorm(10)
Y <- 1 + 5 * hrsemp + rnorm(10)
fit9 <- lm(Y ~ hrsemp)
custom_test(fit9, 'hrsemp', 0, alternative="less")
[1] 0.475

deparse(substitute(coeff)) --例如,请参见 this SO question

然而,你可能会注意到这给了你错误的答案。这是因为您编写的函数不正确。你可能想要这样的东西:

custom_test <- function(data, coeff, alt,
                        alternative = c("two.sided", "greater", "less"),
                        clevel = .95){
    dof <- data$df.residual
    top <- data$coefficients[coeff] - alt
    bottom <- coef(summary(data))[coeff, "Std. Error"]
    stat <- abs(top/bottom)
    if ( alternative == "two.sided" ) {
        return(2 * (1 - pt(stat, dof)))
    } else if ( alternative == "greater" ) {
        return(1 - pt(stat, dof))
    } else if ( alternative == "less" ) {
        return(1 - pt(stat, dof))
    } else {
        stop("Provide a valid alternative hypothesis.", call.=FALSE)
    }
}


custom_test(fit9, 'hrsemp', 0, alternative="less")
      hrsemp 
7.858176e-05 
custom_test(fit9, 'hrsemp', 0, alternative="two.sided")
      hrsemp 
0.0001571635 
coef(summary(fit9))['hrsemp', 'Pr(>|t|)']
[1] 0.0001571635

可以找到许多很好的解释来解释为什么这是正确的计算 here .