R中线性模型的系数和摘要的不同NA动作

Posted 2023-02-15

【中文标题】R中线性模型的系数和摘要的不同NA动作

【英文标题】：Different NA actions for coefficients and summary of linear model in R

【发布时间】：2012-06-11 15:43:19

【问题描述】：

在 R 中，当使用 lm() 时，如果我在对 lm() 的调用中设置 na.action = na.pass，那么在汇总表中，对于任何无法估计的系数都有一个 NA（因为在这种情况下缺少单元格）。

但是，如果我仅使用 summary(myModel)$coefficients 或 coef(summary(myModel)) 从汇总对象中提取系数，则 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

【问题讨论】：

你同意这是一个错误吗？

【参考方案1】：

为什么不直接从拟合模型中提取系数：

> 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 了解更多信息，以及@987654338 @了解更多信息。

【讨论】：

感谢您解释 na.action 设置与此问题无关。从拟合模型中提取系数可能是最后的手段，但我想将置信区间的几列绑定到汇总表。我不只想要估计；我想要标准误差、p 值等，最后附上置信区间。我可以从头开始制作表格，但我认为可能需要更改一些简单的设置以使 coef(summary(myModel)) 和 confint(myModel) 以相同的顺序输出相同数量的行。

@Jdub，你发现了吗？我有完全相同的问题。

这里也一样！同样的问题

@Jdub 这不只是说：summary(model)[coef(model), ] 其中i 的[ 的NA 参数产生一个完全NA 行？我希望这就是您所要求的，因为这是唯一对我有意义的输出。否则，您可能需要更好地描述您正在尝试做的事情。

【参考方案2】：

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))

然后rbind这两个矩阵：

mat_coef <- rbind(summary(myModel)$coefficients,mat_na)

【参考方案3】：

您也可以将汇总拟合表转换为数据框（其中 NA 的变量会丢失）：

fit <- as.data.frame(summary(fit)$coefficients)

然后按名称提取系数：

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

如果“年龄”已被删除，当您尝试从数据帧中提取年龄的 P 值时，您将得到一个 NA