Как эффективно запустить множество логистических регрессий в R и пропустить уравнения, которые выдают ошибки?

В качестве продолжения этого вопрос я хочу запустить сразу несколько уравнений логистической регрессии, а затем отметить, значительно ли группа отличается от контрольной группы. Это решение работает, но оно работает только тогда, когда я не пропускаю значения. Поскольку мои данные состоят из 100 уравнений, они обязательно будут иметь пропущенные значения, поэтому вместо того, чтобы это решение терпело неудачу при обнаружении ошибки, как я могу запрограммировать его на пропуск экземпляров, которые вызывают ошибку?

Вот измененный набор данных, в котором отсутствуют случаи:

library(dplyr)
library(lubridate)
library(broom)

test <- tibble(major = as.factor(c(rep(c("undeclared", "computer science", "english"), 2), "undeclared")),
               time = ymd(c(rep("'2021-01-01", 3), rep("'2020-01-01", 3), rep("'2019-01-01", 1))),
               admit = c(500, 1000, 450, 800, 300, 100, 1000),
               reject = c(1000, 300, 1000, 210, 100, 900, 1500)) %>%
  mutate(total = rowSums(test[ , c("admit", "reject")], na.rm=TRUE)) %>%
  mutate(accept_rate = admit/total)

А вот решение, которое работает, когда в нем есть все случаи (см. набор данных здесь), но когда оно попадает в группу 2019, в которой отсутствуют случаи, оно терпит неудачу:

library(dplyr)
library(lubridate)
library(broom)
library(tidyr)
library(purrr)

test %>% 
  # create year column
  mutate(year = year(time), 
         major = relevel(major, "undeclared")) %>% 
  
  # nest by year
  nest(data = -year) %>% 
  
  # compute regression
  mutate(reg = map(data, ~glm(accept_rate ~ major, data = ., 
                              family = binomial, weights = total, na.action = na.exclude)), 
         
         # use broom::tidy to make a tibble out of model object
         reg_tidy = map(reg, tidy)) %>% 
  
  # get data and regression results back to tibble form
  unnest(c(data, reg_tidy)) %>% 
  filter(term != "(Intercept)") %>%
  
  # create the significant yes/no column
  mutate(significant = ifelse(p.value < 0.05, "Yes", "No")) %>% 
  
  # remove the unnecessary columns
  select(-c(term, estimate, std.error, statistic, p.value, reg))

Я также попытался обернуть решение с помощью пользовательских функций здесь, но я также не смог заставить его работать. Наконец, я также открыт для других идей решения, если оно дает аналогичный результат и устойчиво к этим ошибкам.

Стоит ли изучать PHP в 2023-2024 годах?
Стоит ли изучать PHP в 2023-2024 годах?
Привет всем, сегодня я хочу высказать свои соображения по поводу вопроса, который я уже много раз получал в своем сообществе: "Стоит ли изучать PHP в...
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
Поведение ключевого слова "this" в стрелочной функции в сравнении с нормальной функцией
В JavaScript одним из самых запутанных понятий является поведение ключевого слова "this" в стрелочной и обычной функциях.
Приемы CSS-макетирования - floats и Flexbox
Приемы CSS-макетирования - floats и Flexbox
Здравствуйте, друзья-студенты! Готовы совершенствовать свои навыки веб-дизайна? Сегодня в нашем путешествии мы рассмотрим приемы CSS-верстки - в...
Тестирование функциональных ngrx-эффектов в Angular 16 с помощью Jest
В системе управления состояниями ngrx, совместимой с Angular 16, появились функциональные эффекты. Это здорово и делает код определенно легче для...
Концепция локализации и ее применение в приложениях React ⚡️
Концепция локализации и ее применение в приложениях React ⚡️
Локализация - это процесс адаптации приложения к различным языкам и культурным требованиям. Это позволяет пользователям получить опыт, соответствующий...
Пользовательский скаляр GraphQL
Пользовательский скаляр GraphQL
Листовые узлы системы типов GraphQL называются скалярами. Достигнув скалярного типа, невозможно спуститься дальше по иерархии типов. Скалярный тип...
1
0
66
2
Перейти к ответу Данный вопрос помечен как решенный

Ответы 2

