bootstrap per condition and plot with ggdist.R

library(tidyverse)
library(ggdist) # for stat_histinterval()

############# make the dataset ###################

# 1 sample per subject per condition
data_per_subject_condition = expand_grid(
  subjectID = paste0("S", 1:50),
  independent_variable_A = c("circle", "square"),
  independent_variable_B = 1:3
) %>% 
  mutate(speed = (independent_variable_A=="circle") * (independent_variable_B / 5) + rnorm(n(), 0, 0.5))


# normalize values relative to each subject's performance
baseline = data_per_subject_condition %>% 
  filter(independent_variable_A == "circle") %>% 
  filter(independent_variable_B == 1) %>% 
  select(-independent_variable_A, -independent_variable_B) %>% 
  rename(baseline = speed)

data_per_subject_condition = data_per_subject_condition %>% 
  left_join(baseline) %>% 
  mutate(relative_speed = speed - baseline)

# sanity check (ignores subjects)
ggplot(data_per_subject_condition) + aes(x=independent_variable_B, y=speed, color=independent_variable_A) + geom_smooth()


################ bootstrapped means #################

BOOT_COUNT = 1000

# put into purrr's nested table format and make bootstrap copies
boots = data_per_subject_condition %>% 
  ## reshape dependent variables into separate rows
  # pivot_longer(accuracy:speed, "dependent_variable") %>% 
  # get one tibble per condition and DV
  group_by(independent_variable_A, independent_variable_B) %>% 
  nest() %>% 
  # duplicate each based on number of bootstrap iterations
  expand_grid(boot_index = 1:BOOT_COUNT)

# bootstrapped means: randomly sample each set with replacement and compute the mean
boots = boots %>% 
  #randomly sample with replacement
  mutate(data = map(data, ~slice_sample(., prop = 1, replace = TRUE))) %>% 
  # get the mean of your column of choice
  mutate(relative_speed = map(data, ~mean(.$relative_speed))) %>% 
  mutate(data = NULL) %>% 
  unnest(relative_speed)


################# graph it ######################

boots %>% filter(independent_variable_A != "circle" | independent_variable_B != 1) %>% 
ggplot() +
  aes(x = relative_speed, y = paste(independent_variable_A, independent_variable_B), fill = paste(independent_variable_A, independent_variable_B)) +
  geom_vline(xintercept = 0, linetype="dotted") +
  stat_histinterval(.width=c(.95, .99), breaks=seq(-1, 1, 0.05), point_size=0, interval_size_range = c(1.5,3), justification = -0.1, scale=0.8) +
  scale_fill_brewer(palette = "Set1") +
  theme_classic(20) + theme(axis.ticks.x = element_blank()) +
  guides(fill = "none") +
  labs(title = "Relative Speed", x=NULL, y=NULL)