在R中,线性模型的系数和摘要有不同的NA操作。

11
在R中,当使用lm()时,如果我在调用lm()时设置na.action = na.pass,那么在摘要表中,任何无法估计系数的地方(因为在这种情况下缺少单元格)都会有一个NA。
然而,如果我仅从摘要对象中提取系数,使用summary(myModel)$coefficientscoef(summary(myModel)),那么NA将被省略。
我希望在提取系数时包含NA,就像在打印摘要时包含NA一样。有没有办法做到这一点?
设置options(na.action = na.pass)似乎没有帮助。
以下是示例:
> set.seed(534)
> myGroup1 <- factor(c("a","a","a","a","b","b"))
> myGroup2 <- factor(c("first","second","first","second","first","first"))
> myDepVar <- rnorm(6, 0, 1)
> myModel <- lm(myDepVar ~ myGroup1 + myGroup2 + myGroup1:myGroup2)
> summary(myModel)

Call:
lm(formula = myDepVar ~ myGroup1 + myGroup2 + myGroup1:myGroup2)

Residuals:
       1        2        3        4        5        6 
-0.05813  0.55323  0.05813 -0.55323 -0.12192  0.12192 

Coefficients: (1 not defined because of singularities)
                    Estimate Std. Error t value Pr(>|t|)
(Intercept)         -0.15150    0.23249  -0.652    0.561
myGroup11            0.03927    0.23249   0.169    0.877
myGroup21           -0.37273    0.23249  -1.603    0.207
myGroup11:myGroup21       NA         NA      NA       NA

Residual standard error: 0.465 on 3 degrees of freedom
Multiple R-squared: 0.5605,     Adjusted R-squared: 0.2675 
F-statistic: 1.913 on 2 and 3 DF,  p-value: 0.2914 

> coef(summary(myModel))
               Estimate Std. Error    t value  Pr(>|t|)
(Intercept) -0.15149826  0.2324894 -0.6516352 0.5611052
myGroup11    0.03926774  0.2324894  0.1689012 0.8766203
myGroup21   -0.37273117  0.2324894 -1.6032180 0.2072173

> summary(myModel)$coefficients
               Estimate Std. Error    t value  Pr(>|t|)
(Intercept) -0.15149826  0.2324894 -0.6516352 0.5611052
myGroup11    0.03926774  0.2324894  0.1689012 0.8766203
myGroup21   -0.37273117  0.2324894 -1.6032180 0.2072173

你认为这是一个 bug 吗? - randy
3个回答

3
为什么不直接从拟合的模型中提取系数:
> coef(myModel)
             (Intercept)                myGroup1b 
             -0.48496169              -0.07853547 
          myGroup2second myGroup1b:myGroup2second 
              0.74546233                       NA

这似乎是最简单的选择。

na.action与此无关。请注意,在您的示例中,没有传递na.action = na.pass

na.action是用于处理模型拟合时传递给模型的数据中的NA的全局选项,通常与公式一起使用;它也是一个名为na.action()的函数。 R从data参数和在公式中表达的模型的符号表示中建立所谓的模型框架。此时,任何NA都将被检测到,na.action的默认选项是使用na.omit()通过删除具有任何变量的NA的样本来从数据中删除NA。有其他选择,最有用的是na.exclude(),它会在拟合过程中删除NA,但会将NA添加回适当的位置,如拟合值、残差等。阅读?na.omit?na.action了解更多信息,以及?options了解更多信息。


3
感谢您解释na.action设置对此问题无关。如果万不得已,从拟合模型中提取系数可能会起作用,但我想将几列置于汇总表格以供置信区间使用。我不仅需要估计值;我还需要标准误差、p值等,并在末尾附加置信区间。我可以从头开始制作表格,但我认为可能需要更改某些简单的设置才能使coef(summary(myModel))confint(myModel)按照相同的顺序输出相同数量的行。 - Jdub
@Jdub,你解决了这个问题吗?我也遇到了完全相同的问题。 - half-pass
同样在这里!同样的问题。 - vagabond
@Jdub,这只是一个简单的问题,只需要输入:summary(model)[coef(model), ],其中 [i 参数中的 NA 会生成一行完全为 NA 的结果。我希望这就是你要求的,因为这是唯一有意义的输出。否则,你可能需要更好地描述一下你想做什么。 - AdamO

1

summary.lm的文档中写道:“别名系数在返回对象中被省略,但通过print方法恢复。”似乎没有参数来控制此省略。除了使用@Gavin Simpson建议的coef(summary(myModel))之外,还有另一种解决方法。您可以创建一个矩阵。

nr <- num_regressors - nrow(summary(myModel)$coefficients) ##num_regressors shall be defined previously
nc <- 4
rnames <- names(which(summary(myModel)$aliased))
cnames <- colnames(summary(myModel)$coefficients)
mat_na <- matrix(data = NA,nrow = nr,ncol = nc,
           dimnames = list(rnames,cnames))

然后将这两个矩阵按行合并:
mat_coef <- rbind(summary(myModel)$coefficients,mat_na)

-1
你也可以将摘要适配表格转换为数据框(其中缺失的变量会丢失):
fit <- as.data.frame(summary(fit)$coefficients)

然后按名称提取系数:

fit["age", "Pr(>|z|)"]

如果“age”已被删除,则在尝试从数据框中提取年龄的P值时会得到NA。

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接