daaronr
daaronr

Reputation: 517

How to write a function to use lapply or purrr to broom::tidy a list of (polr) model outputs?

I'm running a list of ordered logit models with different variables etc. I want to convert the output into a tidy tibble to use in ggplot etc. (I also want to save the 'regular model output' so I want to do this separately.)

I want to do this in an automated way, using purrr or lapply or some such, to be able to first 'run all the models' (automating that is another question for later) and then 'tidy all the models', the latter presumably generating a list of tibbles.

I've tried the following, but it throws: Error: No tidy method recognized for this list.

clean_model <- function(polr_results) {
  lapply(polr_results,  
    broom::tidy(polr_results, conf.int = TRUE, exponentiate = TRUE) %>%
      filter(coef.type=="coefficient")  %>% 
      dplyr::arrange(-str_detect(term, 'd2sd'))
      )
}

mtcars_m1 <- mtcars %>% polr(as.factor(cyl) ~ hp , data = ., Hess = TRUE) 
mtcars_m2 <- mtcars %>% polr(as.factor(cyl) ~ hp + qsec , data = ., Hess = TRUE) 

clean_model(c(mtcars_m1, mtcars_m2))

Upvotes: 0

Views: 500

Answers (2)

Julia Silge
Julia Silge

Reputation: 11663

Another way you can do this is with purrr, putting all the different formulae you want into in a dataframe list column:

library(MASS)
library(tidyverse)
library(broom)

formula_dfs <- tibble(formula_id = 1:2,
                      formula = c(as.formula(as.factor(cyl) ~ hp),
                                  as.formula(as.factor(cyl) ~ hp + qsec))) 

formula_dfs
#> # A tibble: 2 x 2
#>   formula_id formula  
#>        <int> <list>   
#> 1          1 <formula>
#> 2          2 <formula>

formula_dfs %>%
  mutate(polr_fit  = map(formula, polr, data = mtcars, Hess = TRUE),
         polr_coef = map(polr_fit, tidy, conf.int = TRUE, exponentiate = TRUE)) %>%
  unnest(polr_coef) %>%
  filter(coef.type=="coefficient")
#> Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
#> Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
#> # A tibble: 3 x 10
#>   formula_id formula   polr_fit term  estimate std.error statistic conf.low
#>        <int> <list>    <list>   <chr>    <dbl>     <dbl>     <dbl>    <dbl>
#> 1          1 <formula> <polr>   hp        1.12    0.0399     2.90     1.06 
#> 2          2 <formula> <polr>   hp        1.13    0.0452     2.72     1.06 
#> 3          2 <formula> <polr>   qsec      1.18    0.369      0.448    0.538
#> # … with 2 more variables: conf.high <dbl>, coef.type <chr>

Created on 2021-05-24 by the reprex package (v2.0.0)

Your regular model output is still there in the polr_fit column.

Upvotes: 1

Ronak Shah
Ronak Shah

Reputation: 389325

Something like this ?

library(broom)
library(tidyverse)

clean_model <- function(polr_results) {
  lapply(polr_results,  function(x) {
    broom::tidy(x, conf.int = TRUE, exponentiate = TRUE) %>%
      filter(coef.type=="coefficient")
  })
}

clean_model(list(mtcars_m1, mtcars_m2))

#[[1]]
# A tibble: 1 x 7
#  term  estimate std.error statistic conf.low conf.high coef.type  
#  <chr>    <dbl>     <dbl>     <dbl>    <dbl>     <dbl> <chr>      
#1 hp        1.12    0.0399      2.90     1.06      1.26 coefficient

#[[2]]
# A tibble: 2 x 7
#  term  estimate std.error statistic conf.low conf.high coef.type  
#  <chr>    <dbl>     <dbl>     <dbl>    <dbl>     <dbl> <chr>      
#1 hp        1.13    0.0452     2.72     1.06       1.29 coefficient
#2 qsec      1.18    0.369      0.448    0.538      2.51 coefficient

Upvotes: 1

Related Questions