Extract draws from the result of a call to emmeans::emmeans() (formerly lsmeans) or emmeans::ref_grid() applied to a Bayesian model.

gather_emmeans_draws(object, value = ".value", ...)

# S3 method for default
gather_emmeans_draws(object, value = ".value", ...)

# S3 method for emm_list
gather_emmeans_draws(object, value = ".value", grid = ".grid", ...)

Arguments

object

An emmGrid object such as returned by emmeans::ref_grid() or emmeans::emmeans().

value

The name of the output column to use to contain the values of draws. Defaults to ".value".

...

Additional arguments passed to the underlying method for the type of object given.

grid

If object is an emmeans::emm_list(), the name of the output column to use to contain the name of the reference grid that a given row corresponds to. Defaults to ".grid".

Value

A tidy data frame of draws. The columns of the reference grid are returned as-is, with an additional column called .value (by default) containing marginal draws. The resulting data frame is grouped by the columns from the reference grid to make use of summary functions like point_interval() straightforward. If object is an emmeans::emm_list(), which contains estimates from different reference grids, an additional column with the default name of ".grid" is added to indicate the reference grid for each row in the output. The name of this column is controlled by the grid argument.

Details

emmeans::emmeans() provides a convenient syntax for generating draws from "estimated marginal means" from a model, and can be applied to various Bayesian models, like rstanarm::stanreg-objects and MCMCglmm::MCMCglmm(). Given a emmeans::ref_grid() object as returned by functions like emmeans::ref_grid() or emmeans::emmeans() applied to a Bayesian model, gather_emmeans_draws returns a tidy format data frame of draws from the marginal posterior distributions generated by emmeans::emmeans().

Author

Matthew Kay

Examples

# \donttest{

library(dplyr)
library(magrittr)
#> 
#> Attaching package: ‘magrittr’
#> The following objects are masked from ‘package:testthat’:
#> 
#>     equals, is_less_than, not
library(brms)
library(emmeans)

# Here's an example dataset with a categorical predictor (`condition`) with several levels:
set.seed(5)
n = 10
n_condition = 5
ABC = tibble(
  condition = rep(c("A","B","C","D","E"), n),
  response = rnorm(n * 5, c(0,1,2,1,-1), 0.5)
)

m = brm(response ~ condition, data = ABC,
  # 1 chain / few iterations just so example runs quickly
  # do not use in practice
  chains = 1, iter = 500)
#> Compiling Stan program...
#> Start sampling
#> 
#> SAMPLING FOR MODEL '42694db623b401afc3f1af679e2e9a25' NOW (CHAIN 1).
#> Chain 1: 
#> Chain 1: Gradient evaluation took 0 seconds
#> Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
#> Chain 1: Adjust your expectations accordingly!
#> Chain 1: 
#> Chain 1: 
#> Chain 1: Iteration:   1 / 500 [  0%]  (Warmup)
#> Chain 1: Iteration:  50 / 500 [ 10%]  (Warmup)
#> Chain 1: Iteration: 100 / 500 [ 20%]  (Warmup)
#> Chain 1: Iteration: 150 / 500 [ 30%]  (Warmup)
#> Chain 1: Iteration: 200 / 500 [ 40%]  (Warmup)
#> Chain 1: Iteration: 250 / 500 [ 50%]  (Warmup)
#> Chain 1: Iteration: 251 / 500 [ 50%]  (Sampling)
#> Chain 1: Iteration: 300 / 500 [ 60%]  (Sampling)
#> Chain 1: Iteration: 350 / 500 [ 70%]  (Sampling)
#> Chain 1: Iteration: 400 / 500 [ 80%]  (Sampling)
#> Chain 1: Iteration: 450 / 500 [ 90%]  (Sampling)
#> Chain 1: Iteration: 500 / 500 [100%]  (Sampling)
#> Chain 1: 
#> Chain 1:  Elapsed Time: 0.009 seconds (Warm-up)
#> Chain 1:                0.007 seconds (Sampling)
#> Chain 1:                0.016 seconds (Total)
#> Chain 1: 

# Once we've fit the model, we can use emmeans() (and functions
# from that package) to get whatever marginal distributions we want.
# For example, we can get marginal means by condition:
m %>%
  emmeans(~ condition) %>%
  gather_emmeans_draws() %>%
  median_qi()
#> # A tibble: 5 x 7
#>   condition .value .lower .upper .width .point .interval
#>   <fct>      <dbl>  <dbl>  <dbl>  <dbl> <chr>  <chr>    
#> 1 A          0.199 -0.170  0.567   0.95 median qi       
#> 2 B          1.00   0.653  1.33    0.95 median qi       
#> 3 C          1.87   1.50   2.23    0.95 median qi       
#> 4 D          1.04   0.606  1.39    0.95 median qi       
#> 5 E         -0.936 -1.28  -0.543   0.95 median qi       

# or we could get pairwise differences:
m %>%
  emmeans(~ condition) %>%
  contrast(method = "pairwise") %>%
  gather_emmeans_draws() %>%
  median_qi()
#> # A tibble: 10 x 7
#>    contrast  .value .lower .upper .width .point .interval
#>    <chr>      <dbl>  <dbl>  <dbl>  <dbl> <chr>  <chr>    
#>  1 A - B    -0.815  -1.28  -0.255   0.95 median qi       
#>  2 A - C    -1.68   -2.17  -1.24    0.95 median qi       
#>  3 A - D    -0.826  -1.38  -0.347   0.95 median qi       
#>  4 A - E     1.14    0.644  1.54    0.95 median qi       
#>  5 B - C    -0.856  -1.33  -0.421   0.95 median qi       
#>  6 B - D    -0.0294 -0.550  0.479   0.95 median qi       
#>  7 B - E     1.95    1.44   2.40    0.95 median qi       
#>  8 C - D     0.819   0.373  1.37    0.95 median qi       
#>  9 C - E     2.81    2.27   3.28    0.95 median qi       
#> 10 D - E     1.95    1.43   2.51    0.95 median qi       

# see the documentation of emmeans() for more examples of types of
# contrasts supported by that packge.

# }