A computational model of responsibility judgments from counterfactual simulations and intention inferences
Sarah A. Wu, Shruti Sridhar, and Tobias Gerstenberg
April 26, 2024
1 Setup
2 Helper functions
2.1 Demographics
print_demographics = function(data) {
# gender
data %>%
group_by(gender) %>%
summarise(n = n()) %>%
print()
# age
print('age:')
mean(data$age, na.rm = T) %>% print()
sd(data$age, na.rm = T) %>% print()
# race
data %>%
group_by(race) %>%
summarise(n = n()) %>%
print()
# ethnicity
data %>%
group_by(ethnicity) %>%
summarise(n = n()) %>%
print()
# time taken
print('time taken (min):')
print(mean(data$time, na.rm = T)/60/1000)
print(sd(data$time, na.rm = T)/60/1000)
}2.2 Scatterplot
scatterplot = function(data, x, x_label, y, y_label, point_fill = 'gray50') {
r = cor(data[[paste0(x, "_mean")]], data[[paste0(y, "_mean")]]) %>%
round(2)
rmse = rmse(data[[paste0(x, "_mean")]], data[[paste0(y, "_mean")]])
g = ggplot(data = data,
aes(x = get(paste0(x, "_mean")),
y = get(paste0(y, "_mean")))) +
geom_abline(intercept = 0, slope = 1,
linetype = 2, size = 0.2) +
# error bars
geom_linerange(size = 0.5,
mapping = aes(ymin = get(paste0(y, "_low")),
ymax = get(paste0(y, "_high"))),
color = 'lightgray') +
geom_errorbarh(mapping = aes(xmin = get(paste0(x, "_low")),
xmax = get(paste0(x, "_high"))),
color = 'lightgray',
height = 0) +
# means
geom_point(shape = 21, size = 3, stroke = 0.1, fill = point_fill) +
geom_text(label = sprintf('RMSE = %.2f\nr = %s', rmse, r),
hjust = 0, # left align
x = 0, y = 95,
size = 5, check_overlap = T) +
scale_x_continuous(name = x_label,
limits = c(0, 100)) +
scale_y_continuous(name = y_label,
limits = c(0, 100)) +
theme(legend.position = "none",
text = element_text(size = 16),
axis.title.x = element_markdown(size = 15),
axis.title.y = element_markdown(size = 15))
return(g)
}2.3 Barplot
Experiment conditions:
conditions = c('cf', 'int', 'effort', 'resp', 'resp_blue', 'resp_red')
condition_labels = c('cf' = 'counterfactual',
'int' = 'intention',
'effort' = 'effort',
'resp' = 'responsibility',
'resp_blue' = 'responsibility (blue)',
'resp_red' = 'responsibility (red)')
condition_colors = c('cf' = 'orchid3',
'int' = 'tan2',
'effort' = 'aquamarine3',
'resp' = 'deepskyblue',
'resp_blue' = 'deepskyblue',
'resp_red' = 'brown1')Explanation features:
features = c('box', 'time', 'red_actions', 'red_mental', 'blue_actions', 'blue_mental')
feature_labels = c('box' = 'box',
'time' = 'time',
'red_actions' = 'red actions',
'red_mental' = 'red mental',
'blue_actions' = 'blue actions',
'blue_mental' = 'blue mental')
feature_colors = c('box' = 'gray90',
'time' = 'gray60',
'red_actions' = 'salmon',
'red_mental' = 'brown1',
'blue_actions' = 'lightskyblue1',
'blue_mental' = 'deepskyblue')barplot = function(data.all, data.models, labels) {
fill_values = condition_colors[names(condition_colors) %in% data.models$condition]
fill_labels = condition_labels[names(condition_labels) %in% data.models$condition]
g = ggplot(data = data.all,
mapping = aes(x = condition,
y = judgment,
fill = condition)) +
stat_summary(fun = "mean",
geom = "bar",
color = "black",
alpha = 0.8,
size = 0.1) +
stat_summary(fun.data = mean_cl_boot,
geom = "linerange",
color = "black") +
geom_point(size = 0.1,
shape = 21,
color = 'gray40',
position = position_jitter(height = 0, width = 0.2),
show.legend = F) +
geom_point(data = data.models,
mapping = aes(y = prediction),
size = 3,
shape = 21,
show.legend = F) +
geom_custom(data = labels %>%
mutate(condition = NA),
mapping = aes(data = grob,
x = -Inf, y = Inf),
grob_fun = function(x) rasterGrob(x,
interpolate = T,
vjust = 0,
hjust = 0,
height = unit(3.8, 'cm'))) +
geom_text(aes(label = trial_label,
x = 0.6, y = 126),
color = 'black',
size = 7,
hjust = 0,
check_overlap = T) +
facet_wrap(~ trial, nrow = 1, scales = 'free_x') +
guides(fill = guide_legend(nrow = 1)) +
scale_fill_manual(values = fill_values,
labels = fill_labels) +
scale_y_continuous(name = 'judgment',
breaks = seq(0, 100, 25),
expand = expansion(mult = 0.01)) +
coord_cartesian(clip = "off",
ylim = c(0, 110)) +
theme(text = element_text(size = 16),
legend.position = 'bottom',
legend.direction = 'horizontal',
legend.text = element_text(
margin = margin(r = 0.5, unit = 'cm')),
axis.title.x = element_blank(),
axis.text.x = element_blank(),
axis.line.x = element_blank(),
axis.ticks.x = element_blank(),
strip.background = element_blank(),
strip.text = element_blank(),
panel.grid.major.y = element_line(),
plot.margin = margin(t = 3.9, l = 0.2,
r = 0.2 , b = 0.2, unit = 'cm'))
print(g)
}3 Experiment 1
data_dir = '../data/experiment1/'
figures_dir = '../figures/experiment1/'
info_dir = '../model/experiments/experiment1/'3.1 Data
3.1.1 Trial info
trials.exp1 = c(1, 9, 14, 16, 17, 21, 23, 24)
info.exp1 = read_csv(paste0(info_dir, 'models.csv'),
show_col_types = F) %>%
rename(cf_model = cf_success_rate,
int_model = numerical_intention,
effort_model = blue_effort) %>%
mutate(outcome = factor(outcome, levels = c('success', 'fail')),
trial_label = c(1, NA, NA, NA, NA, NA, NA, NA,
2, NA, NA, NA, NA, 3, NA, 4,
5, NA, NA, NA, 6, NA, 7, 8),
) %>%
left_join(read_csv(paste0(info_dir, 'heuristics.csv'),
show_col_types = F),
by = 'trial')3.1.2 Counterfactual condition
3.1.2.1 Response data
data.exp1.cf = read_csv(paste0(data_dir, 'counterfactual/trials.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
drop_na(trial) %>%
left_join(info.exp1 %>% select(trial, outcome),
by = 'trial') %>%
mutate(id = factor(id),
plot_id = as.factor(as.integer(id)),
cf_recode = ifelse(outcome == 'success', 100 - cf, cf))Number of participants:
## [1] 503.1.2.2 Participant data
Participants were asked: “What factors influenced how you decided to respond? Do you have any questions or comments regarding the experiment?”
3.1.3 Intention condition
3.1.3.1 Response data
data.exp1.int = read_csv(paste0(data_dir, 'intention/trials.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
drop_na(trial) %>%
left_join(info.exp1 %>% select(trial, outcome),
by = 'trial') %>%
mutate(id = factor(id),
plot_id = as.factor(as.integer(id)),
int_recode = ifelse(outcome == 'success', int, 100 - int))Number of participants:
## [1] 503.1.3.2 Participant data
Participants were asked: “What factors influenced how you decided to respond? Do you have any questions or comments regarding the experiment?”
3.1.4 Effort condition
3.1.4.1 Response data
data.exp1.effort = read_csv(paste0(data_dir, 'effort/trials.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
drop_na(trial) %>%
mutate(id = factor(id),
plot_id = as.factor(as.integer(id)))Number of participants:
## [1] 503.1.4.2 Participant data
Participants were asked: “What factors influenced how you decided to respond? Do you have any questions or comments regarding the experiment?”
participants.exp1.effort = read_csv(paste0(data_dir, 'effort/participants.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
mutate(id = factor(id)) %>%
merge(data.exp1.effort %>%
distinct(id, plot_id),
by = 'id')
participants.exp1.effort %>%
select(plot_id, feedback) %>%
datatable(rownames = F)3.1.5 Responsibility condition
3.1.5.1 Response data
data.exp1.resp = read_csv(paste0(data_dir, 'responsibility/trials.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
drop_na(trial) %>%
mutate(id = factor(id),
plot_id = as.factor(as.integer(id)))Number of participants:
## [1] 503.1.5.2 Participant data
Participants were asked: “What factors influenced how you decided to respond? Do you have any questions or comments regarding the experiment?”
participants.exp1.resp = read_csv(paste0(data_dir, 'responsibility/participants.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
mutate(id = factor(id)) %>%
merge(data.exp1.resp %>%
distinct(id, plot_id),
by = 'id')
participants.exp1.resp %>%
select(plot_id, feedback) %>%
datatable(rownames = F)3.2 Demographics
participants.exp1.cf %>%
rbind(participants.exp1.int) %>%
rbind(participants.exp1.effort) %>%
rbind(participants.exp1.resp) %>%
print_demographics()## # A tibble: 5 × 2
## gender n
## <chr> <int>
## 1 agender 1
## 2 Female 100
## 3 Male 88
## 4 Non-binary 9
## 5 <NA> 2
## [1] "age:"
## [1] 34.26
## [1] 12.9313
## # A tibble: 9 × 2
## race n
## <chr> <int>
## 1 American Indian/Alaska Native 2
## 2 Asian 21
## 3 Black/African American 10
## 4 hispanic 1
## 5 Hispanic 1
## 6 Multiracial 9
## 7 Native Hawaiian/Pacific Islander 1
## 8 White 153
## 9 <NA> 2
## # A tibble: 3 × 2
## ethnicity n
## <chr> <int>
## 1 Hispanic 20
## 2 Non-Hispanic 174
## 3 <NA> 6
## [1] "time taken (min):"
## [1] 11.53095
## [1] 6.7095893.3 Models
data.exp1.means = info.exp1 %>%
left_join(data.exp1.cf %>%
group_by(trial) %>%
summarise(cf_recode = mean(cf_recode)),
by = 'trial') %>%
left_join(data.exp1.int %>%
group_by(trial) %>%
summarise(int_recode = mean(int_recode)),
by = 'trial') %>%
left_join(data.exp1.effort %>%
group_by(trial) %>%
summarise(effort = mean(effort)),
by = 'trial') %>%
mutate(effort_model = effort_model * max(effort) / max(effort_model))
info.exp1$effort_model = data.exp1.means$effort_model
data.exp1.modeling = data.exp1.resp %>%
left_join(data.exp1.means,
by = 'trial')3.3.1 Overall models
Fitting models:
fit.exp1.cf = brm(
formula = resp ~ 1 + cf_recode + (1 | id),
data = data.exp1.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp1.cf"
)
fit.exp1.int = brm(
formula = resp ~ 1 + int_recode + (1 | id),
data = data.exp1.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp1.int"
)
fit.exp1.effort = brm(
formula = resp ~ 1 + effort + (1 | id),
data = data.exp1.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp1.int"
)
fit.exp1.cf_int = brm(
formula = resp ~ 1 + cf_recode + int_recode + (1 | id),
data = data.exp1.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp1.cf_int"
)
fit.exp1.cf_effort = brm(
formula = resp ~ 1 + cf_recode + effort + (1 | id),
data = data.exp1.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp1.cf_effort"
)
fit.exp1.heuristic = brm(
formula = resp ~ 1 + outcome + red_steps + blue_steps + box_steps + (1 | id),
data = data.exp1.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp1.heuristic"
)Comparing models with leave-one-out cross-validation:
fit.exp1.cf = add_criterion(
fit.exp1.cf,
criterion = "loo",
reloo = T,
file = "cache/fit.exp1.cf")
fit.exp1.int = add_criterion(
fit.exp1.int,
criterion = "loo",
reloo = T,
file = "cache/fit.exp1.int")
fit.exp1.cf_int = add_criterion(
fit.exp1.cf_int,
criterion = "loo",
reloo = T,
file = "cache/fit.exp1.cf_int")
fit.exp1.cf_effort = add_criterion(
fit.exp1.cf_effort,
criterion = "loo",
reloo = T,
file = "cache/fit.exp1.cf_effort")
fit.exp1.heuristic = add_criterion(
fit.exp1.heuristic,
criterion = "loo",
reloo = T,
file = "cache/fit.exp1.heuristic")
loo_compare(fit.exp1.cf,
fit.exp1.int,
fit.exp1.cf_int,
fit.exp1.cf_effort,
fit.exp1.heuristic)## elpd_diff se_diff
## fit.exp1.cf_int 0.0 0.0
## fit.exp1.int -162.0 18.5
## fit.exp1.cf -178.6 21.8
## fit.exp1.cf_effort -179.1 21.9
## fit.exp1.heuristic -445.2 28.13.3.2 Individual models
fit.exp1.cf_individual = brm(
formula = resp ~ 1 + cf_recode,
data = data.exp1.modeling %>%
filter(plot_id == 1),
seed = 1,
save_pars = save_pars(all = T),
file = "cache/fit.exp1.cf_individual")
fit.exp1.int_individual = brm(
formula = resp ~ 1 + int_recode,
data = data.exp1.modeling %>%
filter(plot_id == 1),
seed = 1,
save_pars = save_pars(all = T),
file = "cache/fit.exp1.int_individual")
fit.exp1.cf_int_individual = brm(
formula = resp ~ 1 + cf_recode + int_recode,
data = data.exp1.modeling %>%
filter(plot_id == 1),
seed = 1,
save_pars = save_pars(all = T),
file = "cache/fit.exp1.cf_int_individual")
fit.exp1.cf_effort_individual = brm(
formula = resp ~ 1 + cf_recode + effort,
data = data.exp1.modeling %>%
filter(plot_id == 1),
seed = 1,
save_pars = save_pars(all = T),
file = "cache/fit.exp1.cf_effort_individual")
fit.exp1.heuristic_individual = brm(
formula = resp ~ 1 + outcome + red_steps + blue_steps + box_steps,
data = data.exp1.modeling %>%
filter(plot_id == 1),
seed = 1,
save_pars = save_pars(all = T),
file = "cache/fit.exp1.heuristic_individual")
# update model fits for each participant
individual_fits.exp1 = data.exp1.modeling %>%
group_by(plot_id) %>%
nest() %>%
ungroup() %>%
mutate(fit_cf = map(.x = data,
.f = ~ update(fit.exp1.cf_individual,
newdata = .x,
seed = 1)),
fit_int = map(.x = data,
.f = ~ update(fit.exp1.int_individual,
newdata = .x,
seed = 1)),
fit_cf_int = map(.x = data,
.f = ~ update(fit.exp1.cf_int_individual,
newdata = .x,
seed = 1)),
fit_cf_effort = map(.x = data,
.f = ~ update(fit.exp1.cf_effort_individual,
newdata = .x,
seed = 1)),
fit_heuristic = map(.x = data,
.f = ~ update(fit.exp1.heuristic_individual,
newdata = .x,
seed = 1))) %>%
mutate(fit_cf = map(.x = fit_cf,
.f = ~ add_criterion(.x, criterion = "loo",
moment_match = T)),
fit_int = map(.x = fit_int,
.f = ~ add_criterion(.x, criterion = "loo",
moment_match = T)),
fit_cf_int = map(.x = fit_cf_int,
.f = ~ add_criterion(.x, criterion = "loo",
moment_match = T)),
fit_cf_effort = map(.x = fit_cf_effort,
.f = ~ add_criterion(.x, criterion = "loo",
moment_match = T)),
fit_heuristic = map(.x = fit_heuristic,
.f = ~ add_criterion(.x, criterion = "loo",
moment_match = T)),
model_comparison = pmap(.l = list(cf = fit_cf,
int = fit_int,
cf_int = fit_cf_int,
cf_effort = fit_cf_effort,
heuristic = fit_heuristic),
.f = ~ loo_compare(..1, ..2, ..3, ..4, ..5)),
best_model = map_chr(.x = model_comparison,
.f = ~ rownames(.) %>%
.[1]),
best_model = factor(best_model,
levels = c("..1", "..2", "..3", "..4", "..5"),
labels = c("cf",
"int",
"cf_int",
"cf_effort",
"heuristic")))
save(list = c("individual_fits.exp1"),
file = 'cache/individual_fits.exp1.RData')load(file = 'cache/individual_fits.exp1.RData')
individual_fits.exp1 %>%
count(best_model) %>%
arrange(desc(n))## # A tibble: 4 × 2
## best_model n
## <fct> <int>
## 1 cf_int 38
## 2 cf 6
## 3 cf_effort 4
## 4 int 23.4 Figures
data.exp1.all = data.exp1.resp %>%
select(plot_id, trial, judgment = resp) %>%
mutate(condition = 'resp') %>%
rbind(data.exp1.cf %>%
select(plot_id, trial, judgment = cf) %>%
mutate(condition = 'cf')) %>%
rbind(data.exp1.int %>%
select(plot_id, trial, judgment = int) %>%
mutate(condition = 'int')) %>%
rbind(data.exp1.effort %>%
select(plot_id, trial, judgment = effort) %>%
mutate(condition = 'effort')) %>%
mutate(condition = factor(condition,
levels = conditions)) %>%
left_join(info.exp1 %>% select(trial, trial_label),
by = 'trial') %>%
filter(trial %in% trials.exp1)
data.exp1.models = info.exp1 %>%
select(trial, cf_model, int_model, effort_model) %>%
cbind(fit.exp1.cf_int %>%
fitted(newdata = data.exp1.means,
re_formula = NA) %>%
as_tibble() %>%
select(resp_model = Estimate)) %>%
rename_with(~ sub('_model', '', .x)) %>%
filter(trial %in% trials.exp1) %>%
pivot_longer(cols = -trial,
names_to = 'condition',
values_to = 'prediction') %>%
mutate(condition = factor(condition,
levels = conditions))labels.exp1 = data.frame(trial = trials.exp1) %>%
mutate(grob = map(.x = trial,
.f = ~ readPNG(paste0(figures_dir, 'trial_stills/',
.x, '.png'))))3.4.1 Barplot
3.4.2 Simulation models
plot.exp1.cf_model = scatterplot(
data.exp1.cf %>%
group_by(trial) %>%
summarise(cf_mean = mean(cf),
cf_low = smean.cl.boot(cf)[2],
cf_high = smean.cl.boot(cf)[3]) %>%
left_join(info.exp1 %>%
select(trial, cf_model_mean = cf_model),
by = 'trial') %>%
mutate(cf_model_low = cf_model_mean,
cf_model_high = cf_model_mean),
'cf_model', '<span style="color:orchid3">counterfactual simulation model</span>',
'cf', '<span style="color:orchid3">counterfactual judgment</span>'
)
plot.exp1.int_model = scatterplot(
data.exp1.int %>%
group_by(trial) %>%
summarise(int_mean = mean(int),
int_low = smean.cl.boot(int)[2],
int_high = smean.cl.boot(int)[3]) %>%
left_join(info.exp1 %>%
select(trial, int_model_mean = int_model),
by = 'trial') %>%
mutate(int_model_low = int_model_mean,
int_model_high = int_model_mean),
'int_model', '<span style="color:tan2">intention inference model</span>',
'int', '<span style="color:tan2">intention judgment</span>'
)
plot.exp1.effort_model = scatterplot(
data.exp1.effort %>%
group_by(trial) %>%
summarise(effort_mean = mean(effort),
effort_low = smean.cl.boot(effort)[2],
effort_high = smean.cl.boot(effort)[3]) %>%
left_join(info.exp1 %>%
select(trial, effort_model_mean = effort_model),
by = 'trial') %>%
mutate(effort_model_low = effort_model_mean,
effort_model_high = effort_model_mean),
'effort_model', '<span style="color:aquamarine3">effort model</span>',
'effort', '<span style="color:aquamarine3">effort judgment</span>'
)
plot.exp1.cf_model + plot.exp1.int_model + plot.exp1.effort_model +
theme(plot.tag = element_text(size = 24))
3.4.3 Responsibility scatterplots
data.exp1.resp_models = data.exp1.resp %>%
group_by(trial) %>%
summarise(resp_mean = mean(resp),
resp_low = smean.cl.boot(resp)[2],
resp_high = smean.cl.boot(resp)[3]) %>%
cbind(fit.exp1.cf %>%
fitted(newdata = data.exp1.means,
re_formula = NA) %>%
as_tibble() %>%
select(cf_model = Estimate)) %>%
cbind(fit.exp1.int %>%
fitted(newdata = data.exp1.means,
re_formula = NA) %>%
as_tibble() %>%
select(int_model = Estimate)) %>%
cbind(fit.exp1.cf_int %>%
fitted(newdata = data.exp1.means,
re_formula = NA) %>%
as_tibble() %>%
select(cf_int_model = Estimate)) %>%
cbind(fit.exp1.cf_effort %>%
fitted(newdata = data.exp1.means,
re_formula = NA) %>%
as_tibble() %>%
select(cf_effort_model = Estimate)) %>%
cbind(fit.exp1.heuristic %>%
fitted(newdata = data.exp1.means,
re_formula = NA) %>%
as_tibble() %>%
select(heuristic_model = Estimate)) %>%
# clip model predictions between 0 and 100
mutate_at(vars(ends_with('_model')),
~ pmax(pmin(.x, 100), 0))plot.exp1.cf = scatterplot(
data.exp1.resp_models %>%
select(trial, resp_mean, resp_low, resp_high,
cf_model_mean = cf_model,
cf_model_low = cf_model,
cf_model_high = cf_model),
'cf_model', '<span style="color:orchid3">counterfactual</span>',
'resp', '<span style="color:deepskyblue2">responsibility</span>',
)
plot.exp1.int = scatterplot(
data.exp1.resp_models %>%
select(trial, resp_mean, resp_low, resp_high,
int_model_mean = int_model,
int_model_low = int_model,
int_model_high = int_model),
'int_model', '<span style="color:tan2">intention</span>',
'resp', '<span style="color:deepskyblue2">responsibility</span>',
)
plot.exp1.cf_int = scatterplot(
data.exp1.resp_models %>%
select(trial, resp_mean, resp_low, resp_high,
cf_int_model_mean = cf_int_model,
cf_int_model_low = cf_int_model,
cf_int_model_high = cf_int_model),
'cf_int_model', '<span style="color:orchid3">counterfactual</span> + <span style="color:tan2">intention</span>',
'resp', '<span style="color:deepskyblue2">responsibility</span>',
)
plot.exp1.cf_effort = scatterplot(
data.exp1.resp_models %>%
select(trial, resp_mean, resp_low, resp_high,
cf_effort_model_mean = cf_effort_model,
cf_effort_model_low = cf_effort_model,
cf_effort_model_high = cf_effort_model),
'cf_effort_model', '<span style="color:orchid3">counterfactual</span> + <span style="color:aquamarine3">effort</span>',
'resp', '<span style="color:deepskyblue2">responsibility</span>',
)
plot.exp1.heuristic = scatterplot(
data.exp1.resp_models %>%
select(trial, resp_mean, resp_low, resp_high,
heuristic_model_mean = heuristic_model,
heuristic_model_low = heuristic_model,
heuristic_model_high = heuristic_model),
'heuristic_model', 'heuristic model',
'resp', '<span style="color:deepskyblue2">responsibility</span>',
)
plot.exp1.cf + plot.exp1.int + plot.exp1.cf_int + plot.exp1.cf_effort + plot.exp1.heuristic +
plot_layout(nrow = 1) +
plot_annotation(tag_levels = 'A') &
theme(plot.tag = element_text(size = 24))
4 Experiment 2
data_dir = '../data/experiment2/'
figures_dir = '../figures/experiment2/'
info_dir = '../model/experiments/experiment2/'4.1 Data
4.1.1 Trial info
trials.exp2 = c(23, 24, 1, 2, 19, 20, 21, 22)
info.exp2 = read_csv(paste0(info_dir, 'models.csv'),
show_col_types = F) %>%
rename(cf_model = cf_success_rate,
int_model = numerical_intention) %>%
mutate(outcome = factor(outcome, levels = c('success', 'fail')),
trial_label = c('2A', '2B', NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, '3A', '3B', '4A', '4B', '1A', '1B')
) %>%
left_join(read_csv(paste0(info_dir, 'heuristics.csv'),
show_col_types = F),
by = 'trial')4.1.2 Counterfactual condition
4.1.2.1 Response data
data.exp2.cf = read_csv(paste0(data_dir, 'counterfactual/trials.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
drop_na(trial) %>%
left_join(info.exp2 %>% select(trial, outcome),
by = 'trial') %>%
mutate(id = factor(id),
plot_id = as.factor(as.integer(id)),
cf_recode = ifelse(outcome == 'success', 100 - cf, cf))Number of participants:
## [1] 504.1.2.2 Participant data
Participants were asked: “What factors influenced how you decided to respond? Do you have any questions or comments regarding the experiment?”
4.1.3 Intention condition
4.1.3.1 Response data
data.exp2.int = read_csv(paste0(data_dir, 'intention/trials.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
drop_na(trial) %>%
left_join(info.exp2 %>% select(trial, outcome),
by = 'trial') %>%
mutate(id = factor(id),
plot_id = as.factor(as.integer(id)),
int_recode = ifelse(outcome == 'success', int, 100 - int))Number of participants:
## [1] 504.1.3.2 Participant data
Participants were asked: “What factors influenced how you decided to respond? Do you have any questions or comments regarding the experiment?”
4.1.4 Responsibility condition
4.1.4.1 Response data
data.exp2.resp = read_csv(paste0(data_dir, 'responsibility/trials.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
drop_na(trial) %>%
mutate(id = factor(id),
plot_id = as.factor(as.integer(id)))Number of participants:
## [1] 504.1.4.2 Participant data
Participants were asked: “What factors influenced how you decided to respond? Do you have any questions or comments regarding the experiment?”
participants.exp2.resp = read_csv(paste0(data_dir, 'responsibility/participants.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
mutate(id = factor(id)) %>%
merge(data.exp2.resp %>%
distinct(id, plot_id),
by = 'id')
participants.exp2.resp %>%
select(plot_id, feedback) %>%
datatable(rownames = F)4.1.5 Explanation condition
4.1.5.1 Response data
data.exp2.expl_raw = read_csv(paste0(data_dir, 'explanation/trials.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
drop_na(trial) %>%
mutate(id = factor(id),
plot_id = as.factor(as.integer(id)),
explanation = str_to_sentence(explanation))Number of participants:
## [1] 50Coded features:
data.exp2.expl = read_csv(paste0(data_dir, 'explanation/trials_coded.csv'),
show_col_types = F) %>%
pivot_longer(cols = -trial,
names_to = 'feature',
values_to = 'mentioned') %>%
mutate(feature = factor(feature, levels = features)) %>%
left_join(info.exp2 %>% select(trial, trial_label),
by = 'trial')4.1.5.2 Participant data
Participants were asked: “What factors influenced how you decided to respond? Do you have any questions or comments regarding the experiment?”
participants.exp2.expl = read_csv(paste0(data_dir, 'explanation/participants.csv'),
show_col_types = F) %>%
rename(id = workerid) %>%
mutate(id = factor(id)) %>%
merge(data.exp2.expl_raw %>%
distinct(id, plot_id),
by = 'id')
participants.exp2.expl %>%
select(plot_id, feedback) %>%
datatable(rownames = F)4.2 Demographics
participants.exp2.cf %>%
rbind(participants.exp2.int) %>%
rbind(participants.exp2.resp) %>%
rbind(participants.exp2.expl) %>%
print_demographics()## # A tibble: 5 × 2
## gender n
## <chr> <int>
## 1 Female 98
## 2 Male 93
## 3 Non-binary 7
## 4 trans 1
## 5 <NA> 1
## [1] "age:"
## [1] 36.01508
## [1] 12.37972
## # A tibble: 7 × 2
## race n
## <chr> <int>
## 1 American Indian/Alaska Native 3
## 2 Asian 18
## 3 Black/African American 16
## 4 Hispanic 1
## 5 hispanic/latino 1
## 6 Multiracial 13
## 7 White 148
## # A tibble: 3 × 2
## ethnicity n
## <chr> <int>
## 1 Hispanic 18
## 2 Non-Hispanic 172
## 3 <NA> 10
## [1] "time taken (min):"
## [1] 20.77965
## [1] 12.449454.3 Models
4.3.1 Responsibility (blue)
data.exp2.means = info.exp2 %>%
left_join(data.exp2.cf %>%
group_by(trial) %>%
summarise(cf_recode = mean(cf_recode)),
by = 'trial') %>%
left_join(data.exp2.int %>%
group_by(trial) %>%
summarise(int_recode = mean(int_recode)),
by = 'trial')
data.exp2.modeling = data.exp2.resp %>%
select(-resp_red) %>%
left_join(data.exp2.means, by = 'trial')4.3.1.1 Overall models
Fitting models:
fit.exp2.blue.cf = brm(
formula = resp_blue ~ 1 + cf_recode + (1 | id),
data = data.exp2.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp2.blue.cf"
)
fit.exp2.blue.int = brm(
formula = resp_blue ~ 1 + int_recode + (1 | id),
data = data.exp2.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp2.blue.int"
)
fit.exp2.blue.cf_int = brm(
formula = resp_blue ~ 1 + cf_recode + int_recode + (1 | id),
data = data.exp2.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp2.blue.cf_int"
)
fit.exp2.blue.heuristic = brm(
formula = resp_blue ~ 1 + outcome + red_steps + blue_steps + box_steps + (1 | id),
data = data.exp2.modeling,
iter = 2000,
seed = 1,
file = "cache/fit.exp2.blue.heuristic"
)Comparing models with leave-one-out cross-validation:
fit.exp2.blue.cf = add_criterion(
fit.exp2.blue.cf,
criterion = "loo",
reloo = T,
file = "cache/fit.exp2.blue.cf")
fit.exp2.blue.int = add_criterion(
fit.exp2.blue.int,
criterion = "loo",
reloo = T,
file = "cache/fit.exp2.blue.int")
fit.exp2.blue.cf_int = add_criterion(
fit.exp2.blue.cf_int,
criterion = "loo",
reloo = T,
file = "cache/fit.exp2.blue.cf_int")
fit.exp2.blue.heuristic = add_criterion(
fit.exp2.blue.heuristic,
criterion = "loo",
reloo = T,
file = "cache/fit.exp2.blue.heuristic")
loo_compare(fit.exp2.blue.cf,
fit.exp2.blue.int,
fit.exp2.blue.cf_int,
fit.exp2.blue.heuristic)## elpd_diff se_diff
## fit.exp2.blue.cf_int 0.0 0.0
## fit.exp2.blue.int -92.1 13.8
## fit.exp2.blue.cf -142.1 17.5
## fit.exp2.blue.heuristic -350.5 25.44.3.1.2 Individual models
fit.exp2.blue.cf_individual = brm(
formula = resp_blue ~ 1 + cf_recode,
data = data.exp2.modeling %>%
filter(plot_id == 1),
seed = 1,
save_pars = save_pars(all = T),
file = "cache/fit.exp2.blue.cf_individual")
fit.exp2.blue.int_individual = brm(
formula = resp_blue ~ 1 + int_recode,
data = data.exp2.modeling %>%
filter(plot_id == 1),
seed = 1,
save_pars = save_pars(all = T),
file = "cache/fit.exp2.blue.int_individual")
fit.exp2.blue.cf_int_individual = brm(
formula = resp_blue ~ 1 + cf_recode + int_recode,
data = data.exp2.modeling %>%
filter(plot_id == 1),
seed = 1,
save_pars = save_pars(all = T),
file = "cache/fit.exp2.blue.cf_int_individual")
fit.exp2.blue.heuristic_individual = brm(
formula = resp_blue ~ 1 + outcome + red_steps + blue_steps + box_steps,
data = data.exp2.modeling %>%
filter(plot_id == 1),
seed = 1,
save_pars = save_pars(all = T),
file = "cache/fit.exp2.blue.heuristic_individual")
# update model fits for each participant
individual_fits.exp2 = data.exp2.modeling %>%
group_by(plot_id) %>%
nest() %>%
ungroup() %>%
mutate(fit_cf = map(.x = data,
.f = ~ update(fit.exp2.blue.cf_individual,
newdata = .x,
seed = 1)),
fit_int = map(.x = data,
.f = ~ update(fit.exp2.blue.int_individual,
newdata = .x,
seed = 1)),
fit_cf_int = map(.x = data,
.f = ~ update(fit.exp2.blue.cf_int_individual,
newdata = .x,
seed = 1)),
fit_heuristic = map(.x = data,
.f = ~ update(fit.exp2.blue.heuristic_individual,
newdata = .x,
seed = 1))) %>%
mutate(fit_cf = map(.x = fit_cf,
.f = ~ add_criterion(.x, criterion = "loo",
moment_match = T)),
fit_int = map(.x = fit_int,
.f = ~ add_criterion(.x, criterion = "loo",
moment_match = T)),
fit_cf_int = map(.x = fit_cf_int,
.f = ~ add_criterion(.x, criterion = "loo",
moment_match = T)),
fit_heuristic = map(.x = fit_heuristic,
.f = ~ add_criterion(.x, criterion = "loo",
moment_match = T)),
model_comparison = pmap(.l = list(cf = fit_cf,
int = fit_int,
cf_int = fit_cf_int,
heuristic = fit_heuristic),
.f = ~ loo_compare(..1, ..2, ..3, ..4)),
best_model = map_chr(.x = model_comparison,
.f = ~ rownames(.) %>%
.[1]),
best_model = factor(best_model,
levels = c("..1", "..2", "..3", "..4"),
labels = c("cf",
"int",
"cf_int",
"heuristic")))
save(list = c("individual_fits.exp2"),
file = 'cache/individual_fits.exp2.RData')load(file = 'cache/individual_fits.exp2.RData')
individual_fits.exp2 %>%
count(best_model) %>%
arrange(desc(n))## # A tibble: 4 × 2
## best_model n
## <fct> <int>
## 1 cf_int 26
## 2 int 15
## 3 cf 7
## 4 heuristic 24.4 Figures
data.exp2.all = data.exp2.resp %>%
select(-id) %>%
pivot_longer(cols = starts_with('resp'),
names_to = 'condition',
values_to = 'judgment') %>%
rbind(data.exp2.cf %>%
select(plot_id, trial, judgment = cf) %>%
mutate(condition = 'cf')) %>%
rbind(data.exp2.int %>%
select(plot_id, trial, judgment = int) %>%
mutate(condition = 'int')) %>%
mutate(condition = factor(condition, levels = conditions)) %>%
left_join(info.exp2 %>% select(trial, trial_label),
by = 'trial') %>%
filter(trial %in% trials.exp2) %>%
mutate(trial = as.integer(factor(trial, levels = trials.exp2)))
data.exp2.models = info.exp2 %>%
select(trial, cf_model, int_model) %>%
cbind(fit.exp2.blue.cf_int %>%
fitted(newdata = data.exp2.means,
re_formula = NA) %>%
as_tibble() %>%
select(resp_blue_model = Estimate)) %>%
left_join(data.exp2.resp %>%
cbind(fit.exp2.red %>%
fitted(newdata = data.exp2.resp) %>%
as_tibble()) %>%
group_by(trial) %>%
summarise(resp_red_model = mean(Estimate)),
by = 'trial') %>%
rename_with(~ sub('_model', '', .x)) %>%
filter(trial %in% trials.exp2) %>%
pivot_longer(cols = -trial,
names_to = 'condition',
values_to = 'prediction') %>%
mutate(condition = factor(condition, levels = conditions),
trial = as.integer(factor(trial, levels = trials.exp2)))labels.exp2 = data.frame(trial = trials.exp2) %>%
mutate(grob = map(.x = trial,
.f = ~ readPNG(paste0(figures_dir, 'trial_stills/',
.x, '.png'))),
trial = as.integer(factor(trial,
levels = trials.exp2)))4.4.1 Barplot
4.4.2 Simulation models
plot.exp2.cf_model = scatterplot(
data.exp2.cf %>%
group_by(trial) %>%
summarise(cf_mean = mean(cf),
cf_low = smean.cl.boot(cf)[2],
cf_high = smean.cl.boot(cf)[3]) %>%
left_join(info.exp2 %>%
select(trial, cf_model_mean = cf_model),
by = 'trial') %>%
mutate(cf_model_low = cf_model_mean,
cf_model_high = cf_model_mean),
'cf_model', '<span style="color:orchid3">counterfactual simulation model</span>',
'cf', '<span style="color:orchid3">counterfactual judgment</span>'
)
plot.exp2.int_model = scatterplot(
data.exp2.int %>%
group_by(trial) %>%
summarise(int_mean = mean(int),
int_low = smean.cl.boot(int)[2],
int_high = smean.cl.boot(int)[3]) %>%
left_join(info.exp2 %>%
select(trial, int_model_mean = int_model),
by = 'trial') %>%
mutate(int_model_low = int_model_mean,
int_model_high = int_model_mean),
'int_model', '<span style="color:tan2">intention inference model</span>',
'int', '<span style="color:tan2">intention judgment</span>'
)
plot.exp2.cf_model + plot.exp2.int_model +
theme(plot.tag = element_text(size = 24))
4.4.3 Responsibility scatterplots
data.exp2.resp_blue_models = data.exp2.resp %>%
group_by(trial) %>%
summarise(resp_mean = mean(resp_blue),
resp_low = smean.cl.boot(resp_blue)[2],
resp_high = smean.cl.boot(resp_blue)[3]) %>%
cbind(fit.exp2.blue.cf %>%
fitted(newdata = data.exp2.means,
re_formula = NA) %>%
as_tibble() %>%
select(cf_model = Estimate)) %>%
cbind(fit.exp2.blue.int %>%
fitted(newdata = data.exp2.means,
re_formula = NA) %>%
as_tibble() %>%
select(int_model = Estimate)) %>%
cbind(fit.exp2.blue.cf_int %>%
fitted(newdata = data.exp2.means,
re_formula = NA) %>%
as_tibble() %>%
select(cf_int_model = Estimate)) %>%
cbind(fit.exp2.blue.heuristic %>%
fitted(newdata = data.exp2.means,
re_formula = NA) %>%
as_tibble() %>%
select(heuristic_model = Estimate)) %>%
# clip model predictions between 0 and 100
mutate_at(vars(ends_with('_model')),
~ pmax(pmin(.x, 100), 0))# responsibility for blue
plot.exp2.blue.cf = scatterplot(
data.exp2.resp_blue_models %>%
select(trial, resp_mean, resp_low, resp_high,
cf_model_mean = cf_model,
cf_model_low = cf_model,
cf_model_high = cf_model),
'cf_model', '<span style="color:orchid3">counterfactual</span>',
'resp', '<span style="color:deepskyblue2">responsibility (blue)</span>',
)
plot.exp2.blue.int = scatterplot(
data.exp2.resp_blue_models %>%
select(trial, resp_mean, resp_low, resp_high,
int_model_mean = int_model,
int_model_low = int_model,
int_model_high = int_model),
'int_model', '<span style="color:tan2">intention</span>',
'resp', '<span style="color:deepskyblue2">responsibility (blue)</span>',
)
plot.exp2.blue.cf_int = scatterplot(
data.exp2.resp_blue_models %>%
select(trial, resp_mean, resp_low, resp_high,
cf_int_model_mean = cf_int_model,
cf_int_model_low = cf_int_model,
cf_int_model_high = cf_int_model),
'cf_int_model', '<span style="color:orchid3">counterfactual</span> + <span style="color:tan2">intention</span>',
'resp', '<span style="color:deepskyblue2">responsibility (blue)</span>',
)
plot.exp2.blue.heuristic = scatterplot(
data.exp2.resp_blue_models %>%
select(trial, resp_mean, resp_low, resp_high,
heuristic_model_mean = heuristic_model,
heuristic_model_low = heuristic_model,
heuristic_model_high = heuristic_model),
'heuristic_model', 'heuristic model',
'resp', '<span style="color:deepskyblue2">responsibility (blue)</span>',
)
# responsibility for red
plot.exp2.red = scatterplot(
data.exp2.resp %>%
group_by(trial) %>%
summarise(resp_mean = mean(resp_red),
resp_low = smean.cl.boot(resp_red)[2],
resp_high = smean.cl.boot(resp_red)[3]) %>%
left_join(data.exp2.resp %>%
cbind(fit.exp2.red %>%
fitted(newdata = data.exp2.resp) %>%
as_tibble()) %>%
group_by(trial) %>%
summarise(model_mean = mean(Estimate)),
by = 'trial') %>%
mutate(model_low = model_mean,
model_high = model_mean),
'model', '<span style="color:deepskyblue2">responsibility (blue)</span>',
'resp', '<span style="color:brown1">responsibility (red)</span>',
)
plot.exp2.blue.cf + plot.exp2.blue.int + plot.exp2.blue.cf_int + plot.exp2.blue.heuristic + plot.exp2.red +
plot_layout(nrow = 1) +
plot_annotation(tag_levels = 'A') &
theme(plot.tag = element_text(size = 24))
4.4.4 Explanations
plot.exp2.expl_aggregate = ggplot(data = data.exp2.expl,
mapping = aes(x = feature,
y = mentioned,
fill = feature)) +
stat_summary(fun = "mean",
geom = "bar",
color = "black",
size = 0.1)
plot.exp2.expl_subset = ggplot(data = data.exp2.expl %>%
filter(trial %in% c(20, 22)),
mapping = aes(x = feature,
y = mentioned,
fill = feature)) +
stat_summary(fun = "mean",
geom = "bar",
color = "black",
size = 0.1) +
facet_wrap(~ trial_label, nrow = 1)
plot.exp2.expl_aggregate + plot.exp2.expl_subset +
plot_annotation(tag_levels = 'A') +
plot_layout(nrow = 1, widths = c(1, 2),
guides = 'collect') &
scale_fill_manual(values = feature_colors,
labels = feature_labels) &
scale_y_continuous(name = 'frequency',
limits = c(0, 1),
breaks = seq(0, 1, 0.2),
expand = expansion(mult = 0.01)) &
guides(fill = guide_legend(byrow = T)) &
theme(text = element_text(size = 16),
legend.position = 'right',
legend.title = element_blank(),
legend.spacing.y = unit(0.2, 'cm'),
axis.title.x = element_blank(),
axis.text.x = element_blank(),
axis.line.x = element_blank(),
axis.ticks.x = element_blank(),
strip.text = element_text(size = 18),
panel.grid.major.y = element_line(),
plot.tag = element_text(size = 24))
5 Session info
## R version 4.1.3 (2022-03-10)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur/Monterey 10.16
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] grid stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] forcats_0.5.1 stringr_1.4.0 dplyr_1.0.8 purrr_0.3.4
## [5] readr_2.1.2 tidyr_1.2.0 tibble_3.1.6 tidyverse_1.3.2
## [9] ggtext_0.1.2 egg_0.4.5 gridExtra_2.3 png_0.1-7
## [13] patchwork_1.1.1 Hmisc_4.6-0 ggplot2_3.4.3 Formula_1.2-4
## [17] survival_3.2-13 lattice_0.20-45 Metrics_0.1.4 DT_0.22
## [21] brms_2.17.0 Rcpp_1.0.8.3 ggstatsplot_0.9.1 knitr_1.38
##
## loaded via a namespace (and not attached):
## [1] readxl_1.4.0 backports_1.4.1 plyr_1.8.7
## [4] igraph_1.3.0 splines_4.1.3 crosstalk_1.2.0
## [7] TH.data_1.1-0 rstantools_2.2.0 inline_0.3.19
## [10] digest_0.6.29 htmltools_0.5.2 fansi_1.0.3
## [13] magrittr_2.0.3 checkmate_2.0.0 googlesheets4_1.0.0
## [16] paletteer_1.4.0 cluster_2.1.2 tzdb_0.3.0
## [19] modelr_0.1.8 RcppParallel_5.1.5 matrixStats_0.61.0
## [22] vroom_1.5.7 xts_0.12.1 sandwich_3.0-1
## [25] prettyunits_1.1.1 jpeg_0.1-9 colorspace_2.0-3
## [28] rvest_1.0.2 haven_2.4.3 xfun_0.30
## [31] callr_3.7.0 crayon_1.5.1 jsonlite_1.8.0
## [34] zeallot_0.1.0 zoo_1.8-9 glue_1.6.2
## [37] gargle_1.2.0 gtable_0.3.0 emmeans_1.8.2
## [40] statsExpressions_1.3.1 distributional_0.3.0 pkgbuild_1.3.1
## [43] rstan_2.21.3 abind_1.4-5 scales_1.2.0
## [46] mvtnorm_1.1-3 DBI_1.1.2 miniUI_0.1.1.1
## [49] gridtext_0.1.5 xtable_1.8-4 performance_0.9.0
## [52] htmlTable_2.4.0 bit_4.0.4 foreign_0.8-82
## [55] stats4_4.1.3 StanHeaders_2.21.0-7 httr_1.4.2
## [58] datawizard_0.4.0 htmlwidgets_1.5.4 threejs_0.3.3
## [61] RColorBrewer_1.1-3 posterior_1.2.1 ellipsis_0.3.2
## [64] pkgconfig_2.0.3 reshape_0.8.8 loo_2.5.1
## [67] farver_2.1.0 dbplyr_2.1.1 nnet_7.3-17
## [70] sass_0.4.1 utf8_1.2.2 labeling_0.4.2
## [73] tidyselect_1.1.2 rlang_1.1.1 reshape2_1.4.4
## [76] later_1.3.0 cellranger_1.1.0 munsell_0.5.0
## [79] tools_4.1.3 cli_3.6.1 generics_0.1.2
## [82] broom_0.8.0 ggridges_0.5.3 evaluate_0.15
## [85] fastmap_1.1.0 yaml_2.3.5 rematch2_2.1.2
## [88] bit64_4.0.5 fs_1.5.2 processx_3.5.3
## [91] WRS2_1.1-3 nlme_3.1-155 mime_0.12
## [94] xml2_1.3.3 correlation_0.8.0 compiler_4.1.3
## [97] bayesplot_1.9.0 shinythemes_1.2.0 rstudioapi_0.13
## [100] reprex_2.0.1 bslib_0.3.1 stringi_1.7.6
## [103] highr_0.9 ps_1.6.0 parameters_0.17.0
## [106] Brobdingnag_1.2-9 Matrix_1.5-1 markdown_1.1
## [109] shinyjs_2.1.0 tensorA_0.36.2 vctrs_0.6.2
## [112] pillar_1.7.0 lifecycle_1.0.3 mc2d_0.1-21
## [115] jquerylib_0.1.4 bridgesampling_1.1-2 estimability_1.4.1
## [118] data.table_1.14.2 insight_0.17.0 httpuv_1.6.5
## [121] R6_2.5.1 latticeExtra_0.6-29 bookdown_0.26
## [124] promises_1.2.0.1 codetools_0.2-18 colourpicker_1.1.1
## [127] MASS_7.3-55 gtools_3.9.2 assertthat_0.2.1
## [130] withr_2.5.0 shinystan_2.6.0 multcomp_1.4-18
## [133] hms_1.1.1 bayestestR_0.11.5 parallel_4.1.3
## [136] rpart_4.1.16 coda_0.19-4 rmarkdown_2.13
## [139] googledrive_2.0.0 lubridate_1.8.0 shiny_1.7.1
## [142] base64enc_0.1-3 dygraphs_1.1.1.6