使用 dplyr::do 按组拟合多个模型时处理错误的正确方法

Question

使用dplyr::do，可以非常简单地按组拟合多个模型，如下所示：

library(tidyverse)
set.seed(100)
tbl <- tibble(
  group_id = rep(1:3, each = 10),
  y1 = rnorm(30),
  y2 = runif(30),
  x1 = rnorm(30),
  x2 = runif(30)
)

tbl %>%
  group_by(group_id) %>%
  do(
    model1 = lm(y1 ~ x1 + x2, data = .),
    model2 = lm(y2 ~ x1 + x2, data = .)
  )
#> Source: local data frame [3 x 3]
#> Groups: <by row>
#> 
#> # A tibble: 3 x 3
#>   group_id model1   model2  
#> *    <int> <list>   <list>  
#> 1        1 <S3: lm> <S3: lm>
#> 2        2 <S3: lm> <S3: lm>
#> 3        3 <S3: lm> <S3: lm>

这是用于 broom::tidy 和 broom::glance 以按组提取 r.squared 和系数的理想格式。但是，当一组 group_id == 3 具有所有缺失值时会出现问题：

tbl2 <- mutate(tbl, y2 = c(runif(20), rep(NA, 10)))

tbl2 %>%
  group_by(group_id) %>%
  do(
    model1 = lm(y1 ~ x1 + x2, data = .),
    model2 = lm(y2 ~ x1 + x2, data = .)
  )
#> Error in lm.fit(x, y, offset = offset, singular.ok = singular.ok, ...): 0 (non-NA) cases

正如预期的那样，因为 group_id == 3 没有 y2 的非缺失值，model2 无法拟合任何内容。我发现的其他问题建议在拟合之前简单地删除具有 NA 值的行，但是我不想这样做，因为那样我会失去 model1 的成功拟合。我想到的另一种方法是用 try 捕获错误，但我无法只用缺失值替换错误。我尝试了以下使用 purrr::modify_if 的代码的许多变体，但不知道为什么不替换该值（例如，

modify_if(list(1, "a", TRUE), ~ inherits(., "numeric"), `is.na<-`)

工作正常。）您可以看到使用 map 和 inherits 正确地发现了哪个单元格是 class try-error，但将其包裹在 modify_if 让它不再被发现。

tbl2 %>%
  group_by(group_id) %>%
  do(
    model1 = lm(y1 ~ x1 + x2, data = .),
    model2 = try(
      lm(y2 ~ x1 + x2, data = .),
      silent = TRUE
    )
  ) %>%
  ungroup() %>%
  mutate_all(
    function(col) map_lgl(col, function(cell) inherits(cell, "try-error"))
  )
#> # A tibble: 3 x 3
#>   group_id model1 model2
#>   <lgl>    <lgl>  <lgl> 
#> 1 FALSE    FALSE  FALSE 
#> 2 FALSE    FALSE  FALSE 
#> 3 FALSE    FALSE  TRUE

tbl2 %>%
  group_by(group_id) %>%
  do(
    model1 = lm(y1 ~ x1 + x2, data = .),
    model2 = try(
      lm(y2 ~ x1 + x2, data = .),
      silent = TRUE
    )
  ) %>%
  ungroup() %>%
  mutate_at(
    .vars = vars(starts_with("model_")),
    .funs = function(col) {
      modify_if(
        .x = col,
        .p = function(cell) inherits(cell, "try-error"),
        .f = function(cell) unclass(`is.na<-`(cell)))
    }
  )
#> # A tibble: 3 x 3
#>   group_id model1   model2         
#> *    <int> <list>   <list>         
#> 1        1 <S3: lm> <S3: lm>       
#> 2        2 <S3: lm> <S3: lm>       
#> 3        3 <S3: lm> <S3: try-error>

由 reprex package (v0.2.0) 创建于 2018-04-17。

我的实际数据有~80k组和~10个模型供参考。非常感谢任何改进此代码的建议或更好的捕获错误的方法。

Answer 1

我认为这是我找到的处理此问题的最佳方法。与其使用 modify 尝试替换错误模型，不如将它们过滤掉并替换 glance 之后缺失的行。这是因为 glance 无论如何都不能很好地处理格式错误的 lm 输出。

tbl2 %>%
  group_by(group_id) %>%
  do(
    model1 = lm(y1 ~ x1 + x2, data = .),
    model2 = try(
      lm(y2 ~ x1 + x2, data = .),
      silent = TRUE
    )
  ) %>%
  ungroup() %>%
  gather(model, lm, starts_with("model")) %>%
  mutate(error = map_lgl(lm, ~inherits(., "try-error"))) %>%
  filter(error == FALSE) %>%
  rowwise() %>%
  glance(lm) %>%
  ungroup() %>%
  complete(group_id = 1:3, model = c("model1", "model2"))
#> # A tibble: 6 x 14
#>   group_id model  error r.squared adj.r.squared  sigma statistic p.value
#>      <int> <chr>  <lgl>     <dbl>         <dbl>  <dbl>     <dbl>   <dbl>
#> 1        1 model1 FALSE    0.0215       -0.258   0.629    0.0769   0.927
#> 2        1 model2 FALSE    0.107        -0.149   0.329    0.418    0.674
#> 3        2 model1 FALSE    0.208        -0.0184  0.868    0.919    0.442
#> 4        2 model2 FALSE    0.0808       -0.182   0.362    0.308    0.745
#> 5        3 model1 FALSE    0.0707       -0.195   0.738    0.266    0.774
#> 6        3 model2 NA      NA            NA      NA       NA       NA    
#> # ... with 6 more variables: df <int>, logLik <dbl>, AIC <dbl>, BIC <dbl>,
#> #   deviance <dbl>, df.residual <int>

使用 dplyr::do 按组拟合多个模型时处理错误的正确方法

Proper way of dealing with errors when fitting many models by group with dplyr::do

r

lm

dplyr