“mlm”对象的“stats::confint”置信区间错误？答案

【问题标题】：wrong confidence intervals from `stats::confint` for `mlm` objects?“mlm”对象的“stats::confint”置信区间错误？
【发布时间】：2020-02-22 10:33:23
【问题描述】：

我不确定我是否做错了什么，或者这是 R 本身中的 confint 函数中的错误，但我得到的回归估计置信区间不包含估计值。

这里是reprex：

# model (converting all numeric columns in data to z-scores)
mod <-
  stats::lm(
    formula = cbind(mpg, disp) ~ wt,
    data = purrr::modify_if(.x = mtcars, .p = is.numeric, .f = scale)
  )

# tidy dataframe
broom::tidy(mod, conf.int = TRUE)
#> # A tibble: 4 x 8
#>   response term         estimate std.error statistic  p.value conf.low conf.high
#>   <chr>    <chr>           <dbl>     <dbl>     <dbl>    <dbl>    <dbl>     <dbl>
#> 1 mpg      (Intercept)  8.05e-17    0.0893  9.01e-16 1.00e+ 0   -0.182     0.182
#> 2 mpg      wt          -8.68e- 1    0.0908 -9.56e+ 0 1.29e-10   -1.05     -0.682
#> 3 disp     (Intercept) -1.70e-16    0.0826 -2.06e-15 1.00e+ 0   -0.182     0.182
#> 4 disp     wt           8.88e- 1    0.0840  1.06e+ 1 1.22e-11   -1.05     -0.682


confint(mod)
#>                   2.5 %     97.5 %
#> :(Intercept) -0.1824544  0.1824544
#> :wt          -1.0530332 -0.6822855
#> :(Intercept) -0.1824544  0.1824544
#> :wt          -1.0530332 -0.6822855

如果您绘制估计值，这会变得更加明显：

我在这里做错了吗？或者这是预期的行为？

【问题讨论】：

确实，间隔看起来完全错误。见stackoverflow.com/questions/28442141/…

标签： r statistics

【解决方案1】：

confint() 现在可以与 mlm 一起使用：

confint(lm(mpg~wt,data=mtcars))
                2.5 %    97.5 %
(Intercept) 33.450500 41.119753
wt          -6.486308 -4.202635

confint(lm(disp~wt,data=mtcars))
                2.5 %   97.5 %
(Intercept) -204.0914 -58.2054
wt            90.7579 134.1984

confint(lm(cbind(mpg,disp)~wt,data=mtcars))
                       2.5 %     97.5 %
mpg:(Intercept)    33.450500  41.119753
mpg:wt             -6.486308  -4.202635
disp:(Intercept) -204.091436 -58.205395
disp:wt            90.757897 134.198380

问题在于在 dplyr 或 purrr 中使用 scale，我在 scale 后使用 as.data.frame 因为你的所有列都是数字：

confint(lm(cbind(mpg,disp)~wt,data=as.data.frame(scale(mtcars))))

                          2.5 %     97.5 %
    mpg:(Intercept)  -0.1824544  0.1824544
    mpg:wt           -1.0530332 -0.6822855
    disp:(Intercept) -0.1687740  0.1687740
    disp:wt           0.7165054  1.0594545

对于 dplyr 和 purrr，你可以看到变量名也乱了：

library(dplyr)
library(purrr)

confint(lm(cbind(mpg,disp)~wt,data=modify_if(.x = mtcars, .p = is.numeric, .f = scale)))

                  2.5 %     97.5 %
:(Intercept) -0.1824544  0.1824544
:wt          -1.0530332 -0.6822855
:(Intercept) -0.1824544  0.1824544
:wt          -1.0530332 -0.6822855

confint(lm(cbind(mpg,disp)~wt,data=mutate_if(mtcars,is.numeric,scale))

                  2.5 %     97.5 %
:(Intercept) -0.1824544  0.1824544
:wt          -1.0530332 -0.6822855
:(Intercept) -0.1824544  0.1824544
:wt          -1.0530332 -0.6822855

我猜当你进行缩放时，继承的属性会与 confint 中的某些内容混淆（仍在查看代码。）解决此问题的一种方法是将输出强制为向量（使用 c 或 as.数字）：

confint(lm(cbind(mpg,disp)~wt,data=modify_if(mtcars,is.numeric,~c(scale(.x)))))

或

confint(lm(cbind(mpg,disp)~wt,data=mutate_if(mtcars,is.numeric,~c(scale(.x)))))

                      2.5 %     97.5 %
mpg:(Intercept)  -0.1824544  0.1824544
mpg:wt           -1.0530332 -0.6822855
disp:(Intercept) -0.1687740  0.1687740
disp:wt           0.7165054  1.0594545

不幸的是，tidy(..,conf.int=TRUE) 并不能很好地工作，因此要进行绘图，您需要使用统计数据中的 confint 并进行一些操作。

【讨论】：