Skip to content

Instantly share code, notes, and snippets.

@jebyrnes
Last active February 10, 2023 20:29
Show Gist options
  • Save jebyrnes/d1bdea4ad4c8736c4b9e7ac290e8c940 to your computer and use it in GitHub Desktop.
Save jebyrnes/d1bdea4ad4c8736c4b9e7ac290e8c940 to your computer and use it in GitHub Desktop.
Posthoc contrasts with emmeans, tidybayes, and brms
library(tidyverse)
library(emmeans)
library(brms)
library(tidybayes)
warp.brms <- brm(breaks ~ wool * tension, data = warpbreaks)
#get the adjusted means
warp_em <- emmeans (warp.brms, ~ tension | wool)
warp_em
#get all possible contrasts
cont <- contrast(warm_em, "tukey")
cont
#get the posterior draws from the contrasts
cont_posterior <- gather_emmeans_draws(cont)
#plot
ggplot(cont_posterior,
aes(y = contrast, x = .value)) +
geom_halfeyeh() +
facet_wrap(~wool) +
geom_vline(xintercept = 0, color = "red", lty = 2)
ggplot(cont_posterior,
aes(y = contrast, x = .value, fill = wool, group = wool)) +
geom_halfeyeh(alpha = 0.5)+
geom_vline(xintercept = 0, color = "red", lty = 2)
#need to figure out how to make this look all posterior-y
emmip(warp.brms, ~ tension | wool, CIs = TRUE, cov.reduce = range)
@mjskay
Copy link

mjskay commented Nov 3, 2018

This is great!

Here's some quick brainstorming inspired by your the comment on your call to emmip...

warp.brms %>%
  emmeans( ~ tension | wool) %>%
  gather_emmeans_draws() %>%
  ggplot(aes(x = tension, y = .value)) +
  stat_lineribbon() +
  scale_fill_brewer(palette = "Greys") +
  facet_grid(~ wool) +
  theme_light()

image

warp.brms %>%
  emmeans( ~ tension | wool) %>%
  gather_emmeans_draws() %>%
  ggplot(aes(x = tension, y = .value)) +
  geom_eye() +
  stat_summary(aes(group = NA), fun.y = mean, geom = "line") +
  facet_grid(~ wool) +
  theme_light()

image

(reminds me I need to implement the vertical geom_halfeye sometime... or just a generalization of all of the interval+density variants)

warp.brms %>%
  emmeans( ~ tension | wool) %>%
  gather_emmeans_draws() %>%
  ggplot(aes(x = tension, y = .value, group = .draw)) +
  geom_line(alpha = .01) +
  scale_fill_brewer() +
  facet_grid(~ wool) +
  theme_light()

image

warp.brms %>%
  emmeans( ~ tension | wool) %>%
  gather_emmeans_draws() %>%
  ggplot(aes(x = tension, y = .value, fill = wool, color = wool)) +
  stat_lineribbon(alpha = 1/4) +
  theme_light()

image

Too many intervals here probably, it's a little mind-bending...

warp.brms %>%
  emmeans( ~ tension | wool) %>%
  gather_emmeans_draws() %>%
  ggplot(aes(x = tension, y = .value, fill = wool, color = wool)) +
  stat_lineribbon(alpha = 1/3, .width = .95) +
  theme_light()

image

Or you could go the other way and approximate a gradient plot (mind you this brings all the perceptual issues with gradients...)

warp.brms %>%
  emmeans( ~ tension | wool) %>%
  gather_emmeans_draws() %>%
  ggplot(aes(x = tension, y = .value, fill = wool, color = wool)) +
  stat_lineribbon(alpha = 1/41, .width = ppoints(40)) +
  theme_light()

image

Heh, the way stat_lineribbon interleaves layers makes that last one kind of painterly. Neat, but I wouldn't use it in practice :)

warp.brms %>%
  emmeans( ~ tension | wool) %>%
  gather_emmeans_draws() %>%
  ggplot(aes(x = tension, y = .value)) +
  stat_lineribbon(fill = "gray25", color = "gray25", size = .5, alpha = 1/41, .width = ppoints(40)) +
  stat_pointinterval() +
  facet_grid(~ wool) +
  theme_light()

image

warp.brms %>%
  emmeans( ~ tension | wool) %>%
  gather_emmeans_draws() %>%
  ggplot(aes(x = tension, y = .value)) +
  geom_line(aes(group = .draw), alpha = .005) +
  stat_pointinterval() +
  facet_grid(~ wool) +
  theme_light()

image

@JAQuent
Copy link

JAQuent commented Apr 24, 2022

Hey thank you so much for this code! This definitely helped me to get going myself. I am wondering though is there any reason why you used contrast(..., "tukey") instead of pairwise? I am getting the same results for both for my model but I want to be sure that I am not missing anything about using tukey for Bayesian mdoels.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment