Starter code for live coding using OLS, robust SEs, and {marginaleffects} to analysis data from survey experiments

starter-lm.R

# load packages
library(tidyverse)

# create data frame of the potential outcomes
po <- tribble(
  ~Name,        ~Y1, ~Y0,
  "Alex Smith",     7,   5,
  "Jamie Doe",      2,   3,
  "Pat Johnson",    7,   7,
  "Jordan Lee",     5,   4,
  "Taylor Green",   3,   6,
  "Morgan White",   4,   4,
  "Casey Brown",    5,   3,
  "Drew Wilson",    3,   3,
  "Chris Bailey",   4,   2,
  "Sam Rivera",     1,   3,
  "Jesse Kim",      7,   1,
  "Robin Parker",   5,   3
) 


# random assignment
zeros_and_ones <- rep(0:1, length.out = nrow(po))
W <- sample(zeros_and_ones)

# create (fake) observed data set
observed <- data.frame(
  W_num = W, 
  W_fct = ifelse(W == 1, "Treatment", "Control"),
  Y = po$Y1*W + po$Y0*(1 - W)
)

# estimate variance
treat <- filter(observed, W_num == 1)   # separate data frame for treatment group
control <- filter(observed, W_num == 0) # separate data frame for control group
var(treat$Y)/6 + var(control$Y)/6  # formula to estimate variance

# 90% CI: estimate +/- 1.64*SE
estimate <- mean(treat$Y) - mean(control$Y)
var <- var(treat$Y)/6 + var(control$Y)/6
se <- sqrt(var) # SE = sqrt(variance)
lwr90 <- estimate - 1.64*se
upr90 <- estimate + 1.64*se

# how to do this the easy way
# 1. linear regression (means for conditions, with extensions!)
# 2. robust standard errors (HC2) and 90% CIs
# 3. avg_predictions() in {marginaleffects} (to compute the means for each condition)
# 4. hypotheses() from {maringaleffects} (to *compare* the means)
# lm() + avg_predictions() + hypotheses()