欧洲防风草与基 R 的预测标准误差很大

Large standard error of prediction from parsnip vs base R

似乎 predict 产生了一个太大的标准错误。我用 parsnip 模型得到 0.820,但用基本 R 模型得到 0.194。标准误差 0.194 似乎更合理,因为在我的预测上下约 2*0.195 是置信区间的终点。我的 problem/misunderstanding 是什么?

library(parsnip)
library(dplyr)

# example data
mod_dat <- mtcars %>%
  as_tibble() %>%
  mutate(cyl_8 = as.numeric(cyl == 8)) %>%
  select(mpg, cyl_8)

parsnip_mod <- logistic_reg() %>%
  set_engine("glm") %>%
  fit(as.factor(cyl_8) ~ mpg, data = mod_dat)

base_mod <- glm(as.factor(cyl_8) ~ mpg, data = mod_dat, family = "binomial")

parsnip_pred <- tibble(mpg = 18) %>%
  bind_cols(predict(parsnip_mod, new_data = ., type = 'prob'),
            predict(parsnip_mod, new_data = ., type = 'conf_int', std_error = T)) %>%
  select(!ends_with("_0"))

base_pred <- predict(base_mod, tibble(mpg = 18), se.fit = T, type = "response") %>%
  unlist()

# these give the same prediction but different SE
parsnip_pred
#> # A tibble: 1 x 5
#>     mpg .pred_1 .pred_lower_1 .pred_upper_1 .std_error
#>   <dbl>   <dbl>         <dbl>         <dbl>      <dbl>
#> 1    18   0.614         0.230         0.895      0.820
base_pred
#>          fit.1       se.fit.1 residual.scale 
#>      0.6140551      0.1942435      1.0000000

reprex package (v0.3.0)

于 2020-06-04 创建

--编辑--

正如@thelatemail 和@Limey 所说,对基本模型使用 type="link" 将给出 logit 量表 (0.820) 的标准误差。但是,我想要概率标度上的标准误差。 parsnip documentation 中是否有我遗漏的选项?我想使用 parsnip.

@thelatemail 是正确的。来自 predict.glm 的在线文档:

类型
所需的预测类型。默认值在线性预测变量的范围内;备选方案 "response" 在响应变量的范围内。因此,对于默认的二项式模型,默认预测是对数赔率(logit 尺度上的概率)并且 type = "response" 给出预测概率。

默认使用 logit 量表报告,'response' 请求原始概率量表的结果。从 parsnip::predict documentation 中并不明显,我发现它是如何选择 return 结果的尺度的,但很明显它使用的是原始概率尺度。

所以这两种方法都是 return 正确答案,只是使用了不同的尺度。

我不想从@thelatemail 窃取公认的解决方案,因此请他们post 对此做出类似的回答。

正如@thelatemail 所说,您可以使用参数 type="raw", opts=list(se.fit=TRUE, type="response") 获得 parsnip 概率标度上的标准误差。但到那时,您还不如使用基本模型,因为输出完全相同。但是,如果您已经在使用 parsnip 模型并且想要基本模型的标准错误输出,这仍然有用。

library(parsnip)
library(dplyr)

mod_dat <- mtcars %>%
  as_tibble() %>%
  mutate(cyl_8 = as.numeric(cyl == 8)) %>%
  select(mpg, cyl_8)

parsnip_mod <- logistic_reg() %>%
  set_engine("glm") %>%
  fit(as.factor(cyl_8) ~ mpg, data = mod_dat)

base_mod <- glm(as.factor(cyl_8) ~ mpg, data = mod_dat, family = "binomial")

predict(parsnip_mod, tibble(mpg = 18), type="raw",
        opts=list(se.fit=TRUE, type="response")) %>% 
  as_tibble()
#> # A tibble: 1 x 3
#>     fit se.fit residual.scale
#>   <dbl>  <dbl>          <dbl>
#> 1 0.614  0.194              1

predict.glm(base_mod, tibble(mpg = 18), se.fit = T, type="response") %>% 
  as_tibble()
#> # A tibble: 1 x 3
#>     fit se.fit residual.scale
#>   <dbl>  <dbl>          <dbl>
#> 1 0.614  0.194              1

reprex package (v0.3.0)

于 2020 年 6 月 11 日创建