Hey all–

On our call two (three?!) weeks ago:

So what I want to do here is just clarify what I did, especially since I was just following the R models that Katie wrote and saved in dropbox.

First, we’re going to load all the cleaned data

# Initial settings --------------------------------------------------------

library(tidyverse)
library(estimatr)
library(metafor)
library(emmeans)
library(magrittr)
library(here)
library(broom)

# Read data ---------------------------------------------------------------
t_data <- tibble(
  nms = str_c(
    here(),
    "/data/cleaned/"
    ) %>%
    dir() %>% 
    str_subset("study ") %>% 
    str_extract(
      "study \\d"
    ),
  tabs = str_c(
    here(),
    "/data/cleaned/"
    ) %>%
    dir(full.names = T) %>% 
    str_subset("study ") %>% 
    map(readRDS)
  )

The object t_data should look like this:

## # A tibble: 5 × 2
##   nms     tabs                  
##   <chr>   <list>                
## 1 study 1 <tibble [1,827 × 220]>
## 2 study 2 <tibble [1,739 × 46]> 
## 3 study 3 <tibble [3,207 × 52]> 
## 4 study 4 <tibble [2,789 × 63]> 
## 5 study 5 <tibble [1,983 × 69]>

Now we estimate all the models I could find, including many which we’re not going to plot

# Study 3 -----------------------------------------------------------------
t_s3 <- tribble(
  ~mod_name, ~frm, ~tab,
  "perc_acc_ovr",      "out2 ~ group + party + male + age + white + college",  t_data$tabs[[3]],
  "perc_acc_ov_nct",   "out2 ~ group ",                                        t_data$tabs[[3]],
  "perc_acc_reps",     "out2 ~ group  + male + age + white + college",         t_data$tabs[[3]] %>%
                                                                                  filter(party == "rep"),
  "perc_acc_reps_nct", "out2 ~ group",                                         t_data$tabs[[3]] %>% 
                                                                                  filter(party == "rep"),
  "perc_acc_dems",     "out2 ~ group  + male + age + white + college",         t_data$tabs[[3]] %>% 
                                                                                  filter(party == "dem"),
  "perc_acc_dems_nct", "out2 ~ group",                                         t_data$tabs[[3]] %>% 
                                                                                  filter(party == "dem"),
  "change_ovr",        "change ~ group + male + age + white + college",        t_data$tabs[[3]],
  "change_reps",       "change ~ group + male + age + white + college",        t_data$tabs[[3]] %>% 
                                                                                  filter(party == "rep"),
  "change_dems",       "change ~ group + male + age + white + college",        t_data$tabs[[3]] %>% 
                                                                                  filter(party == "dem"),
  # certainty analysis, exploratory
  "cert_ovr",          "cert2 ~ group + party + male + age + white + college", t_data$tabs[[3]],
  "cert_ovr_nct",      "cert2 ~ group",                                        t_data$tabs[[3]],
  "cert_reps",         "cert2 ~ group + male + age + white + college",         t_data$tabs[[3]] %>% 
                                                                                  filter(party == "rep"),
  "cert_reps_nct",     "cert2 ~ group",                                        t_data$tabs[[3]] %>% 
                                                                                  filter(party == "rep"),
  "cert_dems",         "cert2 ~ group + male + age + white + college",         t_data$tabs[[3]] %>% 
                                                                                  filter(party == "dem"),
  "cert_dems_nct",     "cert2 ~ group",                                        t_data$tabs[[3]] %>% 
                                                                                  filter(party == "dem")
  ) %>% 
  mutate(
    mods = frm %>% 
      map2(
        tab, 
        \(i, j)
        
        i %>% 
          as.formula %>% 
          lm_robust(
            data = j
          )
      )
  )

# Study 4 -----------------------------------------------------------------

t_s4  <- tribble(
  ~mod_name, ~frm, ~tab,
  "perc_acc_ovr",      "out2 ~ group + party + male + age + white + college",  t_data$tabs[[4]],
  "perc_acc_ovr_nct",  "out2 ~ group",                                         t_data$tabs[[4]],
  "perc_acc_reps",     "out2 ~ group + male + age + white + college",          t_data$tabs[[4]] %>% 
                                                                                  filter(party == "rep"),
  "perc_acc_reps_nct", "out2 ~ group",                                         t_data$tabs[[4]] %>% 
                                                                                  filter(party == "rep"),
  "perc_acc_dems",     "out2 ~ group + male + age + white + college",          t_data$tabs[[4]] %>% 
                                                                                  filter(party == "dem"),
  "perc_acc_dems_nct", "out2 ~ group",                                         t_data$tabs[[4]] %>% 
                                                                                  filter(party == "dem")
  ) %>% 
  mutate(
    mods = frm %>% 
      map2(
        tab, 
        \(i, j)
        
        i %>% 
          as.formula %>% 
          lm_robust(
            data = j
          )
      )
    )


# Study 5 -----------------------------------------------------------------

t_s5  <- tribble(
  ~mod_name, ~frm, ~tab,
  "perc_acc_ovr",      "out2 ~ group + party + male + age + white + college",  t_data$tabs[[5]],
  "perc_acc_ovr_nct",  "out2 ~ group",                                         t_data$tabs[[5]],
  "perc_acc_reps",     "out2 ~ group + male + age + white + college",          t_data$tabs[[5]] %>%
                                                                                  filter(party == "rep"),
  "perc_acc_reps_nct", "out2 ~ group",                                         t_data$tabs[[5]] %>%
                                                                                  filter(party == "rep"),
  "perc_acc_dems",     "out2 ~ group + male + age + white + college",          t_data$tabs[[5]] %>%
                                                                                  filter(party == "dem"),
  "perc_acc_dems_nct", "out2 ~ group",                                         t_data$tabs[[5]] %>%
                                                                                  filter(party == "dem")
  ) %>% 
  mutate(
    mods = frm %>% 
      map2(
        tab, 
        \(i, j)
        
        i %>% 
          as.formula %>% 
          lm_robust(
            data = j
          )
      )
    )

I guess I should pause here to ensure that Katie is satisfied I didn’t screw anything up.

So the central thing that we plot in this figure are the conditional expected values (ie, the point ranges in the :

The quantities are extracted from the lm_robust objects here

list(
  t_s3, 
  t_s4, 
  t_s5
  ) %>% 
  map2(
    "study " %>% 
      str_c(3:5),
    \(i, j)
    
    i %>% 
      filter(
        mod_name %>% 
          equals("perc_acc_ovr")
      ) %>% 
      pluck("mods", 1) %>% 
      emmeans(~group) %>% 
      tidy %>% 
      select(1:2) %>% 
      mutate(
        study = j
      )
  ) %>% 
  list_rbind %>% 
  spread(
    study, estimate
  ) %>% 
  slice(
    c(
      2:4, 5, 1
    )
  )
## # A tibble: 5 × 4
##   group    `study 3` `study 4` `study 5`
##   <chr>        <dbl>     <dbl>     <dbl>
## 1 control       5.13      4.88      5.17
## 2 early         5.23     NA        NA   
## 3 late          5.33     NA        NA   
## 4 standard     NA         5.01      5.14
## 5 best         NA         2.86      3.11

And then, the differences in the right facet are just different comparisons of these values.

lemme know if we’re on the same page!