vincentarelbundock / modelsummary

Beautiful and customizable model summaries in R.
http://modelsummary.com
Other
912 stars 75 forks source link

support for `emmGrid` #402

Closed strengejacke closed 2 years ago

strengejacke commented 2 years ago

Problem:

@mattansb @vincentarelbundock I just saw the modelsummary cannot cope with emGrid objects, also because model_parameters() returns a context-specific column name for the parameters (thus, it's not changed into term):

m <- glm(am ~ factor(cyl), family = binomial(), data = mtcars)

EList <- emmeans::emmeans(m, pairwise ~ cyl, type = "resp")
E <- emmeans::emmeans(m, ~ cyl, type = "resp")
C <- emmeans::contrast(E, method = "pairwise")

modelsummary::modelsummary(E)
#> Warning in get_gof(models[[j]], vcov_type[[i]], ...): `modelsummary could not extract goodness-of-fit statistics from a model
#> of class "emmGrid". The package tried a sequence of 2 helper functions:
#> 
#> broom::glance(model)
#> performance::model_performance(model)
#> 
#> One of these functions must return a one-row `data.frame`. The `modelsummary` website explains how to summarize unsupported models or add support for new models yourself:
#> 
#> https://vincentarelbundock.github.io/modelsummary/articles/modelsummary.html
#> Error in get_estimates(models[[j]], conf_level = conf_level, vcov = vcov[[i]], : `modelsummary could not extract the required information from a model
#> of class "emmGrid". The package tried a sequence of 2 helper functions to extract
#> estimates:
#> 
#> broom::tidy(model)
#> parameters::parameters(model)
#> 
#> To draw a table, one of these commands must return a `data.frame` with a
#> column named "term". The `modelsummary` website explains how to summarize
#> unsupported models or add support for new models yourself:
#> 
#> https://vincentarelbundock.github.io/modelsummary/articles/modelsummary.html
#> 
#> These errors messages were generated during extraction:
#> `broom::tidy(model)` did not return a data.frame with a `term` column.
#> `parameters::parameters(model)` did not return a data.frame with a `term` column.

parameters::model_parameters(E) |> insight::standardize_names(style = "broom")
#>   cyl  estimate  std.error conf.level   conf.low conf.high  statistic df.error
#> 1   4 0.7272727 0.13428152       0.95 0.41433537 0.9095152  1.4487819      Inf
#> 2   6 0.4285714 0.18704391       0.95 0.14373245 0.7701689 -0.3766642      Inf
#> 3   8 0.1428571 0.09352194       0.95 0.03596066 0.4268261 -2.3459642      Inf
#>      p.value
#> 1 0.14739849
#> 2 0.70642313
#> 3 0.01897793

Created on 2021-11-03 by the reprex package (v2.0.1)

What would you suggest? Fixing insight::standardize_names() here to force the column name to be "term"? However, model_parameters() and broom::tidy() show the same behaviour here, so this is probably something to be solved in modelsummary.

broom::tidy(E)
#> # A tibble: 3 x 7
#>     cyl  prob std.error    df  null z.ratio p.value
#>   <dbl> <dbl>     <dbl> <dbl> <dbl>   <dbl>   <dbl>
#> 1     4 0.727    0.134    Inf   0.5   1.45   0.147 
#> 2     6 0.429    0.187    Inf   0.5  -0.377  0.706 
#> 3     8 0.143    0.0935   Inf   0.5  -2.35   0.0190
bwiernik commented 2 years ago

Can we store the variable with the parameter labels as an attribute?

strengejacke commented 2 years ago

you mean something like attr(x, "parameter_column") <- "cyl"?

bwiernik commented 2 years ago

Yeah. That and one for the column with the coefficient/stat (Coefficient, Median, etc) label would actually also be helpful for see functions too

strengejacke commented 2 years ago

Just saw, we already have this as attribute, "parameter_names":

m <- glm(am ~ factor(cyl), family = binomial(), data = mtcars)

E <- emmeans::emmeans(m, ~ cyl, type = "resp")
parameters::model_parameters(E) |> str()
#> Classes 'parameters_model', 'see_parameters_model' and 'data.frame': 3 obs. of  9 variables:
#>  $ cyl        : num  4 6 8
#>  $ Coefficient: num  0.727 0.429 0.143
#>  $ SE         : num  0.1343 0.187 0.0935
#>  $ CI         : num  0.95 0.95 0.95
#>  $ CI_low     : num  0.414 0.144 0.036
#>  $ CI_high    : num  0.91 0.77 0.427
#>  $ z          : num  1.449 -0.377 -2.346
#>  $ df_error   : num  Inf Inf Inf
#>  $ p          : num  0.147 0.706 0.019
#>  - attr(*, "ci")= num 0.95
#>  - attr(*, "test_statistic")= chr "z-statistic"
#>  - attr(*, "verbose")= logi TRUE
#>  - attr(*, "exponentiate")= logi FALSE
#>  - attr(*, "ordinal_model")= logi FALSE
#>  - attr(*, "linear_model")= logi FALSE
#>  - attr(*, "mixed_model")= logi FALSE
#>  - attr(*, "model_class")= chr "emmGrid"
#>   ..- attr(*, "package")= chr "emmeans"
#>  - attr(*, "bootstrap")= logi FALSE
#>  - attr(*, "iterations")= num 1000
#>  - attr(*, "p_adjust")= chr "none"
#>  - attr(*, "robust_vcov")= logi FALSE
#>  - attr(*, "ignore_group")= logi TRUE
#>  - attr(*, "ran_pars")= logi TRUE
#>  - attr(*, "show_summary")= logi FALSE
#>  - attr(*, "title")= chr ""
#>  - attr(*, "model_formula")= chr "NULL"
#>  - attr(*, "coefficient_name")= chr "Probability"
#>  - attr(*, "zi_coefficient_name")= chr "Log-Odds"
#>  - attr(*, "digits")= num 2
#>  - attr(*, "ci_digits")= num 2
#>  - attr(*, "p_digits")= num 3
#>  - attr(*, "footer_digits")= num 3
#>  - attr(*, "object_name")= chr "E"
#>  - attr(*, "parameter_names")= chr "cyl"

Created on 2021-11-03 by the reprex package (v2.0.1)

vincentarelbundock commented 2 years ago

Thanks for the report!

My initial reaction to this is that the current behavior of both broom and parameters feels inconsistent with their usual practice. For standardization, we need column names to identify the type of information contained in the column, and not the information itself; the data is always presented in "long" format. In this particular case, the column name is itself a piece of information, which is inconsistent with how the two packages treat other models.

I think the behavior should be similar to nnet::multinom, where we get a term column and a y.level (in broom) or response (in parameters) column.

> m <- nnet::multinom(factor(cyl) ~ mpg, data = mtcars)

> broom::tidy(m)

# A tibble: 4 × 6
  y.level term        estimate std.error statistic p.value
  <chr>   <chr>          <dbl>     <dbl>     <dbl>   <dbl>
1 6       (Intercept)    47.3      35.0       1.35  0.177 
2 6       mpg            -2.21      1.64     -1.35  0.178 
3 8       (Intercept)    72.4      37.2       1.95  0.0513
4 8       mpg            -3.58      1.77     -2.02  0.0437

> parameters::model_parameters(m) |> insight::standardize_names(style = "broom")
         term  estimate std.error conf.level    conf.low   conf.high statistic
1 (Intercept) 47.252432 34.975171       0.95 -21.2976435 115.8025065  1.351028
2         mpg -2.205418  1.637963       0.95  -5.4157653   1.0049299 -1.346440
3 (Intercept) 72.440246 37.175162       0.95  -0.4217332 145.3022247  1.948619
4         mpg -3.579991  1.774693       0.95  -7.0583242  -0.1016573 -2.017246
  df.error    p.value response
1       28 0.17668650        6
2       28 0.17816078        6
3       28 0.05134088        8
4       28 0.04366989        8
vincentarelbundock commented 2 years ago

Addendum: The problem with fixing this downstream is that modelsummary treats all columns of the parameters output as statistics. If the only information I have is the output shown above, there is no way to know if cyl is actually a term name or just another statistic that modelsummary should display in the table, like p.value or std.error.

vincentarelbundock commented 2 years ago

I opened an issue on the broom repo: https://github.com/tidymodels/broom/issues/1063

bwiernik commented 2 years ago

Addendum: The problem with fixing this downstream is that modelsummary treats all columns of the parameters output as statistics. If the only information I have is the output shown above, there is no way to know if cyl is actually a term name or just another statistic that modelsummary should display in the table, like p.value or std.error.

It's in the attributes as "parameter name". The issue with emGrid here is that there may be multiple columns if the user computes marginal means across multiple variables

vincentarelbundock commented 2 years ago

The issue with emGrid here is that there may be multiple columns if the user computes marginal means across multiple variables

Thanks. I think I get it now. One other issue that this raises is that modelsummary only supports 2 "levels" of parameters: term and group (the latter could be y.level, response, or whatever). This implies that, at the moment, there isn't a great way to print out marginal means on a 3d (or higher dimensional) grid.

For a single variable, it's pretty easy to write a pre-processing function, but I don't yet have a clear view of a solution for more variables:

library(modelsummary)
library(tidyr)
m <- glm(am ~ factor(cyl), family = binomial(), data = mtcars)
x <- emmeans::emmeans(m, ~ cyl, type = "resp")

preprocess <- function(x) {
    ti <- x |>
        parameters::model_parameters() |>
        insight::standardize_names(style = "broom")
    ti <- pivot_longer(ti, cols = attr(ti, "parameter_names"), 
                       names_to = "term", values_to = "group")
    ti$group <- sprintf("%.0f", ti$group)
    out <- list("tidy" = ti, "glance" = NULL) 
    class(out) <- c("modelsummary_list", "list")
    return(out)
}

modelsummary(preprocess(x), group = term + group ~ model, output = "markdown")

|    |   | Model 1 |
|:---|:--|:-------:|
|cyl |4  |  0.727  |
|    |   | (0.134) |
|    |6  |  0.429  |
|    |   | (0.187) |
|    |8  |  0.143  |
|    |   | (0.094) |
bwiernik commented 2 years ago

The logic in emmeans is that each row is atomic--it's the mean, etc. for that specific combination of variable levels. You could accomplish something similar by concatenating the variables into somthing like cyl = 4, am = 0, gear = 3 with term being something like cyl, am, gear.

The reason behind the separate columns for variables in emmeans, besides it being potentially natural for users, is that it enables filtering/sorting by variables separately and facilitates computing of various contrasts that condition on or marginalize across various variables. If modelsummary stored this information about the individual variables somewhere, perhaps as a list column or list attribute, such things could still be possible.

m <- lm(mpg ~ 0 + factor(cyl) + factor(gear) + factor(am), 
        data = mtcars)

emmeans::emmeans(m, ~ cyl + gear + am) |> 
  emmeans::contrast(simple = c("gear", "cyl"))
#> am = 0:
#>  contrast   estimate   SE df t.ratio p.value
#>  4 3 effect    6.459 1.97 26  3.271  0.0091 
#>  6 3 effect    0.408 1.96 26  0.209  0.8364 
#>  8 3 effect   -3.747 1.40 26 -2.681  0.0283 
#>  4 4 effect    5.456 1.09 26  4.986  0.0003 
#>  6 4 effect   -0.594 1.24 26 -0.481  0.7137 
#>  8 4 effect   -4.750 1.97 26 -2.412  0.0348 
#>  4 5 effect    4.341 1.74 26  2.498  0.0344 
#>  6 5 effect   -1.710 1.57 26 -1.086  0.3694 
#>  8 5 effect   -5.865 1.74 26 -3.380  0.0091 
#> 
#> am = 1:
#>  contrast   estimate   SE df t.ratio p.value
#>  4 3 effect    6.459 1.97 26  3.271  0.0091 
#>  6 3 effect    0.408 1.96 26  0.209  0.8364 
#>  8 3 effect   -3.747 1.40 26 -2.681  0.0283 
#>  4 4 effect    5.456 1.09 26  4.986  0.0003 
#>  6 4 effect   -0.594 1.24 26 -0.481  0.7137 
#>  8 4 effect   -4.750 1.97 26 -2.412  0.0348 
#>  4 5 effect    4.341 1.74 26  2.498  0.0344 
#>  6 5 effect   -1.710 1.57 26 -1.086  0.3694 
#>  8 5 effect   -5.865 1.74 26 -3.380  0.0091 
#> 
#> P value adjustment: fdr method for 9 tests

emmeans::emmeans(m, ~ cyl + gear + am) |> 
  emmeans::contrast(simple = list("gear", "cyl"))
#> $`simple contrasts for gear`
#> cyl = 4, am = 0:
#>  contrast estimate   SE df t.ratio p.value
#>  3 effect   1.0403 1.48 26  0.702  0.7331 
#>  4 effect   0.0374 1.08 26  0.035  0.9726 
#>  5 effect  -1.0777 1.34 26 -0.803  0.7331 
#> 
#> cyl = 6, am = 0:
#>  contrast estimate   SE df t.ratio p.value
#>  3 effect   1.0403 1.48 26  0.702  0.7331 
#>  4 effect   0.0374 1.08 26  0.035  0.9726 
#>  5 effect  -1.0777 1.34 26 -0.803  0.7331 
#> 
#> cyl = 8, am = 0:
#>  contrast estimate   SE df t.ratio p.value
#>  3 effect   1.0403 1.48 26  0.702  0.7331 
#>  4 effect   0.0374 1.08 26  0.035  0.9726 
#>  5 effect  -1.0777 1.34 26 -0.803  0.7331 
#> 
#> cyl = 4, am = 1:
#>  contrast estimate   SE df t.ratio p.value
#>  3 effect   1.0403 1.48 26  0.702  0.7331 
#>  4 effect   0.0374 1.08 26  0.035  0.9726 
#>  5 effect  -1.0777 1.34 26 -0.803  0.7331 
#> 
#> cyl = 6, am = 1:
#>  contrast estimate   SE df t.ratio p.value
#>  3 effect   1.0403 1.48 26  0.702  0.7331 
#>  4 effect   0.0374 1.08 26  0.035  0.9726 
#>  5 effect  -1.0777 1.34 26 -0.803  0.7331 
#> 
#> cyl = 8, am = 1:
#>  contrast estimate   SE df t.ratio p.value
#>  3 effect   1.0403 1.48 26  0.702  0.7331 
#>  4 effect   0.0374 1.08 26  0.035  0.9726 
#>  5 effect  -1.0777 1.34 26 -0.803  0.7331 
#> 
#> P value adjustment: fdr method for 3 tests 
#> 
#> $`simple contrasts for cyl`
#> gear = 3, am = 0:
#>  contrast estimate    SE df t.ratio p.value
#>  4 effect    5.419 0.993 26  5.458  <.0001 
#>  6 effect   -0.632 0.944 26 -0.669  0.5094 
#>  8 effect   -4.787 1.108 26 -4.321  0.0003 
#> 
#> gear = 4, am = 0:
#>  contrast estimate    SE df t.ratio p.value
#>  4 effect    5.419 0.993 26  5.458  <.0001 
#>  6 effect   -0.632 0.944 26 -0.669  0.5094 
#>  8 effect   -4.787 1.108 26 -4.321  0.0003 
#> 
#> gear = 5, am = 0:
#>  contrast estimate    SE df t.ratio p.value
#>  4 effect    5.419 0.993 26  5.458  <.0001 
#>  6 effect   -0.632 0.944 26 -0.669  0.5094 
#>  8 effect   -4.787 1.108 26 -4.321  0.0003 
#> 
#> gear = 3, am = 1:
#>  contrast estimate    SE df t.ratio p.value
#>  4 effect    5.419 0.993 26  5.458  <.0001 
#>  6 effect   -0.632 0.944 26 -0.669  0.5094 
#>  8 effect   -4.787 1.108 26 -4.321  0.0003 
#> 
#> gear = 4, am = 1:
#>  contrast estimate    SE df t.ratio p.value
#>  4 effect    5.419 0.993 26  5.458  <.0001 
#>  6 effect   -0.632 0.944 26 -0.669  0.5094 
#>  8 effect   -4.787 1.108 26 -4.321  0.0003 
#> 
#> gear = 5, am = 1:
#>  contrast estimate    SE df t.ratio p.value
#>  4 effect    5.419 0.993 26  5.458  <.0001 
#>  6 effect   -0.632 0.944 26 -0.669  0.5094 
#>  8 effect   -4.787 1.108 26 -4.321  0.0003 
#> 
#> P value adjustment: fdr method for 3 tests

Created on 2021-11-03 by the reprex package (v2.0.1)

strengejacke commented 2 years ago

The last example works in principle, but the components are not visually separated from each other:

m <- lm(mpg ~ 0 + factor(cyl) + factor(gear) + factor(am), 
        data = mtcars)

emmeans::emmeans(m, ~ cyl + gear + am) |> 
  emmeans::contrast(simple = list("gear", "cyl")) |> 
  parameters::model_parameters() |> 
  insight::standardize_names(style = "broom") |> 
  as.data.frame()
#>                                     term    estimate std.error conf.level
#> 1   contrast [3 effect], cyl [4], am [0]  1.04027384 1.4812030       0.95
#> 2   contrast [4 effect], cyl [4], am [0]  0.03744098 1.0793432       0.95
#> 3   contrast [5 effect], cyl [4], am [0] -1.07771483 1.3419318       0.95
#> 4   contrast [3 effect], cyl [6], am [0]  1.04027384 1.4812030       0.95
#> 5   contrast [4 effect], cyl [6], am [0]  0.03744098 1.0793432       0.95
#> 6   contrast [5 effect], cyl [6], am [0] -1.07771483 1.3419318       0.95
#> 7   contrast [3 effect], cyl [8], am [0]  1.04027384 1.4812030       0.95
#> 8   contrast [4 effect], cyl [8], am [0]  0.03744098 1.0793432       0.95
#> 9   contrast [5 effect], cyl [8], am [0] -1.07771483 1.3419318       0.95
#> 10  contrast [3 effect], cyl [4], am [1]  1.04027384 1.4812030       0.95
#> 11  contrast [4 effect], cyl [4], am [1]  0.03744098 1.0793432       0.95
#> 12  contrast [5 effect], cyl [4], am [1] -1.07771483 1.3419318       0.95
#> 13  contrast [3 effect], cyl [6], am [1]  1.04027384 1.4812030       0.95
#> 14  contrast [4 effect], cyl [6], am [1]  0.03744098 1.0793432       0.95
#> 15  contrast [5 effect], cyl [6], am [1] -1.07771483 1.3419318       0.95
#> 16  contrast [3 effect], cyl [8], am [1]  1.04027384 1.4812030       0.95
#> 17  contrast [4 effect], cyl [8], am [1]  0.03744098 1.0793432       0.95
#> 18  contrast [5 effect], cyl [8], am [1] -1.07771483 1.3419318       0.95
#> 19 contrast [4 effect], gear [3], am [0]  5.41888574 0.9927506       0.95
#> 20 contrast [6 effect], gear [3], am [0] -0.63182247 0.9444253       0.95
#> 21 contrast [8 effect], gear [3], am [0] -4.78706327 1.1078771       0.95
#> 22 contrast [4 effect], gear [4], am [0]  5.41888574 0.9927506       0.95
#> 23 contrast [6 effect], gear [4], am [0] -0.63182247 0.9444253       0.95
#> 24 contrast [8 effect], gear [4], am [0] -4.78706327 1.1078771       0.95
#> 25 contrast [4 effect], gear [5], am [0]  5.41888574 0.9927506       0.95
#> 26 contrast [6 effect], gear [5], am [0] -0.63182247 0.9444253       0.95
#> 27 contrast [8 effect], gear [5], am [0] -4.78706327 1.1078771       0.95
#> 28 contrast [4 effect], gear [3], am [1]  5.41888574 0.9927506       0.95
#> 29 contrast [6 effect], gear [3], am [1] -0.63182247 0.9444253       0.95
#> 30 contrast [8 effect], gear [3], am [1] -4.78706327 1.1078771       0.95
#> 31 contrast [4 effect], gear [4], am [1]  5.41888574 0.9927506       0.95
#> 32 contrast [6 effect], gear [4], am [1] -0.63182247 0.9444253       0.95
#> 33 contrast [8 effect], gear [4], am [1] -4.78706327 1.1078771       0.95
#> 34 contrast [4 effect], gear [5], am [1]  5.41888574 0.9927506       0.95
#> 35 contrast [6 effect], gear [5], am [1] -0.63182247 0.9444253       0.95
#> 36 contrast [8 effect], gear [5], am [1] -4.78706327 1.1078771       0.95
#>     conf.low conf.high   statistic df.error      p.value
#> 1  -2.004383  4.084930  0.70231685       26 7.330831e-01
#> 2  -2.181181  2.256063  0.03468867       26 9.725929e-01
#> 3  -3.836095  1.680666 -0.80310698       26 7.330831e-01
#> 4  -2.004383  4.084930  0.70231685       26 7.330831e-01
#> 5  -2.181181  2.256063  0.03468867       26 9.725929e-01
#> 6  -3.836095  1.680666 -0.80310698       26 7.330831e-01
#> 7  -2.004383  4.084930  0.70231685       26 7.330831e-01
#> 8  -2.181181  2.256063  0.03468867       26 9.725929e-01
#> 9  -3.836095  1.680666 -0.80310698       26 7.330831e-01
#> 10 -2.004383  4.084930  0.70231685       26 7.330831e-01
#> 11 -2.181181  2.256063  0.03468867       26 9.725929e-01
#> 12 -3.836095  1.680666 -0.80310698       26 7.330831e-01
#> 13 -2.004383  4.084930  0.70231685       26 7.330831e-01
#> 14 -2.181181  2.256063  0.03468867       26 9.725929e-01
#> 15 -3.836095  1.680666 -0.80310698       26 7.330831e-01
#> 16 -2.004383  4.084930  0.70231685       26 7.330831e-01
#> 17 -2.181181  2.256063  0.03468867       26 9.725929e-01
#> 18 -3.836095  1.680666 -0.80310698       26 7.330831e-01
#> 19  3.378258  7.459514  5.45845608       26 3.019209e-05
#> 20 -2.573117  1.309472 -0.66900204       26 5.093914e-01
#> 21 -7.064337 -2.509789 -4.32093367       26 3.024482e-04
#> 22  3.378258  7.459514  5.45845608       26 3.019209e-05
#> 23 -2.573117  1.309472 -0.66900204       26 5.093914e-01
#> 24 -7.064337 -2.509789 -4.32093367       26 3.024482e-04
#> 25  3.378258  7.459514  5.45845608       26 3.019209e-05
#> 26 -2.573117  1.309472 -0.66900204       26 5.093914e-01
#> 27 -7.064337 -2.509789 -4.32093367       26 3.024482e-04
#> 28  3.378258  7.459514  5.45845608       26 3.019209e-05
#> 29 -2.573117  1.309472 -0.66900204       26 5.093914e-01
#> 30 -7.064337 -2.509789 -4.32093367       26 3.024482e-04
#> 31  3.378258  7.459514  5.45845608       26 3.019209e-05
#> 32 -2.573117  1.309472 -0.66900204       26 5.093914e-01
#> 33 -7.064337 -2.509789 -4.32093367       26 3.024482e-04
#> 34  3.378258  7.459514  5.45845608       26 3.019209e-05
#> 35 -2.573117  1.309472 -0.66900204       26 5.093914e-01
#> 36 -7.064337 -2.509789 -4.32093367       26 3.024482e-04
#>                    component
#> 1  simple contrasts for gear
#> 2  simple contrasts for gear
#> 3  simple contrasts for gear
#> 4  simple contrasts for gear
#> 5  simple contrasts for gear
#> 6  simple contrasts for gear
#> 7  simple contrasts for gear
#> 8  simple contrasts for gear
#> 9  simple contrasts for gear
#> 10 simple contrasts for gear
#> 11 simple contrasts for gear
#> 12 simple contrasts for gear
#> 13 simple contrasts for gear
#> 14 simple contrasts for gear
#> 15 simple contrasts for gear
#> 16 simple contrasts for gear
#> 17 simple contrasts for gear
#> 18 simple contrasts for gear
#> 19  simple contrasts for cyl
#> 20  simple contrasts for cyl
#> 21  simple contrasts for cyl
#> 22  simple contrasts for cyl
#> 23  simple contrasts for cyl
#> 24  simple contrasts for cyl
#> 25  simple contrasts for cyl
#> 26  simple contrasts for cyl
#> 27  simple contrasts for cyl
#> 28  simple contrasts for cyl
#> 29  simple contrasts for cyl
#> 30  simple contrasts for cyl
#> 31  simple contrasts for cyl
#> 32  simple contrasts for cyl
#> 33  simple contrasts for cyl
#> 34  simple contrasts for cyl
#> 35  simple contrasts for cyl
#> 36  simple contrasts for cyl

emmeans::emmeans(m, ~ cyl + gear + am) |> 
  emmeans::contrast(simple = list("gear", "cyl")) |> 
  modelsummary::modelsummary(output = "markdown")
#> Warning in get_gof(models[[j]], vcov_type[[i]], ...): `modelsummary could not extract goodness-of-fit statistics from a model
#> of class "emm_list". The package tried a sequence of 2 helper functions:
#> 
#> broom::glance(model)
#> performance::model_performance(model)
#> 
#> One of these functions must return a one-row `data.frame`. The `modelsummary` website explains how to summarize unsupported models or add support for new models yourself:
#> 
#> https://vincentarelbundock.github.io/modelsummary/articles/modelsummary.html
Model 1
contrast [3 effect], cyl [4], am [0] 1.040
(1.481)
contrast [4 effect], cyl [4], am [0] 0.037
(1.079)
contrast [5 effect], cyl [4], am [0] -1.078
(1.342)
contrast [3 effect], cyl [6], am [0] 1.040
(1.481)
contrast [4 effect], cyl [6], am [0] 0.037
(1.079)
contrast [5 effect], cyl [6], am [0] -1.078
(1.342)
contrast [3 effect], cyl [8], am [0] 1.040
(1.481)
contrast [4 effect], cyl [8], am [0] 0.037
(1.079)
contrast [5 effect], cyl [8], am [0] -1.078
(1.342)
contrast [3 effect], cyl [4], am [1] 1.040
(1.481)
contrast [4 effect], cyl [4], am [1] 0.037
(1.079)
contrast [5 effect], cyl [4], am [1] -1.078
(1.342)
contrast [3 effect], cyl [6], am [1] 1.040
(1.481)
contrast [4 effect], cyl [6], am [1] 0.037
(1.079)
contrast [5 effect], cyl [6], am [1] -1.078
(1.342)
contrast [3 effect], cyl [8], am [1] 1.040
(1.481)
contrast [4 effect], cyl [8], am [1] 0.037
(1.079)
contrast [5 effect], cyl [8], am [1] -1.078
(1.342)
contrast [4 effect], gear [3], am [0] 5.419
(0.993)
contrast [6 effect], gear [3], am [0] -0.632
(0.944)
contrast [8 effect], gear [3], am [0] -4.787
(1.108)
contrast [4 effect], gear [4], am [0] 5.419
(0.993)
contrast [6 effect], gear [4], am [0] -0.632
(0.944)
contrast [8 effect], gear [4], am [0] -4.787
(1.108)
contrast [4 effect], gear [5], am [0] 5.419
(0.993)
contrast [6 effect], gear [5], am [0] -0.632
(0.944)
contrast [8 effect], gear [5], am [0] -4.787
(1.108)
contrast [4 effect], gear [3], am [1] 5.419
(0.993)
contrast [6 effect], gear [3], am [1] -0.632
(0.944)
contrast [8 effect], gear [3], am [1] -4.787
(1.108)
contrast [4 effect], gear [4], am [1] 5.419
(0.993)
contrast [6 effect], gear [4], am [1] -0.632
(0.944)
contrast [8 effect], gear [4], am [1] -4.787
(1.108)
contrast [4 effect], gear [5], am [1] 5.419
(0.993)
contrast [6 effect], gear [5], am [1] -0.632
(0.944)
contrast [8 effect], gear [5], am [1] -4.787
(1.108)

Created on 2021-11-03 by the reprex package (v2.0.1)

vincentarelbundock commented 2 years ago

Thanks both!

I like @bwiernik 's idea of concatenating into a term column as a stopgap measure. But @strengejacke 's example shows the limits of this, as things become unreadable quickly.

In an ideal world, we would allow arbitrary column names in the group argument formula. This would allow users to do this:

modelsummary(mod, group = gear + cyl + am ~ model)

and get a neat table with visual separate for each combination of gear/cyl/am.

Unfortunately, this would require a major refactor of the group code and -- most dauntingly -- would require me to dive deep in the Base R reshape manual. I'm not sure I have the stamina for this right now :)

That said, I do hope to do something about this in the medium run, so I'll leave this issue open. Thanks again for the report and suggestions!

bwiernik commented 2 years ago

You could use datawizard::reshape_longer() instead, which has syntax more like tidyr::pivot_longer()

Another option would be to convert the the group column into a dataframe column.

vincentarelbundock commented 2 years ago

Thought about this a lot and concluded that there is not a clean "generic" way to handle these objects. The number of dimensions is not pre-determined, so I would have to have a fair amount of emmlist-specific code to handle this reasonably.

Thanks both for your input!