purrr::safely позволяет позаботиться об ошибках. Чтобы обернуть вызов glm внутри purrr::safely, я использую вспомогательную функцию glm_safe. glm_safe возвращает список с двумя элементами, result и error.

Когда ошибки нет, result содержит объект модели, а element — это NULL. В случае ошибки сообщение об ошибке сохраняется в error, а result равно NULL.

Чтобы использовать результаты в вашем конвейере, мы должны извлечь элементы result, которые можно получить с помощью transpose(reg)$result.

library(dplyr)
library(lubridate)
library(broom)
library(tidyr)
library(purrr)

test <- tibble(
  major = as.factor(c(rep(c("undeclared", "computer science", "english"), 2), "undeclared")),
  time = ymd(c(rep("'2021-01-01", 3), rep("'2020-01-01", 3), rep("'2019-01-01", 1))),
  admit = c(500, 1000, 450, 800, 300, 100, 1000),
  reject = c(1000, 300, 1000, 210, 100, 900, 1500)
)

test <- test %>%
  mutate(total = rowSums(test[, c("admit", "reject")], na.rm = TRUE)) %>%
  mutate(accept_rate = admit / total)

glm_safe <- purrr::safely(
  function(x) {
    glm(accept_rate ~ major,
      data = x,
      family = binomial, weights = total, na.action = na.exclude
    )
  }
)

test %>%
  # create year column
  mutate(
    year = year(time),
    major = relevel(major, "undeclared")
  ) %>%
  # nest by year
  nest(data = -year) %>%
  # compute regression
  mutate(reg = map(data, glm_safe),
         reg = transpose(reg)$result) |> 
  mutate(reg_tidy = map(reg, tidy)) %>%
  # get data and regression results back to tibble form
  unnest(c(data, reg_tidy)) %>%
  filter(term != "(Intercept)") %>%
  # create the significant yes/no column
  mutate(significant = ifelse(p.value < 0.05, "Yes", "No")) %>%
  # remove the unnecessary columns
  select(-c(term, estimate, std.error, statistic, p.value, reg))
#> # A tibble: 4 × 8
#>    year major            time       admit reject total accept_rate significant
#>   <dbl> <fct>            <date>     <dbl>  <dbl> <dbl>       <dbl> <chr>      
#> 1  2021 computer science 2021-01-01  1000    300  1300       0.769 Yes        
#> 2  2021 english          2021-01-01   450   1000  1450       0.310 No         
#> 3  2020 computer science 2020-01-01   300    100   400       0.75  No         
#> 4  2020 english          2020-01-01   100    900  1000       0.1   Yes

Спасибо за вашу помощь! Когда я пытаюсь запустить это, я получаю сообщение об ошибке: Ошибка: Проблема со столбцом mutate()reg. ℹ reg = transpose(reg)$result. x Элемент 1 ввода списка не является атомарным вектором

J.Sabree 04.04.2022 23:16

Хм. Я только что внес изменения, чтобы добавить полный репрекс, включая данные. Который, как вы могли видеть, работает без ошибок. Может быть, вы могли бы попробовать еще раз.

stefan 04.04.2022 23:19

да, он все еще дает мне ошибку даже после обновления. Возможно, это что-то с моим компьютером. Тем не менее, спасибо, что нашли время, чтобы представить это!

J.Sabree 04.04.2022 23:35
Ответ принят как подходящий

Чтобы игнорировать ошибки, используйте эту функцию:

get_model <- function(df) {
  tryCatch(
    glm(accept_rate ~ major, data = df, family = binomial, weights = total, na.action = na.exclude),
    error = function(e) NULL, warning=function(w) NULL)
}

Используйте его там, где вы вызываете mutate(reg=map()...):

  # compute regression
  mutate(reg = map(data, get_model), reg_tidy = map(reg, tidy))

Выход:

# A tibble: 4 x 8
   year major            time       admit reject total accept_rate significant
  <dbl> <fct>            <date>     <dbl>  <dbl> <dbl>       <dbl> <chr>      
1  2021 computer science 2021-01-01  1000    300  1300       0.769 Yes        
2  2021 english          2021-01-01   450   1000  1450       0.310 No         
3  2020 computer science 2020-01-01   300    100   400       0.75  No         
4  2020 english          2020-01-01   100    900  1000       0.1   Yes

Другие вопросы по теме