All MMs and AMCEs

library(tidyverse)
library(targets)
library(tidybayes)
library(patchwork)
library(ggforce)
library(scales)
library(glue)
library(gt)
library(gtExtras)
library(here)

# Targets stuff
tar_config_set(store = here('_targets'),
               script = here('_targets.R'))

tar_load(c(grid_treatment_only, level_lookup, feature_lookup))
preds_all <- tar_read(preds_conditional_treatment_only)

invisible(list2env(tar_read(graphic_functions), .GlobalEnv))
invisible(list2env(tar_read(table_functions), .GlobalEnv))

theme_set(theme_ngo())

TODO: Hypotheses

• General basic ones, like only each lever, one at a time
• More specific ones that we preregistered. These conditional ones are tricky because there are so many different possible questions we could explore (at least 576 of them!). To avoid selective p-hacking, we pre-registered a handful that we’re most interested in, and we only explore those.

mms_org <- preds_all %>% 
  group_by(feat_org, .draw) %>% 
  summarize(avg = mean(.epred))

mms_issue <- preds_all %>% 
  group_by(feat_issue, .draw) %>% 
  summarize(avg = mean(.epred))

mms_transp <- preds_all %>% 
  group_by(feat_transp, .draw) %>% 
  summarize(avg = mean(.epred))

mms_acc <- preds_all %>% 
  group_by(feat_acc, .draw) %>% 
  summarize(avg = mean(.epred))

mms_funding <- preds_all %>% 
  group_by(feat_funding, .draw) %>% 
  summarize(avg = mean(.epred))

mms_govt <- preds_all %>% 
  group_by(feat_govt, .draw) %>% 
  summarize(avg = mean(.epred))

amces_org <- mms_org %>% 
  group_by(feat_org) %>% 
  compare_levels(variable = avg, by = feat_org, comparison = "control") %>% 
  ungroup() %>% 
    separate_wider_delim(
    feat_org,
    delim = " - ", 
    names = c("feature_level", "reference_level")
  ) %>% 
  add_row(avg = 0, feature_level = unique(.$reference_level))

amces_issue <- mms_issue %>% 
  group_by(feat_issue) %>% 
  compare_levels(variable = avg, by = feat_issue, comparison = "control") %>% 
  ungroup() %>% 
    separate_wider_delim(
    feat_issue,
    delim = " - ", 
    names = c("feature_level", "reference_level")
  ) %>% 
  add_row(avg = 0, feature_level = unique(.$reference_level))

amces_transp <- mms_transp %>% 
  group_by(feat_transp) %>% 
  compare_levels(variable = avg, by = feat_transp, comparison = "control") %>% 
  ungroup() %>% 
    separate_wider_delim(
    feat_transp,
    delim = " - ", 
    names = c("feature_level", "reference_level")
  ) %>% 
  add_row(avg = 0, feature_level = unique(.$reference_level))

amces_acc <- mms_acc %>% 
  group_by(feat_acc) %>% 
  compare_levels(variable = avg, by = feat_acc, comparison = "control") %>% 
  ungroup() %>% 
    separate_wider_delim(
    feat_acc,
    delim = " - ", 
    names = c("feature_level", "reference_level")
  ) %>% 
  add_row(avg = 0, feature_level = unique(.$reference_level))

amces_funding <- mms_funding %>% 
  group_by(feat_funding) %>% 
  compare_levels(variable = avg, by = feat_funding, comparison = "control") %>% 
  ungroup() %>% 
    separate_wider_delim(
    feat_funding,
    delim = " - ", 
    names = c("feature_level", "reference_level")
  ) %>% 
  add_row(avg = 0, feature_level = unique(.$reference_level))

amces_govt <- mms_govt %>% 
  group_by(feat_govt) %>% 
  compare_levels(variable = avg, by = feat_govt, comparison = "control") %>% 
  ungroup() %>% 
    separate_wider_delim(
    feat_govt,
    delim = " - ", 
    names = c("feature_level", "reference_level")
  ) %>% 
  add_row(avg = 0, feature_level = unique(.$reference_level))

amces_all <- bind_rows(
  lst(amces_org, amces_issue, amces_transp, amces_acc, amces_funding, amces_govt),
  .id = "amce_var"
) %>%
  left_join(select(feature_lookup, amce_var, feature_nice), by = join_by(amce_var)) %>%
  left_join(
    select(level_lookup, contains("_level")), 
    by = join_by(feature_level)
  )

mms_all <- bind_rows(
  lst(
    mms_org = rename(mms_org, feature_level = feat_org), 
    mms_issue = rename(mms_issue, feature_level = feat_issue), 
    mms_transp = rename(mms_transp, feature_level = feat_transp), 
    mms_acc = rename(mms_acc, feature_level = feat_acc), 
    mms_funding = rename(mms_funding, feature_level = feat_funding), 
    mms_govt = rename(mms_govt, feature_level = feat_govt)
  ),
  .id = "mm_var"
) %>%
  left_join(select(feature_lookup, mm_var, feature_nice), by = join_by(mm_var)) %>%
  left_join(
    select(level_lookup, contains("_level")), 
    by = join_by(feature_level)
  )

ggplot(
  mms_all,
  aes(x = avg, y = feature_short_level, fill = feature_nice)
) +
  geom_vline(xintercept = 0.25, color = clrs$prism[8], linetype = "dashed", linewidth = 0.25) +
  stat_halfeye(normalize = "groups") +
  guides(fill = "none") +
  facet_col(facets = "feature_nice", scales = "free_y", space = "free") +
  scale_x_continuous(labels = label_percent()) +
  scale_fill_manual(values = clrs$prism[1:6]) +
  labs(
    x = "Marginal means of probabilities",
    y = NULL,
    title = "Posterior marginal means"
  )

ggplot(
  amces_all,
  aes(x = avg, y = feature_short_level, fill = feature_nice)
) +
  geom_vline(xintercept = 0, color = clrs$prism[8], linetype = "dashed", linewidth = 0.25) +
  stat_halfeye(normalize = "groups") +  # Make the heights of the distributions equal within each facet
  guides(fill = "none") +
  facet_col(facets = "feature_nice", scales = "free_y", space = "free") +
  scale_x_continuous(labels = label_pp) +
  scale_fill_manual(values = clrs$prism[1:6]) +
  labs(
    x = "Percentage point change in probability of choice selection",
    y = NULL,
    title = "Posterior AMCEs"
  )