knitr::opts_chunk$set(echo = TRUE)

# load all packages
library(tidyverse)
library(here) # see https://www.tidyverse.org/blog/2017/12/workflow-vs-script/
library(stringr)
library(openxlsx)
library(numbers)
library(dplyr)
library(lubridate)
library(gridExtra) # for plotting multiple plots together
library(clinfun) # for Jonckheere test
library(kableExtra) # for tables in Rmarkdown
library(tableone)  # for Table 2
library(stargazer) # for regression tables
library(multiwayvcov) #  for robust standard errors
library(lmtest) # for robust standard errors
library(sandwich) # for robust standard errors
here <- here::here # because lubridate also has a here function

# overwrite select and filter functions
select <- dplyr::select
filter <- dplyr::filter

# read main data file
df <- read.csv(here("clean data", "all_batches_public.csv"), 1)
library(Rcpp)
ra_coding <- read.xlsx(here("clean data", "ra_coding.xlsx"), sheet = 1) %>% select(-treatment)

# merge ra_coding file with main df
df <- merge(ra_coding, df, by = "id_otree", all.y = TRUE)
rm(ra_coding)

#### Data prep

# time in minutes, with decimals
df <- df %>% 
  mutate(time_h = hour(hms(time_total)) * 3600,
         time_m = minute(hms(time_total)) * 60,
         time_s = second(hms(time_total)),
         time_ = time_m + time_s,  # dropping time in hours because all with 1-24 hours completion time were finished by the system
         clean_time = time_ / 60,
         clean_time = case_when(invalid  == "auto_advanced" ~ NA_real_,
                                invalid == "duplicate" ~ NA_real_,
                                TRUE ~ clean_time)) 

# clean hint vars
df <- df %>% 
  mutate(follow_hint1 = case_when(hint1 == decision1 ~ 1, is.na(hint1) ~ NA_real_, TRUE ~ 0),
         follow_hint2 = case_when(hint2 == decision2 ~ 1, is.na(hint2) ~ NA_real_, TRUE ~ 0))

# generate deviation from social welfare maximum (group level)
df <- df %>% 
  mutate(deviationwelfare1 = group_payoff1 - maximum1, 
         deviationwelfare2 = group_payoff2 - maximum2)

# dummy for duplicates and bots
df <- df %>% 
  mutate(duplicate = case_when(invalid == "duplicate" ~ 1, TRUE ~ 0),
         bot = case_when(invalid == "auto_advanced" ~ 1, TRUE ~ 0))

# dummy for female (collapsing the non-binary and rather not says into "Male"
df <- df %>% 
  mutate(female = case_when(gender == "Female" ~ 1, is.na(gender) ~ NA_real_, TRUE ~ 0)) 

# chose safe option A
df <- df %>% 
  mutate(choseA1 = case_when(decision1 == "A" ~ 1, TRUE ~ 0),
         choseA2 = case_when(decision2 == "A" ~ 1, TRUE ~ 0))

# bayesian deviation
df <- df %>% 
  mutate(bayesian_original = bayesian,
         bayesian = case_when(bayesian == 999999 ~ NA_real_,
                              TRUE ~ bayesian),
         bayesian_abs = abs(60-bayesian),  # absolute difference
         bayesian_qdr = (60-bayesian)^2)   # quadratic, punishing far-off deviations more

# instructions check
df <- df %>% 
  mutate(instr_check = case_when(instr_opened == 0 ~ 0, # never checked instructions
                                 is.na(instr_opened) ~ NA_real_, # never got to instr check phase
                                 TRUE ~ 1), # people checked instructions at least once
         instr_fail = case_when(is.na(check1) ~ NA_real_, # did not get to instr check phase
                                check1 == 0 & check2 == 0 & check3 == 0 & check4 == 0 ~ 0, # not failed check once
                                TRUE ~ 1)) # failed check at least once))

# attention check
df <- df %>% 
  mutate(attention_fail = case_when(attention == 3 ~ 0,
                                    is.na(attention) ~ NA_real_,
                                    TRUE ~ 1))

# risk aversion and experimenter demand
df <- df %>% 
  mutate(box1 = as.numeric(box1),
         box2 = as.numeric(box2),
         riskaverse = 100 - box1,
         exp_demand = box2 - box1)

# clean coding of RA's
df <- df %>% 
  mutate(agreed = case_when(agreed == "Follow Hint" ~ "Follow hint",
                            agreed == "Other " ~ "Other",
                            agreed == "Safe Option" ~ "Safe option",
                            agreed == "pick highest number" ~ "Pick highest number",
                            agreed == "Risky Option" ~ "Risky option",
                            agreed == "Missunderstanding" ~ "Misunderstanding",
                            agreed == "consider hint + other" ~ "Consider hint + Other",
                            agreed == "Ignore hint + other" ~ "Ignore hint + Other",
                            agreed == "ignore hint + other" ~ "Ignore hint + Other",
                            agreed == "ignore hint + safe option" ~ "Ignore hint + Safe option",
                            agreed == "ignore hint + risky option" ~ "Ignore hint + Risky option",
                            agreed == "Consider hint+ Consider risks/benefits" ~ "Consider hint + Consider risks/benefits",
                            agreed == "Consider hint + risky option" ~ "Consider hint + Risky option",
                            agreed == "Consider hint+ risky option" ~ "Consider hint + Risky option",
                            agreed == "Consider hint + safe option" ~ "Consider hint + Safe option",
                            agreed == "Ignore hint +  Safe option" ~ "Ignore hint + Safe option",
                            agreed == "Ignore hint + consider risks/benefits" ~ "Ignore hint + Consider risks/benefits",
                            agreed == "Ignore hint + safe option" ~ "Ignore hint + Safe option",
                            bot == 1 ~ NA_character_,
                            TRUE ~ agreed
                            ))
# stategy dummies
df <- df %>% 
  mutate(strategy_ignore = case_when(treatment != "Partial Info" ~ NA_real_, 
                                     grepl("ignore", agreed, perl = TRUE, ignore.case = TRUE) ~ 1,
                                     is.na(agreed) ~ NA_real_, 
                                     TRUE ~ 0),
         strategy_consider = case_when(treatment != "Partial Info" ~ NA_real_, 
                                       grepl("consider hint", agreed, perl = TRUE, ignore.case = TRUE) ~ 1,
                                       is.na(agreed) ~ NA_real_, 
                                       TRUE ~ 0),
         strategy_follow = case_when(treatment != "Partial Info" ~ NA_real_,
                                     grepl("follow", agreed, perl = TRUE, ignore.case = TRUE) ~ 1,
                                     is.na(agreed) ~ NA_real_, 
                                     TRUE ~ 0),
         strategy_misunderstanding = case_when(grepl("misunderstanding", agreed, perl = TRUE, ignore.case = TRUE) ~ 1,
                                               is.na(agreed) ~ NA_real_, TRUE ~ 0),
         strategy_altruist = case_when(grepl("altruist", agreed, perl = TRUE, ignore.case = TRUE) ~ 1,
                                               is.na(agreed) ~ NA_real_, TRUE ~ 0),
         risk_report = case_when(grepl("benefit", agreed, perl = TRUE, ignore.case = TRUE) ~ 1,
                                 grepl("option", agreed, perl = TRUE, ignore.case = TRUE) ~ 1, 
                                 TRUE ~ 0),
         expdemand_report = case_when(treatment != "Partial Info" ~ NA_real_, 
                                      grepl("own", agreed, perl = TRUE, ignore.case = TRUE) ~ 1, 
                                      TRUE ~ 0),
         other_report = case_when(grepl("other", agreed, perl = TRUE, ignore.case = TRUE) ~ 1, 
                                  TRUE ~ 0),
         noexp_report = case_when(grepl("explanation", agreed, perl = TRUE, ignore.case = TRUE) ~ 1, 
                        TRUE ~ 0),
         safe_report = case_when(risk_report!= 1 ~ NA_real_,
                                 grepl("safe", agreed, perl = TRUE, ignore.case = TRUE) ~ 1, 
                                 TRUE ~ 0),
         strategy_highestnum = case_when(grepl("number", agreed, perl = TRUE, ignore.case = TRUE) ~ 1, 
                        TRUE ~ 0),
         risky_report = case_when(risk_report!= 1 ~ NA_real_,
                       grepl("risky", agreed, perl = TRUE, ignore.case = TRUE) ~ 1, 
                       TRUE ~ 0), 
         tradeoff_report = case_when(risk_report!= 1 ~ NA_real_,
                       grepl("/", agreed, perl = TRUE, ignore.case = TRUE) ~ 1, 
                       TRUE ~ 0))
                       
# change order of variables
df <- df[c(sort(names(df)))]
firstcols = c("id_otree", "session", "date", "treatment", "group_id")
df <- df[,c(firstcols, setdiff(names(df), firstcols))]
lastcols = c("a_coding", "m_coding")
df <- df[,c(setdiff(names(df), lastcols), lastcols)]
rm(firstcols, lastcols)

In-text statistics

Section 1: Following the recommendation

(page 4)

pct_chose_recommended <- round(prop.table(table(df$follow_hint1))["1"]*100,1)
pct_follow_if_B <- round(prop.table(table(subset(df, hint2 == "B")$follow_hint2))["1"]*100,1)

We find that in general, around 71.5% of subjects choose the object they were recommended; this number is higher whenever the recommended option is less risky, or when it is both individually and collectively advantageous to follow the recommendation.

In cases where individual and group interests conflict, we still observe approximately 38% of subjects follow a recommendation that is not in their individual interest.

Section 4: Group counts

(page 16-17)

groups <- read.csv(here("clean data", "group_counts.csv"), 1) %>% 
  group_by(group_id) %>% 
  filter(row_number()==1) 

groups_total_noinfo <- table(groups$treatment)["No Info"]
groups_total_partial <- table(groups$treatment)["Partial Info"]
groups_total_full <- table(groups$treatment)["Full Info"]

dropall <- groups %>% 
  filter(bot == 0 & duplicate == 0) 

dropall_noinfo <- table(dropall$treatment)["No Info"]
dropall_partial <- table(dropall$treatment)["Partial Info"]
dropall_full <- table(dropall$treatment)["Full Info"]

maingroups <- df %>% 
  group_by(group_id) %>% 
  filter(row_number()==1) 

maingroups_noinfo <- table(maingroups$treatment)["No Info"]
maingroups_partial <- table(maingroups$treatment)["Partial Info"]
maingroups_full <- table(maingroups$treatment)["Full Info"]

In total, we collected data for 650 groups (2600 subjects), with 89 in \(No Info\), 281 in \(Partial Info\), and 280 in \(Full Info\). As pre-registered, we drop all groups consisting of only subjects who did not complete the main part of the experiment (Scenarios 1 and 2). Since dropping groups that contain a combination of some drop-outs and some participants who finished would decrease our power under 65%, then, as pre-registered, we keep these groups for the analysis. For the main analysis we are therefore using data from 603 groups, split 80-262-261 between the No-Partial-Full Info treatments, which is slightly more groups than the pre-registered desired minimum.

Section 4.1: Descriptive Statistics

(page 17)

minutes <- round(mean(df$clean_time, na.rm = TRUE),0)
female <- round(prop.table(table((df$gender)))["Female"] * 100, 1)
students <- round(prop.table(table((df$work_student)))["1"] * 100, 1)
earnings <- round(summary(df$payoff_total)["Mean"],2) + 1.50 # which is the showup-fee

rechecked <- round(prop.table(table((df$instr_check)))["1"] * 100, 1)
passed_all <- 100 - round(prop.table(table((df$instr_fail)))["1"] * 100, 1)
df <- df %>% 
  mutate(check_total = check1 + check2 + check3 + check4)
failed_once <- round(table((df$check_total))["1"]/table((df$instr_fail))["1"] *100, 1)

failed_end <- round(prop.table(table((df$attention_fail)))["1"] * 100, 1)
failed_us <- round(table((df$attention))["0"]/table((df$attention_fail))["1"] *100, 1)

exp_demand <- round(mean(df$exp_demand, na.rm = TRUE), 0)
bayesian_dev <- round(mean(df$bayesian_abs, na.rm = TRUE), 1)
bayesian_dev_pct <- round((100*bayesian_dev)/60, 1)
altruist_pct <- round(prop.table(table((df$redistribution)))["1"] * 100, 1)

# check if misunderstanding (instructions failure) correlates with other observables
vars <- c("clean_time", "riskaverse", "exp_demand", "bayesian_abs", "female", "work_student", "instr_check", "instr_fail", "attention_fail", "redistribution")
corr_check <- df %>%  select(vars) 
# cor(corr_check, use = "complete.obs")
# no correlation with instr_fail larger than 0.10 (largest is with attention fail 0.094)
largest_cor <- round(cor(corr_check$instr_fail, corr_check$attention_fail, use = "complete.obs"),3)

Looking at subject characteristics, around 42.3% of our subjects were women, and 7.6% students. The experiment took approximately 14 minutes, and average earnings were 2.36GBP. Around 9.4% of subjects re-checked the instructions after reading them for the first time, but only 36.7% of subjects passed all comprehension checks at the first attempt. All subjects passed all four comprehension questions before moving on to the experiment, and a majority (57.9%) of subjects who failed at least one comprehension check, failed one check, once. Failing at least one comprehension check does not correlate strongly with other subject characteristics (i.e., largest correlation = 0.097 with failing the attention check). Almost 49% of subjects failed our attention check at the end of the experiment.

99.6% of those subjects failed by answering “United States” which suggest they are inattentive humans (versus computerized agents who would probably randomly select one of the 5 answer options).

Most subjects were moderately risk averse, but were willing to increase their exposure to risk by almost one third (corresponding to about 12 boxes) when prompted to do so in our experimental demand task. Most subjects also made mistakes on the Bayesian task, deviating from the correct answer by 13.4 percentage points (which is 22.3% of the correct answer, 60) on average. In contrast, less than a third of the subjects (around 32%) behaved altruistically, i.e., chose to give money to others in their group rather than keep it for themselves.

Section 4.2: Text below Figure 3

(page 21) Can be found at the code for Figure 3.

Section 4.2: Text below Table 5 (page 23)

(page 23) Can be found at the code for Table 5.

Section 4.2: Misunderstanding and altruism

(page 24)

## MISUNDERSTANDING
# scenario 2
table <- table(subset(df, bot == 0 & treatment == "Partial Info")$instr_fail, subset(df, bot == 0 & treatment == "Partial Info")$decision2)
misunderstanders_optimal2 <- round(prop.table(table,1)["1", "A"]*100,1)
understanders_optimal2 <- round(prop.table(table,1)["0", "A"]*100,1)
fisher_OR <- round(fisher.test(table)$estimate,2)
fisher_p <- round(fisher.test(table)$p.value,3)
# scenario 1
table1 <- table(subset(df, bot == 0 & treatment == "Partial Info")$instr_fail, subset(df, bot == 0 & treatment == "Partial Info")$decision1)
misunderstanders_optimal1 <- round(prop.table(table1,1)["1", "A"]*100,1)
understanders_optimal1 <- round(prop.table(table1,1)["0", "A"]*100,1)
fisher_OR1 <- round(fisher.test(table1)$estimate,2)
fisher_p1 <- round(fisher.test(table1)$p.value,3)

## ALTRUISM
# scenario 2
table_altruists <- table(subset(df, bot == 0 & treatment == "Partial Info")$redistribution, subset(df, bot == 0 & treatment == "Partial Info")$decision2)
altruists_optimal2 <- round(prop.table(table_altruists,1)["1", "A"]*100,1) 
non_altruists_optimal2 <- round(prop.table(table_altruists,1)["0", "A"]*100,1) 
fisher_OR_altr <- round(fisher.test(table_altruists)$estimate,2)
fisher_p_altr <- round(fisher.test(table_altruists)$p.value,3) 
# scenario 1
table_altruists1 <- table(subset(df, bot == 0 & treatment == "Partial Info")$redistribution, subset(df, bot == 0 & treatment == "Partial Info")$follow_hint1)
altruists_optimal1 <- round(prop.table(table_altruists1,1)["1", "1"]*100,1) 
non_altruists_optimal1 <- round(prop.table(table_altruists1,1)["0", "1"]*100,1)
fisher_OR_altr1 <- round(fisher.test(table_altruists1)$estimate,2)
fisher_p_altr1 <- round(fisher.test(table_altruists1)$p.value,3)

The raw data confirm the pattern that those who misunderstood the instructions are less likely than those who understood (76.9% versus 84.8%) to make the optimal decision in Scenario 2 (Fisher’s Exact Test, OR = 1.67, p = 0.005). There is no difference between those who misunderstood (57.0%) and those who did not (53.5%) in optimal decisions in Scenario 1 (Fisher’s Exact Test, OR = 0.87 , p = 0.328).

However, we do not find support for altruism driving optimal behavior in either scenario (Fisher’s Exact Test, \(OR_{scenario1}\) = 1.03, \(p_{scenario1}\) = 0.935, \(OR_{scenario2}\) = 1.13, \(p_{scenario2}\) = 0.487). For this reason we caution against placing too much emphasis on the significant coefficient of altruism in Table 6.

Section 4.3: Self-reported strategies

(page 27)

strategies <- df %>% 
  filter(bot == 0) %>% 
  filter(treatment == "Partial Info") %>% 
  select(starts_with("strategy_"), ends_with("report"))  %>% 
  summarize_all(~mean(. == 1, na.rm = TRUE)) %>% 
  t() %>% 
  as.data.frame() %>% 
  rownames_to_column() %>% 
  rename(percent = V1,
         strategy = rowname) 

strategies10 <- strategies %>% 
  filter(strategy %in% c("risk_report", "strategy_misunderstanding", "strategy_altruist", "expdemand_report", "other_report", "noexp_report")) %>% 
  mutate(label = case_when(strategy == "risk_report" ~ "Risk preferences",
                           strategy == "strategy_misunderstanding" ~ "Mistake",
                           strategy == "strategy_altruist" ~ "Altruist",
                           strategy == "expdemand_report" ~ "Experimenter demand",
                           strategy == "other_report" ~ "Other",
                           strategy == "noexp_report" ~ "No explanation"))

strategies11 <- strategies %>% 
  filter(strategy %in% c("safe_report", "risky_report", "tradeoff_report")) %>% 
  mutate(percent = paste0(round(percent*100,1), "%")) 

pct_tradeoff <- strategies11 %>% filter(strategy == "tradeoff_report") %>% select(percent)
pct_safe <- strategies11 %>% filter(strategy == "safe_report") %>% select(percent)
pct_risky <- strategies11 %>% filter(strategy == "risky_report") %>% select(percent)

strategies12 <- df %>% 
  filter(bot == 0) %>% 
  filter(treatment == "Partial Info") %>% 
  select("strategy_ignore", "strategy_follow", "strategy_consider")  %>% 
  rename("Ignore" = "strategy_ignore",
         "Follow" = "strategy_follow",
         "Consider" = "strategy_consider") %>% 
  summarize_all(~mean(. == 1, na.rm = TRUE)) %>% 
  t() %>% 
  as.data.frame() %>% 
  rownames_to_column() %>% 
  rename(percent = V1,
         strategy = rowname) %>% 
  mutate(percent = paste0(round(percent*100,1), "%")) 

pct_followed <- strategies12 %>% filter(strategy == "Follow") %>% select(percent)
pct_consider <- strategies12 %>% filter(strategy == "Consider") %>% select(percent)
pct_ignore <- strategies12 %>% filter(strategy == "Ignore") %>% select(percent)

pct_no_other <- df %>% 
  filter(treatment == "Partial Info" & strategy_follow == 1) %>% 
  select("strategy_ignore", "strategy_consider", "strategy_misunderstanding", "strategy_altruist") %>% 
  mutate(summed_strategies = as.numeric(strategy_ignore) + as.numeric(strategy_consider) +
           as.numeric(strategy_misunderstanding) + as.numeric(strategy_altruist)) 
  
pct_no_other <- round(prop.table(table(pct_no_other$summed_strategies))[1],3)*100




#########

strategies <- df %>% 
  filter(bot == 0) %>% 
  filter(treatment == "Partial Info") %>% 
  select(starts_with("strategy_"), ends_with("report"))  %>% 
  summarize_all(~mean(. == 1, na.rm = TRUE)) %>% 
  t() %>% 
  as.data.frame() %>% 
  rownames_to_column() %>% 
  rename(percent = V1,
         strategy = rowname) %>% 
  mutate(percent = paste0(round(percent*100,1), "%"))  

pct_riskprefs <- strategies %>% filter(strategy == "risk_report") %>% select(percent)
pct_tradeoff <- strategies %>% filter(strategy == "tradeoff_report") %>% select(percent)
pct_other <- strategies %>% filter(strategy == "other_report") %>% select(percent)
pct_safe <- strategies %>% filter(strategy == "safe_report") %>% select(percent)
pct_risky <- strategies %>% filter(strategy == "risky_report") %>% select(percent)
pct_altruist <- strategies %>% filter(strategy == "strategy_altruist") %>% select(percent)
pct_expdemand <- strategies %>% filter(strategy == "expdemand_report") %>% select(percent)
pct_noexp <- strategies %>% filter(strategy == "noexp_report") %>% select(percent)
pct_misunderstanding <- strategies %>% filter(strategy == "strategy_misunderstanding") %>% select(percent)

pct_follow <- strategies %>% filter(strategy == "strategy_follow") %>% select(percent)
pct_ignore <- strategies %>% filter(strategy == "strategy_ignore") %>% select(percent)
pct_consider <- strategies %>% filter(strategy == "strategy_consider") %>% select(percent)

p_bayesian <- round(as.numeric(wilcox.test(bayesian ~ risk_report, data = subset(df, bot == 0 & treatment == "Partial Info"))[3]),3)

majority <- df %>% 
  filter(treatment == "Partial Info" & bot == 0 & strategy_follow == 1) %>% 
  select(agreed) %>% 
  group_by(agreed) %>% 
  tally() %>% 
  mutate(pct = n / sum(n)*100) %>% 
  filter(agreed == "Follow hint") %>% 
  select(pct) 

majority <- round(majority,1)

## check experimenter demand in more detail (how many 'own decision' comments?)

exp_demand_check <- df %>% 
  filter(bot == 0) %>% 
  filter(treatment == "Partial Info") %>% 
  select(id_otree, expdemand_report, strategy) %>% 
  filter(expdemand_report == 1)

## all three are negative experimenter demand: 
# I didn't really care about the recommendations. I made my own decisions based on what I wanted.
# My own decision
# I didn't use the recommendation, I just picked the one I wanted.

Overall, 62.5% of the subjects explained their decision depended on their risk preferences; of these, the majority reported they weighed the risks and benefits of the two options (57.7%), while about a third (37.5%) indicated that they simply chose the safe option. Only a small share of subjects indicated they simply chose the risky option (4.8%). Recall from Section 2.4 our theoretical prediction that risk preferences matter only for (sufficiently) prior-biased agents; however, we find no support in the data that these subjects who considered their risk preferences are any better at the Bayesian updating task than the other subjects (Wilcoxon rank sum test, \(p\) = 0.754). We therefore point to these self-reports as symptomatic of a particular type of mistake: either not updating correctly, or not even realizing that one should update following the recommendation.

The second most common category of reasons given was “other” (10.9%); this included explanations such as “choosing the best option” without explaining what “best” meant, or relying on “gut feeling”. Other reasons provided were relatively uncommon, such as appeals to altruistic motives (0.7%) or experimenter demand (0.3%). 1.8% of the subjects failed to provide an explanation for their choice.

Finally, only (0.9%) of subjects provided an explanation that could be directly identified as a reasoning mistake or confusion. Unsurprisingly, all of these subjects also failed at least one question in our comprehension test, with Q4 being the most common mistake. However, their responses do not provide us enough detail to say anything further about what specifically they do not understand.

On top of these categories, the research assistants also indicated whether the subject mentioned they followed their recommendation (26.6%), took it into consideration (21.3%), or ignored it (49.7%). In cases the subjects did not mention the recommendation at all, the research assistants coded such responses as ignoring the recommendation.

Of those subjects who mentioned they followed the recommendation they were given, the majority (82.3) did not mention any other reason for their choice.

Figures

Figure 2: Treatment effects

df_summary <- df %>% 
  group_by(treatment) %>% 
  summarize(mean_payoff1 = mean(group_payoff1), 
            sd_payoff1 = sd(group_payoff1), 
            n = n()) %>% 
  mutate(high = mean_payoff1 + qt(0.975, n-1)*(sd_payoff1/sqrt(n)), 
         low = mean_payoff1 - qt(0.975, n-1)*(sd_payoff1/sqrt(n)))

left <- ggplot(df_summary, aes(x = treatment, y = mean_payoff1, fill = treatment)) +
  geom_col(position = position_dodge()) +
  geom_errorbar(aes(ymin = low, ymax = high), 
                width = 0.2, position = position_dodge(0.9)) +
  geom_hline(yintercept=1348.48, linetype = 2, size = 1) +
  xlab("Treatment") + ylab("Mean Payoff") +
  scale_fill_manual(values = c("#A0A0A0", "#D9D9D9", "#505050"), guide = FALSE) +
  theme_bw() +
  scale_x_discrete(limits = c("No Info", "Partial Info", "Full Info")) +
  scale_y_continuous(breaks = seq(1000, 1350, by = 50)) +
  coord_cartesian(ylim = c(1000, 1380)) + labs(y = "", x = "\n (a) Scenario 1", title = "Aggregate Social Welfare")

df_summary <- df %>% 
  group_by(treatment) %>% 
  summarize(mean_payoff2 = mean(group_payoff2), 
            sd_payoff2 = sd(group_payoff2), 
            n = n()) %>% 
  mutate(high = mean_payoff2 + qt(0.975, n-1)*(sd_payoff2/sqrt(n)), 
         low = mean_payoff2 - qt(0.975, n-1)*(sd_payoff2/sqrt(n)))

right <- ggplot(df_summary, aes(x = treatment, y = mean_payoff2, fill = treatment)) +
  geom_col(position = position_dodge()) +
  geom_errorbar(aes(ymin = low, ymax = high), 
                width = 0.2, position = position_dodge(0.9)) +
  geom_hline(yintercept = 1148.48, linetype = 2, size = 1) +
  xlab("Treatment") + ylab("Mean Payoff") +
  scale_fill_manual(values = c("#A0A0A0", "#D9D9D9", "#505050"), guide = FALSE) +
  theme_bw() +
  scale_x_discrete(limits = c("No Info", "Partial Info", "Full Info")) +
  scale_y_continuous(breaks = seq(800, 1150, by = 50)) +
  coord_cartesian(ylim = c(800, 1180)) + labs(y = "", x = "\n (b) Scenario 2", title = "Aggregate Social Welfare")


g <- grid.arrange(left, right, nrow=1)
The figure plots the aggregate social welfare reached in Scenarios 1 and 2 with 95% confidence intervals. The dashed line indicates the theoretical maximum aggregate welfare (first best) that can be achieved in expectation.

The figure plots the aggregate social welfare reached in Scenarios 1 and 2 with 95% confidence intervals. The dashed line indicates the theoretical maximum aggregate welfare (first best) that can be achieved in expectation.

ggsave(here("output", "figures", "treatment_effects.pdf"), g, device = 'pdf',
       width = 20, height = 10, units = "cm")
rm(df_summary, left, right)

Figure 3: Reshufflings

df_groups <- df %>%
  group_by(group_id, treatment) %>%
  summarise(across(starts_with("payoff1"), mean)) %>% 
  mutate_at(vars(matches("^payoff1")), funs(4 * .)) %>% 
  group_by(treatment) %>% 
  summarise(across(starts_with("payoff1"), mean))  %>% 
  rename(payoff1_alt1 = payoff1) %>% 
  pivot_longer(cols = starts_with("payoff1_alt"), names_to = "time", 
               names_prefix = "payoff1_alt", values_to = "payoff1_alt") %>% 
  mutate(time = as.numeric(time)) %>% 
  mutate(treatment = ordered(treatment, levels = c("No Info", "Partial Info", "Full Info")))

left2 <-  ggplot(df_groups, aes(x = time, y = payoff1_alt)) +
  geom_hline(yintercept=1195, linetype = 2, size = 0.5) +
  geom_hline(yintercept=1262, linetype = 2, size = 0.5) +
  geom_hline(yintercept=1278, linetype = 2, size = 0.5) +
  geom_bar(aes(fill = treatment), position = "dodge", stat = "identity") + 
  facet_grid(~treatment, scales = "free", switch = "both") +
  xlab("Treatment") + ylab("Mean Payoff") +
  scale_fill_manual(values = c("#D9D9D9", "#505050", "#A0A0A0")) +
  theme_minimal() +
  scale_y_continuous(breaks = seq(1000, 1400, by = 100)) +
  coord_cartesian(ylim = c(1000, 1400)) + labs(y = "", x = "\n (a) Scenario 1", title = "") +
  theme(axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        legend.position = "none")

df_groups <- df %>%
  group_by(group_id, treatment) %>%
  summarise(across(starts_with("payoff2"), mean)) %>% 
  mutate_at(vars(matches("^payoff2")), funs(4 * .)) %>% 
  group_by(treatment) %>% 
  summarise(across(starts_with("payoff2"), mean))  %>% 
  rename(payoff2_alt1 = payoff2) %>% 
  pivot_longer(cols = starts_with("payoff2_alt"), names_to = "time", 
               names_prefix = "payoff2_alt", values_to = "payoff2_alt") %>% 
  mutate(time = as.numeric(time)) %>% 
  mutate(treatment = ordered(treatment, levels = c("No Info", "Partial Info", "Full Info")))

right2 <-  ggplot(df_groups, aes(x = time, y = payoff2_alt)) +
  geom_hline(yintercept=1077, linetype = 2, size = 0.5) +
  geom_hline(yintercept=1031, linetype = 2, size = 0.5) +
  geom_hline(yintercept=1010, linetype = 2, size = 0.5) +
  geom_bar(aes(fill = treatment), position = "dodge", stat = "identity") + 
  facet_grid(~treatment, scales = "free", switch = "both") +
  xlab("Treatment") + ylab("Mean Payoff") +
  scale_fill_manual(values = c("#D9D9D9", "#505050", "#A0A0A0")) +
  theme_minimal() +
  scale_y_continuous(breaks = seq(800, 1200, by = 100)) +
  coord_cartesian(ylim = c(800, 1200)) + labs(y = "", x = "\n (b) Scenario 2", title = "") +
  theme(axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        legend.position = "none")


g <- grid.arrange(left2, right2, nrow=1)
The figure plots the aggregate social welfare reached in Scenarios 1 and 2 under different subject reshufflings within groups. The dashed lines indicate the average social welfare levels achieved in our default ordering of subjects.

The figure plots the aggregate social welfare reached in Scenarios 1 and 2 under different subject reshufflings within groups. The dashed lines indicate the average social welfare levels achieved in our default ordering of subjects.

ggsave(here("output", "figures", "reshufflings.pdf"), g, device = 'pdf',
       width = 20, height = 10, units = "cm")
rm(left2, right2, df_groups)
tests_scenario2 <- df %>%
  filter(treatment != "Full Info") %>%
  select(starts_with("payoff2_alt"), treatment) %>%
  summarize(across(starts_with("payoff2_alt"), ~wilcox.test(.x ~ treatment)$p.value)) %>% 
  t() %>% 
  as.data.frame() %>% 
  rename(p_value = V1) %>% 
  mutate(significant = case_when(p_value < 0.05 ~ 1, TRUE ~ 0))

num_tests <- tests_scenario2 %>% tally() 
num_sig <- tests_scenario2 %>% filter(significant == 1) %>% tally()
pct_sig <- round(num_sig/(num_tests) *100,1)

Our results are generally stable with one exception: in the majority of alternative subject reshufflings under No Info in Scenario 2 the social welfare is lower than in our original dataset. In the 23 additional pairwise (Mann-Whitney-U) tests where we compare the aggregate social welfare in Scenario 2 under partial vs. no information, we find a significant difference in 4 reshufflings (17.4%). This is more than would be expected by pure chance (5%), and so we see it as a qualification of our result from Table 2, suggesting that there might be groups that slightly benefit from partial information as compared to no information even in Scenario 2.

Tables

Table 1: Summary Statistics

### group level
# aggregate welfare
group_summary <- maingroups %>% 
  group_by(treatment) %>% 
  summarise(mean1 = mean(group_payoff1),
         sd1 = sd(group_payoff1),
         mean2 = mean(group_payoff2),
         sd2 = sd(group_payoff2)) %>% 
  mutate_if(is.numeric, round, 0) %>% 
  mutate(means = paste0(mean1, " || ", mean2),
         sds = paste0(sd1, " || ", sd2)) %>% 
  select(means, sds, treatment) %>% 
  pivot_longer(cols = -treatment, names_to = "Treatment") %>% 
  pivot_wider(names_from = treatment, values_from = value)  %>% 
  mutate(rowname = c("Aggregate welfare", "SD")) %>% 
  select(-Treatment)
group_summary <- group_summary[c("rowname", "Full Info", "Partial Info", "No Info")]

# deviation of welfare
group_summary_dev <- maingroups %>% 
  mutate(payoff_perc_max1=100*group_payoff1/maximum1,
         payoff_perc_max2=100*group_payoff2/maximum2) %>% 
  group_by(treatment) %>% 
  summarise(mean1 = mean(payoff_perc_max1),
       sd1 = sd(payoff_perc_max1),
       mean2 = mean(payoff_perc_max2),
       sd2 = sd (payoff_perc_max2)) %>% 
  mutate_if(is.numeric, round, 2) %>% 
  mutate(means = paste0(mean1, " || ", mean2),
         sds = paste0(sd1, " || ", sd2)) %>% 
  select(means, sds, treatment) %>% 
  pivot_longer(cols = -treatment, names_to = "Treatment") %>% 
  pivot_wider(names_from = treatment, values_from = value)  %>% 
  mutate(rowname = c("... as % of max", "SD")) %>% 
  select(-Treatment) 
group_summary_dev <- group_summary_dev[c("rowname", "Full Info", "Partial Info", "No Info")]

# choosing A
chooseA1 <- df %>% 
  filter(bot == 0) %>% 
  group_by(treatment) %>% 
  count(choseA1) %>%
  summarise(choseA1 = choseA1,
            propA1 = round(prop.table(n)*100,2)) %>% 
  filter(choseA1 == 1)
chooseA2 <- df %>% 
  filter(bot == 0) %>% 
  group_by(treatment) %>% 
  count(choseA2) %>%
  summarise(choseA2 = choseA2,
            propA2 = round(prop.table(n)*100,2)) %>% 
  filter(choseA2 == 1)
choseA <- cbind(chooseA1, chooseA2) %>% 
  select(propA1, propA2, treatment...1) %>% 
  rename(treatment = treatment...1) %>%
  mutate("% choosing A" = paste0(as.character(propA1), " || ", as.character(propA2))) %>%
  select("% choosing A", treatment) %>% 
  pivot_longer(cols = -treatment, names_to = "Treatment") %>% 
  pivot_wider(names_from = treatment, values_from = value) %>% 
  rename(rowname = Treatment)
rm(chooseA2, chooseA1)
choseA <- choseA[c("rowname", "Full Info", "Partial Info", "No Info")]

# following Recommendation
follow1 <- round(prop.table(table((subset(df, bot == 0)$follow_hint1)))["1"] * 100, 2)
follow2 <- round(prop.table(table((subset(df, bot == 0)$follow_hint2)))["1"] * 100, 2)
follow <- choseA %>% 
  mutate("Full Info" = "", 
         "No Info" = "", 
         "Partial Info" = paste0(follow1, " || ", follow2),
         rowname = "% following recommendation")
rm(follow1, follow2)

table2 <- rbind(group_summary, group_summary_dev, choseA, follow)

### individual level
vars <- c("clean_time", "riskaverse", "exp_demand", "bayesian_abs", "female", "work_student", "instr_check", "instr_fail", "attention_fail", "redistribution")
table_out <- CreateTableOne(
  data = subset(df, bot == 0), 
  strata = "treatment", 
  vars = vars,
  factorVars = c("female", "work_student", "instr_check", "instr_fail", "attention_fail", "redistribution"),
  test = FALSE)

table_out_stats <- print(table_out, varLabels = TRUE, format = "p", catDigits = 2, printToggle = FALSE)
table_out_stats <- gsub(" \\(.*\\)", "", table_out_stats)
table_out_stats <- as.data.frame(table_out_stats)[2:11,]

table_out_stats <- cbind(table_out_stats[,1], table_out_stats[,3], table_out_stats[,2])
rownames(table_out_stats) <- c("Completion time (min)", "Risk aversion", "Experimenter demand", "Bayesian deviation",
                               "% female", "% student status", "% instructions check", "% instructions failure", "% attention failure", "% altruist")
colnames(table_out_stats) <- c("Full Info", "Partial Info", "No Info")

table2part2 <- as.data.frame(table_out_stats) %>% 
  rownames_to_column()

### merge tables
table2 <- rbind(table2, table2part2)

kable(table2, align = c("lccc")) %>% 
  add_header_above(c(" " = 1, "Scn 1 || Scn 2" = 1, "Scn 1 || Scn 2" = 1, "Scn 1 || Scn 2" = 1)) %>% 
  add_header_above(c(" " = 1, "Full Info" = 1, "Partial Info" = 1, "No Info" = 1)) %>% 
  pack_rows("Main outcomes:", start_row = 1, end_row = 6) %>% 
  pack_rows("Subject characteristics:", start_row = 7, end_row = 16) %>% 
  footnote(general ="Main outcomes of interest are reported for every within- and between-subject treatment (Scenarios 1 and 2, and Full/Partial/No Info). Aggregate welfare refers to group outcomes, whereas shares of subjects  choosing a specific strategy are calculated on an individual level. Individual-level statistics do not include choices determined by the computer for subjects who dropped out. Subject characteristics are reported for every between-subject treatment on an individual level. Completion time refers to the total number of minutes a subject took to complete the experiment, and is reported only for subjects who completed the experiment. Share of instructions check, instructions failure, and attention check refers to the share of participants who went back to re-read the instructions, who failed the instructions comprehension quiz, or who failed the attention check, respectively. Risk aversion refers to (100 - the number of boxes collected on the BRET task), such that higher values mean higher risk aversion. Experimenter demand is measured as the total increase in the number of boxes collected on the BRET task with explicit experimenter instructions as compared to our first task without any nudges about 'right' behavior. The deviation from Bayesian updating on our cookie task is listed as the absolute difference in probability reported and the Bayesian estimate. Finally, we report the percentage of subjects who choose to sacrifice their own payoff and instead benefit others in their group (altruist).")
Full Info
Partial Info
No Info
Scn 1 || Scn 2
Scn 1 || Scn 2
Scn 1 || Scn 2
rowname Full Info Partial Info No Info
Main outcomes:
Aggregate welfare 1278 &#124;&#124; 1077 1262 &#124;&#124; 1031 1195 &#124;&#124; 1010
SD 182 &#124;&#124; 186 168 &#124;&#124; 184 183 &#124;&#124; 168
… as % of max 93.92 &#124;&#124; 92.47 93.15 &#124;&#124; 89.02 89.7 &#124;&#124; 87.89
SD 8.15 &#124;&#124; 10.65 9.12 &#124;&#124; 12.88 11.09 &#124;&#124; 11.43
% choosing A 59.68 &#124;&#124; 76.24 55.78 &#124;&#124; 78.79 66.42 &#124;&#124; 86.35
% following recommendation 74.13 &#124;&#124; 67.02
Subject characteristics:
Completion time (min) 13.26 14.20 14.14
Risk aversion 58.14 59.22 60.70
Experimenter demand 11.73 11.65 13.52
Bayesian deviation 13.00 13.50 14.82
% female 42.28 43.12 39.61
% student status 6.13 9.10 7.06
% instructions check 10.41 8.13 10.20
% instructions failure 60.57 65.76 63.92
% attention failure 47.75 49.73 50.59
% altruist 28.06 29.90 28.78
Note:
Main outcomes of interest are reported for every within- and between-subject treatment (Scenarios 1 and 2, and Full/Partial/No Info). Aggregate welfare refers to group outcomes, whereas shares of subjects choosing a specific strategy are calculated on an individual level. Individual-level statistics do not include choices determined by the computer for subjects who dropped out. Subject characteristics are reported for every between-subject treatment on an individual level. Completion time refers to the total number of minutes a subject took to complete the experiment, and is reported only for subjects who completed the experiment. Share of instructions check, instructions failure, and attention check refers to the share of participants who went back to re-read the instructions, who failed the instructions comprehension quiz, or who failed the attention check, respectively. Risk aversion refers to (100 - the number of boxes collected on the BRET task), such that higher values mean higher risk aversion. Experimenter demand is measured as the total increase in the number of boxes collected on the BRET task with explicit experimenter instructions as compared to our first task without any nudges about ‘right’ behavior. The deviation from Bayesian updating on our cookie task is listed as the absolute difference in probability reported and the Bayesian estimate. Finally, we report the percentage of subjects who choose to sacrifice their own payoff and instead benefit others in their group (altruist).

Table 2: Treatment Effects on Social Welfare

library(DescTools)
maingroups <- maingroups %>% 
  mutate(treatment_trend = ordered(treatment, levels = c("Partial Info", "Full Info", "No Info")),
         treatment_trend2 = ordered(treatment, levels = c("Full Info", "Partial Info", "No Info")))

set.seed(1234)
jt1 <- JonckheereTerpstraTest(maingroups$group_payoff1, maingroups$treatment_trend, 
                       alternative = "decreasing", nperm = 100000, exact = TRUE)$p.value

set.seed(1234)
jt2 <- JonckheereTerpstraTest(maingroups$group_payoff2, maingroups$treatment_trend2, 
                       alternative = "decreasing", nperm = 100000, exact = TRUE)$p.value

mww1.1 <- wilcox.test(group_payoff1 ~ treatment, data = subset(maingroups, treatment %in% c("Full Info", "Partial Info")))$p.value
mww1.2 <- wilcox.test(group_payoff1 ~ treatment, data = subset(maingroups, treatment %in% c("Partial Info", "No Info")))$p.value
mww1.3 <- wilcox.test(group_payoff1 ~ treatment, data = subset(maingroups, treatment %in% c("Full Info", "No Info")))$p.value
mww2.1 <- wilcox.test(group_payoff2 ~ treatment, data = subset(maingroups, treatment %in% c("Full Info", "Partial Info")))$p.value
mww2.2 <- wilcox.test(group_payoff2 ~ treatment, data = subset(maingroups, treatment %in% c("Partial Info", "No Info")))$p.value
mww2.3 <- wilcox.test(group_payoff2 ~ treatment, data = subset(maingroups, treatment %in% c("Full Info", "No Info")))$p.value


p_values <- c(mww1.1, mww1.2, mww1.3, mww2.1, mww2.2, mww2.3)

# manually recoded this https://github.com/BITSS/IDBMarch2018/blob/master/4-MultipleTesting/fdr_sharpened_qvalues.do to an R function
fdr_sharpened_qvalues <- function(p_values) {
  total_pvals <- length(p_values)
  
  # Sort the p-values in ascending order and generate a variable that codes each p-value's rank
  sorted_pvals <- sort(p_values)
  ranks <- rank(p_values)
  
  # Initialize the variables
  q_val <- 1
  bky06_qval <- rep(1, total_pvals)
  
  while (q_val > 0) {
    # First Stage
    q_val_adj <- q_val / (1 + q_val)
    fdr_temp1 <- q_val_adj * ranks / total_pvals
    reject_temp1 <- ifelse(fdr_temp1 >= p_values, 1, 0)
    reject_rank1 <- reject_temp1 * ranks
    total_rejected1 <- max(reject_rank1, na.rm = TRUE)
    
    # Second Stage
    q_val_2st <- q_val_adj * (total_pvals / (total_pvals - total_rejected1))
    fdr_temp2 <- q_val_2st * ranks / total_pvals
    reject_temp2 <- ifelse(fdr_temp2 >= p_values, 1, 0)
    reject_rank2 <- reject_temp2 * ranks
    total_rejected2 <- max(reject_rank2, na.rm = TRUE)
    
    # Update the BKY (2006) sharpened q-values
    bky06_qval <- ifelse(ranks <= total_rejected2 & !is.na(ranks), q_val, bky06_qval)
    
    # Reduce q_val by 0.001 and repeat the loop
    q_val <- q_val - 0.001
  }
  
  # Return the BKY (2006) sharpened q-values
  return(bky06_qval)
}

bky_adjusted <- round(fdr_sharpened_qvalues(p_values),3)

# Custom function to add stars based on p-values
add_stars <- function(x) {
  rounded <- ifelse(round(as.numeric(x), 3) == 0, "0.000", round(x, 3))
  stars <- ifelse(x < 0, "", ifelse(x < 0.001, "***", ifelse(x < 0.01, "**", ifelse(x < 0.05, "*", ""))))
  paste0(rounded, stars)
}

table3 <- data.frame(jt = c(jt1, jt2),
                     mww1 = c(mww1.1, mww2.1),
                     mww2 = c(mww1.2, mww2.2),
                     mww3 = c(mww1.3, mww2.3)) %>% 
  mutate_all(add_stars)

table3top <- table3[1,]
table3bottom <- table3[2,]

# add sharpened q-values
# paste square brackets
brackets <- function(x){
  ifelse(x == "", paste0(x), paste0("[", x, "]"))
}

q1 <- as.data.frame(cbind("", t(bky_adjusted[1:3]))) %>% 
  mutate_all(brackets)
q2 <- as.data.frame(cbind("", t(bky_adjusted[4:6])))%>% 
  mutate_all(brackets)

names(q1) <- names(table3)
names(q2) <- names(table3)

table3 <- rbind(table3top, q1, table3bottom, q2, c(603, 523, 342, 341))
names(table3) <- NULL
rownames(table3) <- NULL

kable(table3, align=c("cccc")) %>% 
  add_header_above(header = c("All comparisons" = 1, "Full vs. Partial" = 1, "Partial vs. No" = 1, "Full vs. No" = 1)) %>% 
  add_header_above(header = c("(Jonckheere-Terpstra)" = 1, "(Mann-Whitney-U)" = 3)) %>% 
  pack_rows("H1: Partial > Full > No", start_row = 1, end_row = 3) %>% 
  pack_rows("H2: Full > Partial > No", start_row = 3, end_row = 5) %>% 
  pack_rows("N(groups)", start_row = 5, end_row = 5) %>% 
  footnote(general = "The first column lists p-values from the Jonckheere-Terpstra trend test for the ordered aggregate social welfare levels, and columns 2-4 list p-values for two-sided pairwise comparisons using the Mann-Whitney-U test. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the six pairwise tests (Anderson, 2008) are in brackets.")
(Jonckheere-Terpstra)
(Mann-Whitney-U)
All comparisons
Full vs. Partial
Partial vs. No
Full vs. No
H1: Partial > Full > No
0.082 0.218 0.002** 0.000***
[0.096] [0.005] [0.002]
H2: Full > Partial > No
0.000*** 0.005** 0.381 0.005**
[0.005] [0.146] [0.005]
N(groups)
603 523 342 341
Note:
The first column lists p-values from the Jonckheere-Terpstra trend test for the ordered aggregate social welfare levels, and columns 2-4 list p-values for two-sided pairwise comparisons using the Mann-Whitney-U test. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the six pairwise tests (Anderson, 2008) are in brackets.

Table 3: Reaching first best

# H1: W* = Partial

firstbest1 <- maingroups %>% 
  filter(treatment == "Partial Info") 
firstbest1_p <- round(wilcox.test(firstbest1$deviationwelfare1, mu = 0)$p.value,3)
firstbest1_n <- nrow(firstbest1)

# H2: W* = Full
firstbest2 <- maingroups %>% 
  filter(treatment == "Full Info") 
firstbest2_n <- nrow(firstbest2)
firstbest2_p <- round(wilcox.test(firstbest2$deviationwelfare2, mu = 0)$p.value,3)

# add stars
table4 <- data.frame(scenario1 = c(firstbest1_p, firstbest1_n),
                     scenario2 = c(firstbest2_p, firstbest2_n)) %>% 
  mutate_all(add_stars)

# add sharpened q-values
bky_adjusted <- round(fdr_sharpened_qvalues(c(firstbest1_p, firstbest2_p)),3)
q1 <- as.data.frame(cbind(t(bky_adjusted[1:2]))) %>% 
  mutate_all(brackets)

names(q1) <- names(table4)
table4top <- table4[1,]
table4bottom <- table4[2,]

table4 <- rbind(table4top, q1, table4bottom)
names(table4) <- NULL
rownames(table4) <- NULL

kable(table4, align=c("cccc")) %>% 
  add_header_above(header = c("H1: W* = Partial" = 1, "H2: W* > Full" = 1)) %>% 
  add_header_above(header = c("Scenario 1" = 1, "Scenario 2" = 1)) %>% 
  pack_rows("N(groups)", start_row = 3, end_row = 3) %>% 
  footnote(general = "The table lists the p-values for the one-sample two-sided Wilcoxon sign rank test, comparing the treatment that was
hypothesized to equal (Scenario 1) or fail to reach (Scenario 2) the aggregate social welfare optimum to the theoretical first best. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the two tests (Anderson, 2008) are in brackets.") 
Scenario 1
Scenario 2
H1: W* = Partial
H2: W* > Full
0.000*** 0.000***
[0.001] [0.001]
N(groups)
262 261
Note:
The table lists the p-values for the one-sample two-sided Wilcoxon sign rank test, comparing the treatment that was
hypothesized to equal (Scenario 1) or fail to reach (Scenario 2) the aggregate social welfare optimum to the theoretical first best. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the two tests (Anderson, 2008) are in brackets.

Table 4: Following the Recommendation

follow1 <- subset(df, bot == 0)$follow_hint1
follow1_p <- round(wilcox.test(follow1, mu = 1)$p.value,3)
follow1_pct <- round(prop.table(table(follow1))[2]*100,2)

follow1a <- subset(df, bot == 0 & hint1 == "A")$follow_hint1
follow1a_p <- round(wilcox.test(follow1a, mu = 1)$p.value,3)
follow1a_pct <- round(prop.table(table(follow1a))[2]*100,2)

follow1b <- subset(df, bot == 0 & hint1 == "B")$follow_hint1
follow1b_p <- round(wilcox.test(follow1b, mu = 1)$p.value,3)
follow1b_pct <- round(prop.table(table(follow1b))[2]*100,2)

follow2a <- subset(df, bot == 0 & hint2 == "A")$follow_hint2
follow2a_p <- round(wilcox.test(follow2a, mu = 1)$p.value,3)
follow2a_pct <- round(prop.table(table(follow2a))[2]*100,2)

follow2b <- subset(df, bot == 0 & hint2 == "B")$follow_hint2
follow2b_p <- round(wilcox.test(follow2b, mu = 1)$p.value,3)
follow2b_pct <- round(prop.table(table(follow2b))[2]*100,2)

follow2 <- subset(df, bot == 0)$follow_hint2
follow2_pct <- round(prop.table(table(follow2))[2]*100,2)

follow12_pct <- follow1_pct - follow2_pct
follow12_p <- round(wilcox.test(follow1, follow2)$p.value,3)
follow12a_pct <- follow1a_pct - follow2a_pct
follow12a_p <- round(wilcox.test(follow1a, follow2a)$p.value,3)
follow12b_pct <- follow1b_pct - follow2b_pct
follow12b_p <- round(wilcox.test(follow1b, follow2b)$p.value,3)

p_values <- c(follow1_p, follow1a_p, follow1b_p, follow2a_p, follow2b_p, follow12_p, follow12a_p, follow12b_p)
bky_adjusted <- format(fdr_sharpened_qvalues(p_values),digits = 3, nsmall =0)


table5 <- data.frame(scenario1 = c(follow1_pct, follow1_p, follow1a_pct, follow1a_p, follow1b_pct, follow1b_p),
                     scenario2 = c(NA_real_, NA_real_, follow2a_pct, follow2a_p, follow2b_pct, follow2b_p),
                     scenario12 = c(NA_real_, NA_real_, follow12a_pct, follow12a_p, follow12b_pct, follow12b_p)) %>% 
  mutate_all(add_stars) 

table5top <- table5[1:2,]
table5middle <- table5[3:4,]
table5bottom <- table5[5:6,]

q1 <- as.data.frame(cbind(bky_adjusted[1], NA_real_, NA_real_)) %>% 
  mutate_all(brackets) 
q1[2] <- " "
q1[3] <- " "
table5top[2,2] <-  " "
table5top[2,3] <-  " "
table5top[1,3] <-  " "
table5top[1,2] <-  " "

q2 <- as.data.frame(cbind(t(bky_adjusted[3:5]))) %>% 
  mutate_all(brackets)
q3 <- as.data.frame(cbind(t(bky_adjusted[6:8]))) %>% 
  mutate_all(brackets)

names(q1) <- names(table5)
names(q2) <- names(table5)
names(q3) <- names(table5)

table5 <- rbind(table5top, q1, 
                table5middle, q2,
                table5bottom, q3, 
                c(943, 943, 943))
names(table5) <- NULL
rownames(table5) <- NULL

table5[7,3] <- format(round(as.numeric(table5[7,3]),2), nsmall = 2)


kable(table5, align=c("ccc")) %>% 
  add_header_above(header = c("H3: always follow" = 1, "H3: only follow A" = 1, " " = 1)) %>% 
  add_header_above(header = c("Scenario 1" = 1, "Scenario 2" = 1, "Scenario 1 vs. 2" = 1)) %>% 
  pack_rows("% follow", start_row = 1, end_row = 3) %>% 
  pack_rows("% follow A", start_row = 4, end_row = 6) %>% 
  pack_rows("% follow B", start_row = 7, end_row = 9) %>% 
  pack_rows("N (subjects)", start_row = 10, end_row = 10) %>% 
  footnote(general = "The table lists the average share of subjects who follow the Recommendation in the Partial Info treatment and the associated p-values for comparisons w.r.t. theoretical benchmarks (columns 1 and 2) and between scenarios (column 3). We do not run the aggregate test for hint following (row 1) in scenario 2 since the theoretical prediction depends on the realized share of subjects receiving each recommendation; we therefore only test conditional on realized recommendations (rows 2 and 3). * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the seven tests (Anderson, 2008) are in square brackets.")
Scenario 1
Scenario 2
Scenario 1 vs. 2
H3: always follow
H3: only follow A
% follow
74.13
0.000***
[0.001]
% follow A
79.62 95.38 -15.76
0.000*** 0.000*** 0.000***
[0.001] [0.001] [0.001]
% follow B
68.52 38.12 30.40
0.000*** 0.000*** 0.000***
[0.001] [0.001] [0.001]
N (subjects)
943 943 943
Note:
The table lists the average share of subjects who follow the Recommendation in the Partial Info treatment and the associated p-values for comparisons w.r.t. theoretical benchmarks (columns 1 and 2) and between scenarios (column 3). We do not run the aggregate test for hint following (row 1) in scenario 2 since the theoretical prediction depends on the realized share of subjects receiving each recommendation; we therefore only test conditional on realized recommendations (rows 2 and 3). * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the seven tests (Anderson, 2008) are in square brackets.

Table 5: Strategies played

actionprofiles <- df %>% 
  filter(bot == 0) %>% 
  mutate(doubleaction = paste0(decision1, decision2),
         doublehint = paste0(hint1, hint2, " "),
         Ahigher1 = case_when(valuationsb1 == 100 ~ "A",
                              valuationsb1 == 300 ~ "=",
                              TRUE ~ "B"),
         Ahigher2 = case_when(valuationsb2 == 0 | valuationsb2 == 200 ~ "A",
                              TRUE ~ "B"),
         Ahigher = paste0(Ahigher1, Ahigher2),
         Ahigher = ordered(Ahigher, levels = c("AA", "AB", "BA", "BB", "=A", "=B")))

# No Info
part1 <- actionprofiles %>% 
  filter(treatment == "No Info") %>%
  select(doubleaction) %>% 
  table() %>% 
  as.data.frame() %>% 
  rename(rowname = "doubleaction", " " = Freq) %>% 
  column_to_rownames() %>% t()

averse_neutral <- round((part1[,"AA"]+ part1[,"BA"])/sum(part1) * 100,1)

# Partial Info: Recommendations Received
part2 <- actionprofiles %>% 
  filter(treatment == "Partial Info") %>% 
  select(doublehint, doubleaction) %>% 
  table() %>% 
  as.data.frame.matrix()

partial_optimal <- round(((sum(part2[1:2,1]) + sum(part2[3:4,3])) / sum(part2))*100,1)

# Full Info: Highest expected payoff
part3 <- actionprofiles %>% 
  filter(treatment == "Full Info") %>% 
  select(Ahigher, doubleaction) %>% 
  table() %>% 
  as.data.frame.matrix()

full_optimal <- round((sum(diag(as.matrix(part3[1:4,]))) + part3[5,1] + part3[5,3] + part3[6,2] + part3[6,4])/sum(part3)*100, 2)

n_subjects <- actionprofiles %>% 
  select(doubleaction) %>% 
  table() %>% 
  as.data.frame() %>% 
  rename(rowname = "doubleaction", "  " = Freq) %>% 
  column_to_rownames() %>% t()

table6 <- rbind(part1, part2, part3, n_subjects)
rm(part1, part2, n_subjects)

kable(table6, align = c("lcccc")) %>% 
  add_header_above(c(" " = 1, "Actions taken" = 4)) %>% 
  pack_rows("No Info:", start_row = 1, end_row = 1) %>% 
  pack_rows("Partial Info: Recommendations Received", start_row = 2, end_row = 5) %>% 
  pack_rows("Full Info: Highest expected payoff", start_row = 6, end_row = 11) %>% 
  pack_rows("N (subjects)", start_row = 12, end_row = 12) 
Actions taken
AA AB BA BB
No Info:
164 16 70 21
Partial Info: Recommendations Received
AA 225 1 43 9
AB 96 57 35 10
BA 64 2 122 10
BB 66 15 92 96
Full Info: Highest expected payoff
AA 113 16 6 6
AB 7 29 1 1
BA 15 5 128 8
BB 1 4 6 36
=A 266 18 119 18
=B 32 49 15 31
N (subjects)
1049 212 637 246

text below Table 5

As Table 5 shows, the majority of subjects behave in accordance with the theoretical predictions of Bayesian, self-interested profit maximizers. In the Full Info treatment, most (82.9%) subjects choose the individually payoff-maximizing object. In contrast, only 56.7% of subjects do so in the Partial Info treatment. In No Info, the individually optimal payoff-maximizing option depends on risk preferences: a full 86.3% of subjects play strategies consistent either with risk aversion or risk neutrality (AA or BA), which is largely consistent with past experimental research on risk preferences

Table 6: Optimal Choices across Scenarios: Partial Info

df.d2 <- df %>% 
  mutate(recommendedA = case_when(hint1 == "A"~ 1, TRUE ~ 0)) %>% 
  filter(treatment == "Partial Info" & bot == 0) %>% 
  filter(!is.na(riskaverse) & !is.na(bayesian_abs) & !is.na(attention_fail) & !is.na(instr_fail) & !is.na(redistribution) & !is.na(exp_demand))


############################### scenario 1

# Fit the logistic regression model using `glm`
model1a <- glm(follow_hint1 ~ recommendedA, data = df.d2, family = binomial())
model1 <- glm(follow_hint1 ~ recommendedA + riskaverse * bayesian_abs + attention_fail + instr_fail + redistribution + exp_demand, data = df.d2, family = binomial())

# Compute cluster-robust standard errors using `coeftest` and `vcovCL` functions
cluster_se <- sqrt(diag(vcovCL(model1, cluster = ~ id_otree, type = "HC1")))
cluster_sea <- sqrt(diag(vcovCL(model1a, cluster = ~ id_otree, type = "HC1")))
coef_tests1 <- coeftest(model1, vcov = sandwich)[,4]
coef_tests1a <- coeftest(model1a, vcov = sandwich)[,4]
# compute odds ratios by exponentiating
odds_ratios <- exp(coef(model1))
odds_ratiosa <- exp(coef(model1a))

############################### scenario 2

df.d2.2 <- df %>% 
  mutate(recommendedA = case_when(hint2 == "A"~ 1, TRUE ~ 0)) %>% 
  filter(treatment == "Partial Info" & bot == 0) %>% 
  filter(!is.na(riskaverse) & !is.na(bayesian_abs) & !is.na(attention_fail) & !is.na(instr_fail) & !is.na(redistribution) & !is.na(exp_demand))


model2a <- glm(choseA2 ~ recommendedA, data = df.d2.2, family = binomial())
model2 <- glm(choseA2 ~ recommendedA + riskaverse * bayesian_abs + attention_fail + instr_fail + redistribution + exp_demand, data = df.d2.2, family = binomial())

odds_ratios2 <- exp(coef(model2))
odds_ratios2a <- exp(coef(model2a))
cluster_se2 <- sqrt(diag(vcovCL(model2, cluster = ~ id_otree, type = "HC1")))
cluster_se2a <- sqrt(diag(vcovCL(model2a, cluster = ~ id_otree, type = "HC1")))
coef_tests2 <- coeftest(model2, vcov = sandwich)[,4]
coef_tests2a <- coeftest(model2a, vcov = sandwich)[,4]

# models 2 with comprehension questions

model2a_check <- glm(choseA2 ~ recommendedA +check1 + check2 + check3 + check4, data = df.d2.2, family = binomial())
model2_check <- glm(choseA2 ~ recommendedA + riskaverse * bayesian_abs + attention_fail +  + redistribution + exp_demand+check1 + check2 + check3 + check4, data = df.d2.2, family = binomial())

odds_ratios2_check <- exp(coef(model2_check))
odds_ratios2a_check <- exp(coef(model2a_check))
cluster_se2_check <- sqrt(diag(vcovCL(model2_check, cluster = ~ id_otree, type = "HC1")))
cluster_se2a_check <- sqrt(diag(vcovCL(model2a_check, cluster = ~ id_otree, type = "HC1")))
coef_tests2_check <- coeftest(model2_check, vcov = sandwich)[,4]
coef_tests2a_check <- coeftest(model2a_check, vcov = sandwich)[,4]

# in the next code chunk I'm running this, but not including it (otherwise it gets knitted twice)

# stargazer_output <- stargazer(model1a, model1, model2a, model2, model2a_check, model2_check,
# coef =list(odds_ratiosa, odds_ratios, odds_ratios2a, odds_ratios2, odds_ratios2a_check, odds_ratios2_check),
# p = list(coef_tests1a, coef_tests1, coef_tests2a, coef_tests2, coef_tests2a_check, coef_tests2_check),
# se = list(cluster_sea, cluster_se, cluster_se2a, cluster_se2, cluster_se2a_check, cluster_se2_check),
# type = 'text',
# omit.table.layout = 'sn',
# order = c(1:3,8,4,5,6,7), # to get the interaction term to the right place
# star.cutoffs = c(0.05, 0.01, 0.001),
# header = FALSE
#           )
### this is a bit of manual work to rework the stargazer output to a dataframe, such that we can add some rows at the bottom

# Create an empty dataframe to store variable names and coefficients
output_df <- data.frame(Variable = character(), Coefficient_1 = character(), Coefficient_2 = character(), 
                        Coefficient_3 = character(), Coefficient_4 =character(), 
                        Coefficient_5 = character(), Coefficient_6 =character(),
                        stringsAsFactors = FALSE)

# Extract variable names and coefficients using regular expressions
for (i in 1:length(stargazer_output)) {
  line <- stargazer_output[i]
  
  # Skip the top rows with "Dependent variable"
  if (grepl("Dependent variable:", line)) {
    next
  }
  
  # Check if the line matches the pattern for variable names and coefficients
  if (grepl("^[^=]+\\s+[^\\s]+\\s+[^\\s]+\\s+[^\\s]+", line)) {
    # Extract variable name and coefficients using regex
    line_parts <- trimws(strsplit(line, "\\s+")[[1]])
    var_name <- line_parts[1]
    coefficient_1 <- line_parts[2]
    coefficient_2 <- line_parts[3]
    coefficient_3 <- line_parts[4]
    coefficient_4 <- line_parts[5]
    coefficient_5 <- line_parts[6]
    coefficient_6 <- line_parts[7]
    
    
    # Append variable name and coefficients to the dataframe
    output_df <- rbind(output_df, data.frame(Variable = var_name, Coefficient_1 = coefficient_1, Coefficient_2 = coefficient_2, 
                                             Coefficient_3 = coefficient_3, Coefficient_4 = coefficient_4, 
                                             Coefficient_5 = coefficient_5, Coefficient_6 = coefficient_6,
                                             stringsAsFactors = FALSE))
  }
}

output_df <- output_df %>% 
  filter(!is.na(Coefficient_2)) 

output_df <- output_df[2:26,] # to remove the constant
output_df[4:17,7] <- output_df[4:17,4]
output_df[4:17,5] <- output_df[4:17,3]
output_df[4:17,3] <- output_df[4:17,2]
output_df[17:26,6] <- output_df[17:26,2]
output_df[17:26,7] <- output_df[17:26,3]

output_df[4:26,2] <- c(" ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", "", "")
output_df[4:26,4] <- c(" ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", "", "")
output_df[4:17,6] <- c(" ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ")
output_df[18:26,3] <- c("", "", "", "", "", "", "", "", "")
output_df[18:26,5] <- c("", "", "", "", "", "", "", "", "")
output_df[12:13,7] <- c("", "")

output_df <- output_df[1:25,]

names(output_df) <- c("var", "col1", "col2", "col3", "col4", "col5", "col6") 

output_df <- output_df %>% 
  mutate(var = case_when(var == "recommendedA" ~ "Recommended A", 
                         var == "riskaverse" ~ "Risk averse",
                         var == "riskaverse:bayesian_abs" ~ "Risk averse X non-Bayesian",
                         var == "bayesian_abs" ~ "Non Bayesian",
                         var == "attention_fail" ~ "Inattention",
                         var == "instr_fail" ~ "# attempts comprehension Qs",
                         var == "redistribution" ~ "Altruist",
                         var == "exp_demand" ~ "Experimenter demand",
                         var == "check1" ~ "Failed comprehension Q1",
                         var == "check2" ~ "Failed comprehension Q2",
                         var == "check3" ~ "Failed comprehension Q3",
                         var == "check4" ~ "Failed comprehension Q4",
                         TRUE ~ var))
names(output_df) <- c(" ", "col1", "col2", "col3", "col4", "col5", "col6") 

total <- data.frame(" " = "N (total)", 
                    col1 = 869,
                    col2 = 869,
                    col3 = 869,
                    col4 = 869,
                    col5 = 869, 
                    col6 = 869)
names(total) <- c(" ", "col1", "col2", "col3", "col4", "col5", "col6") 

tabled2 <- rbind(output_df, total) 
tabled2 <- tabled2[2:26,]
rownames(tabled2) <- NULL
names(tabled2) <- NULL

kable(tabled2, align = c("lcccccc")) %>%  # use "format = 'latex', booktabs = TRUE" if you want to output latex 
  add_header_above(c(" " = 1,  "Scenario 1" = 2, "Scenario 2" = 4)) %>% 
  pack_rows(" ", start_row = 25, end_row = 25) %>% 
  footnote(general = "The table shows the odds ratios for making the theoretically optimal choices in the Partial Info treatment. Optimal strategies follow the Recommendation in the first Scenario, and select A in the second Scenario. Only subjects who completed all ancillary tasks are included. Simple logit estimation. Robust standard errors are in parentheses. Clustering on individual level. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001")
Scenario 1
Scenario 2
Recommended A 1.745*** 1.697*** 14.391*** 16.065*** 15.564*** 16.725***
(0.158) (0.159) (0.267) (0.266) (0.276) (0.274)
Risk averse 1.005 1.015 1.014
(0.006) (0.008) (0.008)
Non Bayesian 1.025 1.043 1.042
(0.016) (0.027) (0.028)
Risk averse X non-Bayesian 1.000 0.999 0.999
(0.0003) (0.0004) (0.0004)
Inattention 0.949 0.904 0.941
(0.159) (0.194) (0.197)
# attempts comprehension Qs 1.328 0.526**
(0.163) (0.215)
Altruist 1.013 0.637* 0.629*
(0.173) (0.207) (0.209)
Experimenter demand 1.006 1.001 1.000
(0.005) (0.006) (0.005)
Failed comprehension Q1 1.034 1.080
(0.192) (0.159)
Failed comprehension Q2 0.904 0.885
(0.190) (0.195)
Failed comprehension Q3 0.655** 0.634**
(0.138) (0.141)
Failed comprehension Q4 0.632* 0.659*
(0.188) (0.195)
N (total) 869 869 869 869 869 869
Note:
The table shows the odds ratios for making the theoretically optimal choices in the Partial Info treatment. Optimal strategies follow the Recommendation in the first Scenario, and select A in the second Scenario. Only subjects who completed all ancillary tasks are included. Simple logit estimation. Robust standard errors are in parentheses. Clustering on individual level. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001

Table 7: Choices of subjects who misunderstand in Partial Info

#### how often is each question failed?
partial_info_only <- df %>% filter(treatment == "Partial Info")
failed_at_least_one <- round(prop.table(table(df$instr_fail))[2]*100,1)

failed1pct <- round((1-prop.table(table(partial_info_only$check1))[1])*100,1) # 12.2%
failed2pct <- round((1-prop.table(table(partial_info_only$check2))[1])*100,1) # 23.9%
failed3pct <- round((1-prop.table(table(partial_info_only$check3))[1])*100,1) # 36.9%
failed4pct <- round((1-prop.table(table(partial_info_only$check4))[1])*100,1) # 29.0%

#### those who fail each of the four, the descriptive split of sc1 and sc2 choices
failed1 <- partial_info_only %>% filter(check1 != 0) %>% select(decision1, decision2)
failed1sc1 <- round(prop.table(table(failed1$decision1))[1]*100,1) # 49.6% chose A in scenario 1
failed1sc2 <- round(prop.table(table(failed1$decision2))[1]*100,1) # 82.3% chose A in scenario 2

failed2 <- partial_info_only %>% filter(check2 != 0) %>% select(decision1, decision2)
failed2sc1 <- round(prop.table(table(failed2$decision1))[1]*100,1) # 59.3 chose A in scenario 1
failed2sc2 <- round(prop.table(table(failed2$decision2))[1]*100,1) # 78.7% chose A in scenario 2

failed3 <- partial_info_only %>% filter(check3 != 0) %>% select(decision1, decision2)
failed3sc1 <- round(prop.table(table(failed3$decision1))[1]*100,1) # 57.8% chose A in scenario 1
failed3sc2 <- round(prop.table(table(failed3$decision2))[1]*100,1) # 74.8% chose A in scenario 2

failed4 <- partial_info_only %>% filter(check4 != 0) %>% select(decision1, decision2)
failed4sc1 <- round(prop.table(table(failed4$decision1))[1]*100,1) # 56.0% chose A in scenario 1
failed4sc2 <- round(prop.table(table(failed4$decision2))[1]*100,1) # 74.6% chose A in scenario 2

# compared to 
sc1 <- round(prop.table(table(partial_info_only$decision1))[1]*100,1) # 56% overall in partial info
sc2 <- round(prop.table(table(partial_info_only$decision2))[1]*100,1) # 76.2% overall in partial info
total <- round(prop.table(table(df$instr_fail))[2]*100,1)
########## turn into figure

# stitch numbers into one table
# put it after table 6

newtable <- data.frame(" " = c("Q1", "Q2", "Q3", "Q4", "Overall"),
                       total = c(failed1pct, failed2pct, failed3pct, failed4pct, total),
                       scenario1 = c(failed1sc1, failed2sc1, failed3sc1, failed4sc1, sc1),
                       scenario2 = c(failed1sc2, failed2sc2, failed3sc2, failed4sc2, sc2)
)

colnames(newtable) <- c(" ", " ", "chose A in scenario 1", "chose A in scenario 2")

kable(newtable, type = "text") %>% 
      add_header_above(c(" " =1, "wrong at least once" = 1, "of those with a wrong answer..." =2 ))
wrong at least once
of those with a wrong answer…
chose A in scenario 1 chose A in scenario 2
Q1 12.2 49.6 82.3
Q2 23.9 59.3 78.7
Q3 36.9 57.8 74.8
Q4 29.0 56.0 74.6
Overall 63.3 56.0 76.2

Appendix (Robustness)

Figure A1: Simulations No Info

Set up bots with oTree command line bots. In oTree models.py (lines 78-96) set player type to get the required simulation (e.g. alwayssafe0 means 0% alwaysA and 100% bayesian). oTree code of this experiment can be found at here. Stored results of the oTree bots can be found here.

        if self.round_number == 1:
            for p in players:
                # SETTING TYPES FOR BOTS:
                types = []
                for g in self.get_groups():
                    if g.treatment == "Full Info":
                        types = ['highest'] * 95 + ['protest'] * 5
                    elif g.treatment == "Partial Info":
                        types = ['bayesian'] * 100 + ['over-complier'] * 0 + ['incorrect-complier'] * 0 + ['defier'] * 0 + ['alwaysA'] * 0 + ['alwaysB'] * 0
                    elif g.treatment == "OR + Full Info":
                        types = ['highest'] * 0 + ['over-complier'] * 100
                    elif g.treatment == "No Info":
                        types = ['highest_ev'] * 100 + ['safe'] * 0 + ['gamble'] * 0
                p.participant.vars["player_type"] = random.choice(types)

The tests.py file defines what the player types do.

class PlayerBot(Bot):

    def play_round(self):
        num = self.player.round_number

        if self.round_number < 3:
            player_type = self.participant.vars["player_type"]
            if self.group.treatment == "Full Info":
                if self.participant.vars['valuations'][num][0] > self.participant.vars['valuations'][num][1]:
                    if player_type == 'highest':
                        yield Submission(pages.Decision, dict(pref='A', opened=0, player_type=player_type),
                                         check_html=False)
                    else:
                        yield Submission(pages.Decision, dict(pref='B', opened=0, player_type=player_type),
                                         check_html=False)
                elif self.participant.vars['valuations'][num][0] == self.participant.vars['valuations'][num][1]:
                    yield Submission(pages.Decision, dict(pref=random.choice(['A', 'B']), opened=0,
                                                          player_type=player_type), check_html=False)
                else:
                    if player_type == 'highest':
                        yield Submission(pages.Decision, dict(pref='B', opened=0, player_type=player_type),
                                         check_html=False)
                    else:
                        yield Submission(pages.Decision, dict(pref='A', opened=0, player_type=player_type),
                                         check_html=False)
            elif self.group.treatment == "Partial Info":
                if player_type == 'over-complier':
                    yield Submission(pages.Decision, dict(pref=self.player.participant.vars["hint"][num], opened=0,
                                                          player_type=player_type), check_html=False)
                elif player_type == 'defier':
                    if self.player.participant.vars["hint"][num] == 'A':
                        pref = 'B'
                    else:
                        pref = 'A'
                    yield Submission(pages.Decision, dict(pref=pref, opened=0,
                                                          player_type=player_type), check_html=False)
                elif player_type == "alwaysA":
                    yield Submission(pages.Decision, dict(pref='A', opened=0, player_type=player_type),
                                     check_html=False)
                elif player_type == "alwaysB":
                    yield Submission(pages.Decision, dict(pref='B', opened=0, player_type=player_type),
                                     check_html=False)
                elif player_type == "bayesian" and num == 1:
                    yield Submission(pages.Decision, dict(pref=self.player.participant.vars["hint"][num], opened=0,
                                                          player_type=player_type), check_html=False)
                elif player_type == "bayesian" and num == 2:
                    yield Submission(pages.Decision, dict(pref='A', opened=0, player_type=player_type),
                                     check_html=False)
                elif player_type == "incorrect-complier" and num == 1:
                    yield Submission(pages.Decision, dict(pref='A', opened=0,
                                                          player_type=player_type), check_html=False)
                elif player_type == "incorrect-complier" and num == 2:
                    yield Submission(pages.Decision, dict(pref=self.player.participant.vars["hint"][num], opened=0,
                                                          player_type=player_type), check_html=False)
            elif self.group.treatment == "OR + Full Info":
                if player_type == 'over-complier':
                    yield Submission(pages.Decision, dict(pref=self.player.participant.vars["hint"][num], opened=0,
                                                          player_type=player_type), check_html=False)
                elif player_type == 'highest':
                    if self.participant.vars['valuations'][num][0] > self.participant.vars['valuations'][num][1]:
                        yield Submission(pages.Decision, dict(pref='A', opened=0, player_type=player_type),
                                         check_html=False)
                    elif self.participant.vars['valuations'][num][0] == self.participant.vars['valuations'][num][1]:
                        yield Submission(pages.Decision, dict(pref=random.choice(['A', 'B']), opened=0,
                                                              player_type=player_type), check_html=False)
                    else:
                        yield Submission(pages.Decision, dict(pref='B', opened=0, player_type=player_type),
                                         check_html=False)
            elif self.group.treatment == "No Info":
                if player_type == 'highest_ev' and num == 1:
                    yield Submission(pages.Decision, dict(pref=random.choice(['A', 'B']), opened=0,
                                                          player_type=player_type), check_html=False)
                elif player_type == 'highest_ev' and num == 2:
                    yield Submission(pages.Decision, dict(pref='A', opened=0, player_type=player_type),
                                     check_html=False)
                elif player_type == 'safe':
                    yield Submission(pages.Decision, dict(pref='A', opened=0, player_type=player_type),
                                     check_html=False)
                elif player_type == 'gamble':
                    yield Submission(pages.Decision, dict(pref='B', opened=0, player_type=player_type),
                                     check_html=False)

In oTree pages.py (lines 452-472) comment out all other pages except Decision (otherwise the bots will also come up with simulated answer to all other questions, like demographics, and the bots run very slow.)

page_sequence = [
    # Welcome,
    # GeneralInstructions,
    # Instructions,
    # PracticeDecision,
    # ResultsPractice,
    # Check,
    # Check2,
    # Check3,
    # Check4,
    Decision,
    # ExplainTimeout,
    # Strategy,
    # CookieQuestion,
    # InstructionsBRET,
    # BRET,
    # Demographics,
    # Redistribution,
    # Payment,
    # Thanks
]

Run in oTree terminal: otree test matching2 10000 --export=matching/cleaning_data_in_r/bots/alwayssafe0 where matching2 refers to the Partial Info treatment, 10000 refers to the number of bots, and the final line defines where the bots should be stored.

knitr::opts_chunk$set(echo = TRUE)

alwayssafe0 <- read.csv(here("bots", "alwayssafe0", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "alwayssafe", pct_bayesian = 100)
alwayssafe20 <- read.csv(here("bots", "alwayssafe20", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "alwayssafe", pct_bayesian = 80)
alwayssafe40 <- read.csv(here("bots", "alwayssafe40", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "alwayssafe", pct_bayesian = 60)
alwayssafe60 <- read.csv(here("bots", "alwayssafe60", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "alwayssafe", pct_bayesian = 40)
alwayssafe80 <- read.csv(here("bots", "alwayssafe80", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "alwayssafe", pct_bayesian = 20)
alwayssafe100 <- read.csv(here("bots", "alwayssafe100", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "alwayssafe", pct_bayesian = 0)

defiers0 <- read.csv(here("bots", "defiers0", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "defiers", pct_bayesian = 100)
defiers20 <- read.csv(here("bots", "defiers20", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "defiers", pct_bayesian = 80)
defiers40 <- read.csv(here("bots", "defiers40", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "defiers", pct_bayesian = 60)
defiers60 <- read.csv(here("bots", "defiers60", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "defiers", pct_bayesian = 40)
defiers80 <- read.csv(here("bots", "defiers80", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "defiers", pct_bayesian = 20)
defiers100 <- read.csv(here("bots", "defiers100", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "defiers", pct_bayesian = 0)

overcompliers0 <- read.csv(here("bots", "overcomplier0", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "overcompliers", pct_bayesian = 100)
overcompliers20 <- read.csv(here("bots", "overcomplier20", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "overcompliers", pct_bayesian = 80)
overcompliers40 <- read.csv(here("bots", "overcomplier40", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "overcompliers", pct_bayesian = 60)
overcompliers60 <- read.csv(here("bots", "overcomplier60", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "overcompliers", pct_bayesian = 40)
overcompliers80 <- read.csv(here("bots", "overcomplier80", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "overcompliers", pct_bayesian = 20)
overcompliers100 <- read.csv(here("bots", "overcomplier100", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "overcompliers", pct_bayesian = 0)

# in oTree terminal:
#  types = ['highest_ev'] * 0 + ['safe'] * 100
# otree test matching1 10000 --export=matching/cleaning_data_in_r/bots/no_info_random0

no_info_random0 <- read.csv(here("bots", "no_info_random0", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "noinfo random", pct_bayesian = 0)
no_info_random20 <- read.csv(here("bots", "no_info_random20", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "noinfo random", pct_bayesian = 20)
no_info_random40 <- read.csv(here("bots", "no_info_random40", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "noinfo random", pct_bayesian = 40)
no_info_random60 <- read.csv(here("bots", "no_info_random60", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "noinfo random", pct_bayesian = 60)
no_info_random80 <- read.csv(here("bots", "no_info_random80", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "noinfo random", pct_bayesian = 80)
no_info_random100 <- read.csv(here("bots", "no_info_random100", "all_apps_wide.csv"), sep = ',', header=TRUE) %>% mutate(simulation = "noinfo random", pct_bayesian = 100)

choices <- rbind(defiers0, defiers20, defiers40, defiers60, defiers80, defiers100,
                 overcompliers0, overcompliers20, overcompliers40, overcompliers60, overcompliers80, overcompliers100,
                 alwayssafe0, alwayssafe20, alwayssafe40, alwayssafe60, alwayssafe80, alwayssafe100,
                 no_info_random0, no_info_random20, no_info_random40, no_info_random60, no_info_random80, no_info_random100
                 )

rm(defiers0, defiers20, defiers40, defiers60, defiers80, defiers100,
   overcompliers0, overcompliers20, overcompliers40, overcompliers60, overcompliers80, overcompliers100,
   alwayssafe0, alwayssafe20, alwayssafe40, alwayssafe60, alwayssafe80, alwayssafe100,
   no_info_random0, no_info_random20, no_info_random40, no_info_random60, no_info_random80, no_info_random100
   )

choices <- choices  %>%
  select(starts_with(c("matching", "participant", "session", "time", "date", "simulation", "pct_bayesian")))

thenames <- names(choices)
function.clean_varnames <-  function(text){
  if(grepl("player", text, perl=TRUE | grepl("matching", text, perl=TRUE))){
    if(grepl("matching", text, perl=TRUE)){text <- sub(pattern = "matching.", replacement = "", x = text)}
    if(grepl(".player", text, perl=TRUE)){text <- sub(pattern = ".player", replacement = "", x = text)} 
    text <- strsplit(text, '[.]')
    text <- text[[1]][c(2,1)]
    text <- paste(text, collapse='_')
  }
  else {text <- text}} 
thenames <- lapply(thenames, function.clean_varnames)
names(choices) <- thenames 
rm(thenames, function.clean_varnames)

choices <- choices %>% 
  rename(treatment = matching.3.group.treatment,
         id_in_group = id_in_group_1,
         decision1 = pref_1,
         decision2 = pref_2,
         rank1 = position1_1,
         rank2 = position2_1,
         allocation1 = allocation_1,
         allocation2 = allocation_2,
         success1 = successful_1,
         success2 = successful_2,
         payoff1 = payoff_1,
         payoff2 = payoff_2,
         group_id = matching.1.group.group_id) %>% 
  mutate(treatment = factor(treatment, levels = c("No Info", "Full Info", "Partial Info")),
         rank0 = 3, # by definition
         valuationsb0 = as.numeric(substr(valuations0_1,6,7)),
         valuationsb1 = as.numeric(case_when(substr(valuations1_1,7,9) == "0]" ~ "0",
                                  TRUE ~ substr(valuations1_1,7,9))),
         valuationsb2 = as.numeric(case_when(substr(valuations2_1,7,9) == "0]" ~ "0",
                                  TRUE ~ substr(valuations2_1,7,9))),
         hint0 = case_when(treatment == "Partial Info" ~ substr(hint_1, 3,3), TRUE ~ NA_character_),
         hint1 = case_when(treatment == "Partial Info" ~ substr(hint_1, 8,8), TRUE ~ NA_character_),
         hint2 = case_when(treatment == "Partial Info" ~ substr(hint_1, 13,13), TRUE ~ NA_character_)) %>%
  group_by(group_id) %>% 
  mutate(group_payoff1 = sum(payoff1),
         group_payoff2 = sum(payoff2),
         group_completed = case_when(group_payoff1 == 0 ~ 0,
                                     TRUE ~ 1),
         group_completed = factor(group_completed, levels = c(0,1), labels = c("Not yet", "Completed"))) %>%
  select(c("simulation", "pct_bayesian", 
    "group_id", "treatment", "group_completed", "group_payoff1", "group_payoff2", "id_in_group", 
    "rank1", "hint1", "decision1", "allocation1", "success1", "valuationsb1", "payoff1", 
    "rank2", "hint2", "decision2", "allocation2", "success2", "valuationsb2", "payoff2",  
    )) 

choices <- choices[order(choices$group_id, choices$id_in_group),] 

Now the bots are loaded, we need to calculate the payoffs in R:

choices1 <- choices %>% 
  mutate(decisionfullinfo1 = case_when(valuationsb1 > 300 ~ "B",
                                       valuationsb1 == 300 ~ sample(c("A", "B"), 1), # if both are the same, select randomly
                                       TRUE ~ "A"),
         decisionfullinfo2 = case_when(valuationsb2 > 300 ~ "B",
                                       valuationsb2 == 300 ~ sample(c("A", "B"), 1), # if both are the same, select randomly
                                       TRUE ~ "A"),
         allocationfullinfo1 = case_when(rank1 < 3 ~ decisionfullinfo1,
                                         TRUE ~ ""), # if you are a high rank you always get what you choose)
         allocationfullinfo2 = case_when(rank2 < 3 ~ decisionfullinfo2,
                                         TRUE ~ "") 
  )


choices2 <- choices1 %>% 
  group_by(group_id) %>% 
  summarise(string = paste0(paste0(allocationfullinfo1, collapse = '')),
            string02 = paste0(paste0(allocationfullinfo2, collapse = ''))) %>% 
  ungroup() 

choices3 <- merge(choices1, choices2, by = "group_id") 
rm(choices1, choices2)
choices3 <- choices3 %>% 
  mutate(allocationfullinfo1 = case_when(rank1 > 2 & string == "AA" ~ "B",
                                         rank1 > 2 & string == "BB" ~ "A",
                                         rank1 == 3 & string == "AB" ~ decisionfullinfo1,
                                         rank1 == 3 & string == "BA" ~ decisionfullinfo1,
                                         TRUE ~ allocationfullinfo1),
         allocationfullinfo2 = case_when(rank2 > 2 & string02 == "AA" ~ "B",
                                         rank2 > 2 & string02 == "BB" ~ "A",
                                         rank2 == 3 & string02 == "AB" ~ decisionfullinfo2,
                                         rank2 == 3 & string02 == "BA" ~ decisionfullinfo2,
                                         TRUE ~ allocationfullinfo2))

choicescheck <- choices3 %>% select(group_id, rank1, string, allocationfullinfo1, decisionfullinfo1, valuationsb1,
                                    rank2, string02, allocationfullinfo2, decisionfullinfo2, valuationsb2)
# only rank 4 has no allocation here

choices4 <- choices3 %>% 
  group_by(group_id) %>% 
  summarise(string2 = paste0(paste0(allocationfullinfo1, collapse = '')),
            string22 = paste0(paste0(allocationfullinfo2, collapse = ''))) %>% 
  ungroup() 

choices5 <- merge(choices4, choices3, by = "group_id")

choices6 <- choices5 %>% select(rank1, string, string2, allocationfullinfo1, group_id, decisionfullinfo1, valuationsb1,
                                rank2, string02, string22, allocationfullinfo2, decisionfullinfo2, valuationsb2) 

choices6 <- choices6 %>% 
  mutate(allocationfullinfo1 = case_when(rank1 == 4 & string2 == "AAB"  ~ "B",
                                         rank1 == 4 & string2 == "ABA"  ~ "B",
                                         rank1 == 4 & string2 == "BAA"  ~ "B",
                                         rank1 == 4 & string2 == "BBA"  ~ "A",
                                         rank1 == 4 & string2 == "BAB"  ~ "A",
                                         rank1 == 4 & string2 == "ABB"  ~ "A",
                                         TRUE ~ allocationfullinfo1),
         payofffullinfo1 = case_when(allocationfullinfo1 == "A" ~ 300,
                                     allocationfullinfo1 == "B" ~ valuationsb1),
         allocationfullinfo2 = case_when(rank2 == 4 & string22 == "AAB"  ~ "B",
                                         rank2 == 4 & string22 == "ABA"  ~ "B",
                                         rank2 == 4 & string22 == "BAA"  ~ "B",
                                         rank2 == 4 & string22 == "BBA"  ~ "A",
                                         rank2 == 4 & string22 == "BAB"  ~ "A",
                                         rank2 == 4 & string22 == "ABB"  ~ "A",
                                         TRUE ~ allocationfullinfo2),
         payofffullinfo2 = case_when(allocationfullinfo2 == "A" ~ 300,
                                     allocationfullinfo2 == "B" ~ valuationsb2)) %>% 
  group_by(group_id) %>% 
  mutate(grouppayfullinfo1 = sum(payofffullinfo1),
         grouppayfullinfo2 = sum(payofffullinfo2)) %>% 
  select(grouppayfullinfo1, grouppayfullinfo2, group_id) %>% 
  slice(1)


choices <- merge(choices, choices6, by = "group_id")
rm(choices1, choices3, choices4, choices5, choices6, choicescheck)


max1 <- choices %>% 
  group_by(group_id) %>% 
  arrange(desc(valuationsb1), .by_group = TRUE) %>% 
  do(head(.,2)) %>% 
  mutate(maximum1 = sum(valuationsb1) + 2 * 300) %>% 
  filter(row_number()==1) %>% 
  ungroup() %>% 
  select(group_id, maximum1) 

max2 <- choices %>% 
  group_by(group_id) %>% 
  arrange(desc(valuationsb2), .by_group = TRUE) %>% 
  do(head(.,2)) %>% 
  mutate(maximum2 = sum(valuationsb2) + 2 * 300) %>%
  filter(row_number()==1) %>% 
  ungroup() %>% 
  select(group_id, maximum2) 

min1 <- choices %>% 
  group_by(group_id) %>% 
  arrange(desc(-valuationsb1), .by_group = TRUE) %>% 
  do(head(.,2)) %>% 
  mutate(minimum1 = sum(valuationsb1) + 2 * 300) %>% 
  filter(row_number()==1) %>% 
  ungroup() %>% 
  select(group_id, minimum1) 

min2 <- choices %>% 
  group_by(group_id) %>% 
  arrange(desc(-valuationsb2), .by_group = TRUE) %>% 
  do(head(.,2)) %>% 
  mutate(minimum2 = sum(valuationsb2) + 2 * 300) %>%
  filter(row_number()==1) %>% 
  ungroup() %>% 
  select(group_id, minimum2) 

choices <- merge(choices, max1, by = "group_id")
choices <- merge(choices, max2, by = "group_id")
choices <- merge(choices, min1, by = "group_id")
choices <- merge(choices, min2, by = "group_id")
rm(max1, max2, min1, min2)

group_payoffs <- choices %>% 
  filter(treatment != "No Info") %>% 
  filter(id_in_group == 1) %>% 
  select(treatment, group_payoff1, maximum1, group_payoff2, maximum2, minimum1, minimum2, group_id, simulation, pct_bayesian, grouppayfullinfo1, grouppayfullinfo2
         ) %>% 
  mutate(maximum = "Maximum")
group_payoffs_ni <- choices %>%
  filter(treatment == "No Info") %>%
  filter(id_in_group == 1) %>%
  select(treatment, group_payoff1, maximum1, group_payoff2, maximum2, minimum1, minimum2, group_id, simulation, pct_bayesian, grouppayfullinfo1, grouppayfullinfo2) %>%
  mutate(maximum = "Maximum")


plot_ni1 <- group_payoffs_ni %>%
  group_by(pct_bayesian, simulation) %>%
  summarise(payoff1 = mean(group_payoff1, na.rm = TRUE),
            payoff2 = mean(group_payoff2, na.rm = TRUE),
            payoff1_sd = sd(group_payoff1, na.rm = TRUE),
            payoff2_sd = sd(group_payoff2, na.rm = TRUE)) %>%
  mutate(group = "Partial + No Info",
         pct_bayesian = pct_bayesian/100,
         pct_type = 1 - pct_bayesian)

plot_ni2 <- group_payoffs_ni %>%
  group_by(pct_bayesian, simulation) %>%
  summarise(payoff1 = mean(grouppayfullinfo1, na.rm = TRUE),
            payoff2 = mean(grouppayfullinfo2, na.rm = TRUE),
            payoff1_sd = sd(grouppayfullinfo1, na.rm = TRUE),
            payoff2_sd = sd(grouppayfullinfo2, na.rm = TRUE)) %>%
  mutate(group = "Full Info",
         pct_bayesian = pct_bayesian/100,
         pct_type = 1 - pct_bayesian)

plot_ni3 <- group_payoffs_ni %>%
  group_by(pct_bayesian, simulation) %>%
  summarise(payoff1 = mean(maximum1, na.rm = TRUE),
            payoff2 = mean(maximum2, na.rm = TRUE),
            payoff1_sd = sd(maximum1, na.rm = TRUE),
            payoff2_sd = sd(maximum2, na.rm = TRUE)) %>%
  mutate(group = "Maximum",
         pct_bayesian = pct_bayesian/100,
         pct_type = 1 - pct_bayesian)

plot_ni4 <- group_payoffs_ni %>%
  group_by(pct_bayesian, simulation) %>%
  summarise(payoff1 = mean(minimum1, na.rm = TRUE),
            payoff2 = mean(minimum2, na.rm = TRUE),
            payoff1_sd = sd(minimum1, na.rm = TRUE),
            payoff2_sd = sd(minimum2, na.rm = TRUE)) %>%
  mutate(group = "Minimum",
         pct_bayesian = pct_bayesian/100,
         pct_type = 1 - pct_bayesian)

plot_ni <- rbind(plot_ni1, plot_ni2, plot_ni3, plot_ni4) 
plot_ni$group <- ordered(plot_ni$group, levels = c("Maximum", "Full Info", "Partial + No Info", "Minimum"))  # to change the order
plot_ni <- plot_ni[order(plot_ni$group, plot_ni$pct_type),] %>%
  ungroup() %>%
  select(-pct_bayesian) %>% 
  mutate(simulation = "Always safe option A")

write.xlsx(as.data.frame(plot_ni), here("output", "simulations", "simulations_noinfo.xlsx"), row.names=FALSE)

cbp1 <- c("#999999", "#E69F00", "#CC79A7", "#009E73")
shapes <- c(24, 22, 21, 25)

decision1ni <- ggplot(data=plot_ni, aes(x=pct_type, y=payoff1, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1) + theme_minimal() +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes  
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 1",
       subtitle = "A = 300, B = [100, 300, 500]",
       color = "", fill = "", linetype = "", shape = "") +
  theme(legend.position = "none", axis.title.y = element_text(size=12), plot.title = element_text(size = 14, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "italic"))


decision2ni <- ggplot(data=plot_ni, aes(x=pct_type, y=payoff2, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1)  + theme_minimal() +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 2",
       color = "Simulation", fill = "Simulation", linetype = "Simulation", shape = "Simulation") +
  theme(legend.position = "right")

# Create user-defined function, which extracts legends from ggplots
extract_legend <- function(my_ggp) {
  step1 <- ggplot_gtable(ggplot_build(my_ggp))
  step2 <- which(sapply(step1$grobs, function(x) x$name) == "guide-box")
  step3 <- step1$grobs[[step2]]
  return(step3)
}

legendni <- extract_legend(decision2ni)

decision2ni <- ggplot(data=plot_ni, aes(x=pct_type, y=payoff2, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1)  + theme_minimal() +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 2",
       subtitle = "A = 300, B = [0, 200, 400]",
       color = "", fill = "", linetype = "", shape = "") +
  theme(legend.position = "none", axis.title.y = element_text(size=12), plot.title = element_text(size = 14, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "italic"))

g <- grid.arrange(decision1ni, decision2ni, legendni, nrow=1, widths= c(2,2,1))

ggsave(here("output", "simulations", "simulations_noinfo.pdf"), g, device = 'pdf',
       width = 20, height = 10, units = "cm")

Figure C2: Simulations always safe

Export the figures:

plot1 <- group_payoffs %>% 
  group_by(pct_bayesian, simulation) %>% 
  summarise(payoff1 = mean(maximum1, na.rm = TRUE),
            payoff2 = mean(maximum2, na.rm = TRUE),
            payoff1_sd = sd(maximum1, na.rm = TRUE),
            payoff2_sd = sd(maximum2, na.rm = TRUE)) %>% 
  mutate(group = "Maximum" , 
         pct_bayesian = pct_bayesian/100,
         pct_type = 1 - pct_bayesian)  

plot1a <- group_payoffs %>% 
  group_by(pct_bayesian, simulation) %>% 
  summarise(payoff1 = mean(minimum1, na.rm = TRUE),
            payoff2 = mean(minimum2, na.rm = TRUE),
            payoff1_sd = sd(minimum1, na.rm = TRUE),
            payoff2_sd = sd(minimum2, na.rm = TRUE)) %>% 
  mutate(group = "Minimum" , 
         pct_bayesian = pct_bayesian/100,
         pct_type = 1 - pct_bayesian)  

plot2 <- group_payoffs %>% 
  group_by(pct_bayesian, simulation) %>% 
  summarise(payoff1 = mean(group_payoff1, na.rm = TRUE),
            payoff2 = mean(group_payoff2, na.rm = TRUE),
            payoff1_sd = sd(group_payoff1, na.rm = TRUE),
            payoff2_sd = sd(group_payoff2, na.rm = TRUE)) %>% 
  mutate(group = "Partial Info", 
         pct_bayesian = pct_bayesian/100,
         pct_type = 1 - pct_bayesian)  

plot3 <- group_payoffs %>% 
  group_by(pct_bayesian, simulation) %>% 
  summarise(payoff1 = mean(grouppayfullinfo1, na.rm = TRUE),
            payoff2 = mean(grouppayfullinfo2, na.rm = TRUE),
            payoff1_sd = sd(grouppayfullinfo1, na.rm = TRUE),
            payoff2_sd = sd(grouppayfullinfo2, na.rm = TRUE)) %>% 
  mutate(group = "Full Info",
         pct_bayesian = pct_bayesian/100,
         pct_type = 1 - pct_bayesian)  

plot <- rbind(plot1, plot1a, plot2, plot3) %>% 
  mutate(simulation = case_when(simulation == "defiers" ~ "Always opposite to Recommendation",
                                simulation == "overcompliers" ~ "Always follow Recommendation",
                                simulation == "alwayssafe" ~ "Always safe option A"))
plot$group <- factor(plot$group, levels = c("Maximum", "Full Info", "Partial Info", "Minimum"))  # to change the order
plot <- plot[order(plot$simulation, plot$group, plot$pct_type),] %>% 
  ungroup() %>% 
  select(-pct_bayesian)

write.xlsx(as.data.frame(plot), here("output", "simulations", "simulations.xlsx"), row.names=FALSE)

plot_original <- plot
rm(plot1, plot1a, plot2, plot3)

# The palette with grey:
cbp1 <- c("#999999", "#E69F00", "#56B4E9", "#009E73")
shapes <- c(24, 22, 23, 25)
  
decision1_a <- ggplot(data=subset(plot_original, simulation == "Always follow Recommendation"), 
                           aes(x=pct_type, y=payoff1, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1)  + theme_minimal() +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 1", 
       subtitle = "A = 300, B = [100, 300, 500]",
       color = "", fill = "", linetype = "", shape = "") +
  theme(legend.position = "none", axis.title.y = element_text(size=12), plot.title = element_text(size = 14, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "italic"))

decision1_b <- ggplot(data=subset(plot_original, simulation == "Always opposite to Recommendation"), 
                           aes(x=pct_type, y=payoff1, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1)  + theme_minimal() +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 1", 
       subtitle = "A = 300, B = [100, 300, 500]",
       color = "", fill = "", linetype = "", shape = "") +
  theme(legend.position = "none", axis.title.y = element_text(size=12), plot.title = element_text(size = 14, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "italic"))

decision1_c <- ggplot(data=subset(plot_original, simulation == "Always safe option A"), 
                           aes(x=pct_type, y=payoff1, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1)  + theme_minimal() +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 1", 
       subtitle = "A = 300, B = [100, 300, 500]",
       color = "", fill = "", linetype = "", shape = "") +
  theme(legend.position = "none", axis.title.y = element_text(size=12), plot.title = element_text(size = 14, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "italic"))

decision2_a <- ggplot(data=subset(plot_original, simulation == "Always follow Recommendation"), 
                      aes(x=pct_type, y=payoff2, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1)  + theme_minimal() +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 2",
       subtitle = "A = 300, B = [0, 200, 400]",
       color = "", fill = "", linetype = "", shape = "") +
  theme(legend.position = "none", axis.title.y = element_text(size=12), plot.title = element_text(size = 14, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "italic"))

decision2_b <- ggplot(data=subset(plot_original, simulation == "Always opposite to Recommendation"), 
                      aes(x=pct_type, y=payoff2, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1)  + theme_minimal() +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 2",
       subtitle = "A = 300, B = [0, 200, 400]",
       color = "", fill = "", linetype = "", shape = "") +
  theme(legend.position = "none", axis.title.y = element_text(size=12), plot.title = element_text(size = 14, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "italic"))

decision2_c <- ggplot(data=subset(plot_original, simulation == "Always safe option A"), 
                      aes(x=pct_type, y=payoff2, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1)  + theme_minimal() +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 2",
       subtitle = "A = 300, B = [0, 200, 400]",
       color = "", fill = "", linetype = "", shape = "") +
  theme(legend.position = "none", axis.title.y = element_text(size=12), plot.title = element_text(size = 14, face = "bold"),
        plot.subtitle = element_text(size = 10, face = "italic"))

decision2 <- ggplot(data=plot_original, aes(x=pct_type, y=payoff2, color=group, fill = group)) +
  geom_line(aes(linetype = group)) + 
  geom_point(aes(shape = group)) +
  facet_wrap(~simulation, ncol = 1)  + theme_minimal() +
  scale_y_continuous(limits = c(800, 1400)) +
  scale_x_continuous(labels = scales::percent_format(accuracy = 1L),
                     limits = c(0,1),
                     breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1)) +
  scale_colour_manual(values=cbp1) +
  scale_fill_manual(values=cbp1) +
  scale_shape_manual(values=shapes) + # to the the requested shapes
  labs(x = "\n percentage of each type", y = "\nGroup payoff", title = "Scenario 2",
       color = "Simulation", fill = "Simulation", linetype = "Simulation", shape = "Simulation") +
  theme(legend.position = "right")


legend <- extract_legend(decision2)

g3 <- grid.arrange(decision1_c, decision2_c, legend, nrow=1, widths= c(2,2,1))

ggsave(here("output", "simulations", "simulations3.pdf"), g3, device = 'pdf',
       width = 20, height = 9, units = "cm")

Figure C4: Simulations always follow

g1 <- grid.arrange(decision1_a, decision2_a, legend, nrow=1, widths= c(2,2,1))

ggsave(here("output", "simulations", "simulations1.pdf"), g1, device = 'pdf',
       width = 20, height = 9, units = "cm")

Figure C5: Simulations opposite to Recommendation

g2 <- grid.arrange(decision1_b, decision2_b, legend, nrow=1, widths= c(2,2,1))

ggsave(here("output", "simulations", "simulations2.pdf"), g2, device = 'pdf',
       width = 20, height = 9, units = "cm")

Appendix F: Simulations for Policy Analysis

   simulations2 <- df %>% 
    group_by(group_id) %>% 
    arrange(desc(valuationsb1), desc(redistribution), .by_group = TRUE) %>% 
    filter(treatment == "Partial Info") %>% 
    summarise(altruist_num = sum(redistribution),
              valuationsb1_ = paste0(paste0(paste0(paste0(valuationsb1, collapse = ' ')))),
              altruists_ = paste0(paste0(paste0(paste0(redistribution, collapse = '')))),
              hints1_original = as.character(paste0(paste0(paste0(paste0(hint1, collapse = ''))))),
              num_in_groupb1 = row_number(),
              id_otree = id_otree
    ) %>% 
    mutate(valuationb1_1 = as.numeric(substr(valuationsb1_, 1, 4)),
           valuationb1_2 = as.numeric(substr(valuationsb1_, 5, 8)),
           valuationb1_3 = as.numeric(substr(valuationsb1_, 9, 12)),
           valuationb1_4 = as.numeric(substr(valuationsb1_, 13, 16)),
           altruist1_1 = substr(altruists_, 1, 1),
           altruist1_2 = substr(altruists_, 2, 2),
           altruist1_3 = substr(altruists_, 3, 3),
           altruist1_4 = substr(altruists_, 4, 4))%>% 
    mutate(ties1 = case_when(valuationb1_1 == valuationb1_4 ~ paste0("all equal ", as.character(valuationb1_2)),
                            valuationb1_1 == valuationb1_2 & valuationb1_2 == valuationb1_3 ~ paste0("first three ", as.character(valuationb1_2)),
                            valuationb1_2 == valuationb1_3 & valuationb1_3 == valuationb1_4 ~ paste0("last three ", as.character(valuationb1_2)),
                            valuationb1_2 == valuationb1_3 ~ paste0("middle two ", as.character(valuationb1_3)), 
                            TRUE ~ "no ties")) # no relevant ties (last or first two can be tied but we have two seats)
  
# check ties:
# table(simulations2$valuationsb1_, simulations2$ties1) 

simulations2 <- simulations2 %>%
  select(-starts_with("valuation"), -altruists_)
  
simulations22 <- df %>%
  group_by(group_id) %>%
  arrange(desc(valuationsb2), desc(redistribution), .by_group = TRUE) %>%
  filter(treatment == "Partial Info") %>%
  summarise(
    valuationsb2_ = paste0(paste0(paste0(paste0(valuationsb2, collapse = ' ')))),
    altruists_ = paste0(paste0(paste0(paste0(redistribution, collapse = '')))),
    hints2_original = as.character(paste0(paste0(paste0(paste0(hint2, collapse = ''))))),
    num_in_group_b2 = row_number(),
    id_otree = id_otree
  ) %>%
  mutate(
    valuationb2_1 = as.numeric(strsplit(valuationsb2_," ")[[1]][1]),
    valuationb2_2 = as.numeric(strsplit(valuationsb2_," ")[[1]][2]),
    valuationb2_3 = as.numeric(strsplit(valuationsb2_," ")[[1]][3]),
    valuationb2_4 = as.numeric(strsplit(valuationsb2_," ")[[1]][4]),
    altruist_1 = substr(altruists_, 1, 1),
    altruist_2 = substr(altruists_, 2, 2),
    altruist_3 = substr(altruists_, 3, 3),
    altruist_4 = substr(altruists_, 4, 4)
  ) %>%
  mutate(
    ties2 = case_when(
      valuationb2_1 == valuationb2_4 ~ paste0("all equal ", as.character(valuationb2_2)),
      valuationb2_1 == valuationb2_2 & valuationb2_2 == valuationb2_3 ~ paste0("first three ", as.character(valuationb2_2)),
      valuationb2_2 == valuationb2_3 & valuationb2_3 == valuationb2_4 ~ paste0("last three ", as.character(valuationb2_2)),
      valuationb2_2 == valuationb2_3 ~ paste0("middle two ", as.character(valuationb2_3)), 
      TRUE ~ "no ties"
    )
  ) 


# check ties:
# table(simulations22$valuationsb2_, simulations22$ties2) 

simulations22 <- simulations22 %>%
  select(-starts_with("valuation"))
  
simulations <- merge(simulations2, simulations22, by = "id_otree") 
  
rm(simulations22, simulations2)
  
set.seed(1234) # set seed because some random samples are picked
simulations <- merge(simulations, df, by = "id_otree")  %>% 
    mutate(group_id = group_id.x) %>% 
    select(-group_id.x, -group_id.y) %>% 
    group_by(group_id) %>% 
    mutate(hints1_alternative = case_when(# first the case with only 1 altruist, player 4 is the other to take the A hint, because B leads to 500 which is higher
                                          ties1 == "first three 500" & altruist_num == 1 & altruist1_1 == 1 ~ "ABBA", # altruist1 can take the A hint 
                                          ties1 == "first three 500" & altruist_num == 1 & altruist1_2 == 1 ~ "BABA", # altruist2 can take the A hint
                                          ties1 == "first three 500" & altruist_num == 1 & altruist1_3 == 1 ~ "BBAA", # altruist3 can take the A hint
                                          # if we have two altruists, randomly pick one of the two (plus player 4)
                                          ties1 == "first three 500" & altruist_num == 2 & altruist1_1 == 1 & altruist1_2 == 1 ~ sample(c("ABBA", "BABA"),1),
                                          ties1 == "first three 500" & altruist_num == 2 & altruist1_1 == 1 & altruist1_3 == 1 ~ sample(c("ABBA", "BBAA"),1),
                                          ties1 == "first three 500" & altruist_num == 2 & altruist1_2 == 1 & altruist1_3 == 1 ~ sample(c("BABA", "BBAA"),1),
                                          # if all three are altruists, randomly pick one of the three (plus player 4)
                                          ties1 == "first three 500" & altruist_num == 3 & altruist_4 == 0 ~ sample(c("ABBA", "BABA", "BBAA"), 1),
                                          # first the case with only 1 altruist, player 1 for sure gets the B hint, so we need to divide the second B hint among the last three
                                          # if the first one is the altruist, the original hints suffice
                                          ties1 == "last three 100" & altruist_num == 1 & altruist1_2 == 1 ~ "BBAA",
                                          ties1 == "last three 100" & altruist_num == 1 & altruist1_3 == 1 ~ "BABA",
                                          ties1 == "last three 100" & altruist_num == 1 & altruist1_4 == 1 ~ "BAAB",
                                          # with 2 altruists including player 1, the altruist gets the B
                                          ties1 == "last three 100" & altruist_num == 2 & altruist1_1 == 1 & altruist1_2 == 1 ~ "BBAA",
                                          ties1 == "last three 100" & altruist_num == 2 & altruist1_1 == 1 & altruist1_3 == 1 ~ "BABA",
                                          ties1 == "last three 100" & altruist_num == 2 & altruist1_1 == 1 & altruist1_4 == 1 ~ "BAAB",
                                          # with 2 altruists excluding player 1, randomly pick one of the altruists for B
                                          ties1 == "last three 100" & altruist_num == 2 & altruist1_2 == 1 & altruist1_3 == 1 ~ sample(c("BBAA", "BABA"),1),
                                          ties1 == "last three 100" & altruist_num == 2 & altruist1_2 == 1 & altruist1_4 == 1 ~ sample(c("BBAA", "BAAB"),1),
                                          ties1 == "last three 100" & altruist_num == 2 & altruist1_4 == 1 & altruist1_3 == 1 ~ sample(c("BABA", "BAAB"),1),
                                          # with 3 altruists including player 1, randomly pick one of the altruists for B
                                          ties1 == "last three 100" & altruist_num == 3 & altruist1_2 == 0 ~ sample(c("BABA", "BAAB"),1),
                                          ties1 == "last three 100" & altruist_num == 3 & altruist1_3 == 0 ~ sample(c("BBAA", "BAAB"),1),
                                          ties1 == "last three 100" & altruist_num == 3 & altruist1_4 == 0 ~ sample(c("BBAA", "BABA"),1),
                                          # with 3 altruists excluding player 1, the original hint is fine
                                          # first three 300 does not matter because 300 in B is equal to 300 in A, so original hints are fine
                                          # all equal 300 does not matter for the same reason
                                          # last three 300 does not matter either for the same reason
                                          # middle two 300 does not matter either for the same reason
                                          TRUE ~ hints1_original), # when no ties or everyone altruist, original hints suffice),
           hints1_different = case_when(hints1_original == hints1_alternative ~ "equal", TRUE ~ "new"),
           hints2_alternative = case_when(altruist_num == 4 ~ hints2_original,
                                          altruist_num == 0 ~ hints2_original,
                                          # all equal 200 means the better option is A (300) so the altruists take the hit by B (keep reading next lines)
                                          # first the case with only 1 altruist, still need to randomly pick the other
                                          ties2 == "all equal 200" & altruist_num == 1 & altruist_1 == 1 ~ sample(c("BBAA", "BABA", "BAAB"), 1),
                                          ties2 == "all equal 200" & altruist_num == 1 & altruist_2 == 1 ~ sample(c("BBAA", "ABBA", "ABAB"), 1),
                                          ties2 == "all equal 200" & altruist_num == 1 & altruist_3 == 1 ~ sample(c("AABB", "ABBA", "BABA"), 1),
                                          ties2 == "all equal 200" & altruist_num == 1 & altruist_4 == 1 ~ sample(c("AABB", "ABAB", "BAAB"), 1),   
                                          # with exactly two altruists we can divide according to altruism
                                          ties2 == "all equal 200" & altruist_num == 2 & altruist_1 == 1 & altruist_2 == 1 ~ "BBAA",
                                          ties2 == "all equal 200" & altruist_num == 2 & altruist_1 == 1 & altruist_3 == 1 ~ "BABA",
                                          ties2 == "all equal 200" & altruist_num == 2 & altruist_1 == 1 & altruist_4 == 1 ~ "BAAB",
                                          ties2 == "all equal 200" & altruist_num == 2 & altruist_2 == 1 & altruist_3 == 1 ~ "ABBA",
                                          ties2 == "all equal 200" & altruist_num == 2 & altruist_2 == 1 & altruist_4 == 1 ~ "ABAB",
                                          ties2 == "all equal 200" & altruist_num == 2 & altruist_3 == 1 & altruist_4 == 1 ~ "AABB",
                                          # with one altruist too many, we need to randomly pick one of the three to get the good A
                                          ties2 == "all equal 200" & altruist_num == 3 & altruist_1 == 0 ~ sample(c("ABBA", "ABAB", "AABB"),1),
                                          ties2 == "all equal 200" & altruist_num == 3 & altruist_2 == 0 ~ sample(c("AABB", "BABA", "BAAB"),1),
                                          ties2 == "all equal 200" & altruist_num == 3 & altruist_3 == 0 ~ sample(c("ABAB", "BAAB", "BBAA"),1),
                                          ties2 == "all equal 200" & altruist_num == 3 & altruist_4 == 0 ~ sample(c("BBAA", "BABA", "ABBA"),1),
                                          # if all 4 are altruists, we can stick to the original hints
                                          # all equal 400 means the better option is 400 = B (versus A = 300) so the altruists can take the hit by A
                                          ties2 == "all equal 400" & altruist_num == 1 & altruist_1  == 1 ~ sample(c("AABB", "ABAB", "ABBA"), 1),
                                          ties2 == "all equal 400" & altruist_num == 1 & altruist_2  == 1 ~ sample(c("AABB", "BAAB", "BABA"), 1),
                                          ties2 == "all equal 400" & altruist_num == 1 & altruist_3  == 1 ~ sample(c("BBAA", "BAAB", "ABAB"), 1),
                                          ties2 == "all equal 400" & altruist_num == 1 & altruist_4  == 1 ~ sample(c("BBAA", "BABA", "ABBA"), 1),  
                                          # if we have two altruists we can divide according to altruism
                                          ties2 == "all equal 400" & altruist_num == 2 & altruist_1 == 1 & altruist_2 == 1 ~ "AABB",
                                          ties2 == "all equal 400" & altruist_num == 2 & altruist_1 == 1 & altruist_3 == 1 ~ "ABAB",
                                          ties2 == "all equal 400" & altruist_num == 2 & altruist_1 == 1 & altruist_4 == 1 ~ "ABBA",
                                          ties2 == "all equal 400" & altruist_num == 2 & altruist_2 == 1 & altruist_3 == 1 ~ "BAAB",
                                          ties2 == "all equal 400" & altruist_num == 2 & altruist_2 == 1 & altruist_4 == 1 ~ "BABA",
                                          ties2 == "all equal 400" & altruist_num == 2 & altruist_3 == 1 & altruist_4 == 1 ~ "BBAA",
                                          # if we have three altruists, randomly pick two altruists to receive A
                                          ties2 == "all equal 400" & altruist_num == 3 & altruist_1 == 0 ~ sample(c("BAAB", "BABA", "BBAA"),1),
                                          ties2 == "all equal 400" & altruist_num == 3 & altruist_2 == 0 ~ sample(c("BBAA", "ABAB", "ABBA"),1),
                                          ties2 == "all equal 400" & altruist_num == 3 & altruist_3 == 0 ~ sample(c("BABA", "ABBA", "AABB"),1),
                                          ties2 == "all equal 400" & altruist_num == 3 & altruist_4 == 0 ~ sample(c("AABB", "ABAB", "BAAB"),1),
                                          # first three 400, altruist among those can take the A hint
                                          ties2 == "first three 400" & altruist_num == 1 & altruist_1 == 1 ~ "ABBA", 
                                          ties2 == "first three 400" & altruist_num == 1 & altruist_2 == 1 ~ "BABA", 
                                          ties2 == "first three 400" & altruist_num == 1 & altruist_3 == 1 ~ "BBAA", 
                                          # when number 4 is an altruist, this person gets the A anyway, so randomly pick A from the first 3, original hints suffice
                                          # with two altruists in the first three, non-altruist in first three gets B (400), randomly pick one other
                                          ties2 == "first three 400" & altruist_num == 2 & altruist_4 == 0 & altruist_1 == 0 ~ sample(c("BBAA", "BABA"),1), 
                                          ties2 == "first three 400" & altruist_num == 2 & altruist_4 == 0 & altruist_2 == 0 ~ sample(c("BBAA", "ABBA"),1), 
                                          ties2 == "first three 400" & altruist_num == 2 & altruist_4 == 0 & altruist_3 == 0 ~ sample(c("ABBA", "BABA"),1),  
                                          # with two altruists one of which is number 4, the altruist gets the A
                                          ties2 == "first three 400" & altruist_num == 2 & altruist_4 == 1 & altruist_1 == 1 ~ "ABBA", 
                                          ties2 == "first three 400" & altruist_num == 2 & altruist_4 == 1 & altruist_2 == 1 ~ "BABA", 
                                          ties2 == "first three 400" & altruist_num == 2 & altruist_4 == 1 & altruist_3 == 1 ~ "BBAA",  
                                          # with three altruists in the first three, does not matter
                                          # with three altruists including number four, randomly pick one of the two altruists to get the A
                                          ties2 == "first three 400" & altruist_num == 3 & altruist_4 == 1 & altruist_1 == 0 ~ sample(c("BBAA", "BABA"),1), 
                                          ties2 == "first three 400" & altruist_num == 3 & altruist_4 == 1 & altruist_2 == 0 ~ sample(c("BBAA", "ABBA"),1), 
                                          ties2 == "first three 400" & altruist_num == 3 & altruist_4 == 1 & altruist_3 == 0 ~ sample(c("ABBA", "BABA"),1),  
                                          # with four altruists original hints suffice
                                          # last three 0, first person gets a B (200 or 400) and altruist gets the other B (0)
                                          # if one altruist it only matters if this is not number 1
                                          ties2 == "last three 0" & altruist_num == 1 & altruist_2 == 1 ~ "BBAA",
                                          ties2 == "last three 0" & altruist_num == 1 & altruist_3 == 1 ~ "BABA",
                                          ties2 == "last three 0" & altruist_num == 1 & altruist_4 == 1 ~ "BAAB",
                                          # if two altruists of which one is number 1, altruist gets B
                                          ties2 == "last three 0" & altruist_num == 2 & altruist_1 == 1 & altruist_2 == 1 ~"BBAA",
                                          ties2 == "last three 0" & altruist_num == 2 & altruist_1 == 1 & altruist_3 == 1 ~"BABA",
                                          ties2 == "last three 0" & altruist_num == 2 & altruist_1 == 1 & altruist_4 == 1 ~"BAAB",
                                          # if two altruists in the last three, randomly pick one of them for B
                                          ties2 == "last three 0" & altruist_num == 2 & altruist_1 == 0 & altruist_2 == 0 ~ sample(c("BABA", "BAAB"),1),
                                          ties2 == "last three 0" & altruist_num == 2 & altruist_1 == 0 & altruist_3 == 0 ~ sample(c("BBAA", "BAAB"),1),
                                          ties2 == "last three 0" & altruist_num == 2 & altruist_1 == 0 & altruist_4 == 0 ~ sample(c("BABA", "BBAA"),1),
                                          # if three altruists in the last three, it does not matter
                                          # if three altruists including number 1, randomly pick one of them for B
                                          ties2 == "last three 0" & altruist_num == 3 & altruist_1 == 1 & altruist_2 == 0 ~ sample(c("BABA", "BAAB"),1),
                                          ties2 == "last three 0" & altruist_num == 3 & altruist_1 == 1 & altruist_3 == 0 ~ sample(c("BBAA", "BAAB"),1),
                                          ties2 == "last three 0" & altruist_num == 3 & altruist_1 == 1 & altruist_4 == 0 ~ sample(c("BABA", "BBAA"),1),
                                          # if four altruists, original hints suffice
                                          # first three 200, last one is 0. altruist + on random other will take the B hint, and last one always gets A
                                          ties2 == "first three 200" & altruist_num == 1 & altruist_1 == 1 ~ sample(c("BBAA", "BABA"),1), 
                                          ties2 == "first three 200" & altruist_num == 1 & altruist_2 == 1 ~ sample(c("BBAA", "ABBA"),1), 
                                          ties2 == "first three 200" & altruist_num == 1 & altruist_3 == 1 ~ sample(c("ABBA", "BABA"),1), 
                                          # if the fourth one is the altruist, original hints suffice
                                          # with two altruists in the top 3, both get the B
                                          ## NOT AN ISSUE, SOMETIMES IMPROVES WELFARE ASSUMING HINTS ARE FOLLOWED (group id a45aa)
                                          ties2 == "first three 200" & altruist_num == 2 & altruist_4 == 0 & altruist_1 == 0 ~ "ABBA",
                                          ties2 == "first three 200" & altruist_num == 2 & altruist_4 == 0 & altruist_2 == 0 ~ "BABA",
                                          ties2 == "first three 200" & altruist_num == 2 & altruist_4 == 0 & altruist_3 == 0 ~ "BBAA",
                                          # with two altruists one being number four, randomly pick one other (and the altruist) to get B
                                          ties2 == "first three 200" & altruist_num == 2 & altruist_4 == 1 & altruist_1 == 1 ~ sample(c("BBAA", "BABA"),1), 
                                          ties2 == "first three 200" & altruist_num == 2 & altruist_4 == 1 & altruist_2 == 1 ~ sample(c("BBAA", "ABBA"),1), 
                                          ties2 == "first three 200" & altruist_num == 2 & altruist_4 == 1 & altruist_3 == 1 ~ sample(c("ABBA", "BABA"),1), 
                                          # with three altruists in the first three, original hints suffice
                                          # with three altruists one of which number 4, both altruists get B
                                          ties2 == "first three 200" & altruist_num == 3 & altruist_1 == 0 ~ "ABBA",
                                          ties2 == "first three 200" & altruist_num == 3 & altruist_2 == 0 ~ "BABA",
                                          ties2 == "first three 200" & altruist_num == 3 & altruist_3 == 0 ~ "BBAA",
                                          # all four altruists original hints suffice
                                          # last three 200, meaning number 1 is 400 and gets B for sure, altruists get the B
                                          # if altruist is number 1 original hints suffice
                                          ties2 == "last three 200" & altruist_num == 1 & altruist_2 == 1 ~ "BBAA",
                                          ties2 == "last three 200" & altruist_num == 1 & altruist_3 == 1 ~ "BABA",
                                          ties2 == "last three 200" & altruist_num == 1 & altruist_4 == 1 ~ "BAAB",
                                          # with two altruists including player 1, altruist 2/3/4 gets the B
                                          ties2 == "last three 200" & altruist_num == 2 & altruist_1 == 1 & altruist_2 == 1 ~ "BBAA",
                                          ties2 == "last three 200" & altruist_num == 2 & altruist_1 == 1 & altruist_3 == 1 ~ "BABA",
                                          ties2 == "last three 200" & altruist_num == 2 & altruist_1 == 1 & altruist_4 == 1 ~ "BAAB",
                                          # with two altruists excluding player 1, non-altruist gets A, randomly pick another A among altruists
                                          ties2 == "last three 200" & altruist_num == 2 & altruist_1 == 0 & altruist_2 == 0 ~ sample(c("BABA", "BAAB"),1),
                                          ties2 == "last three 200" & altruist_num == 2 & altruist_1 == 0 & altruist_3 == 0 ~ sample(c("BBAA", "BAAB"),1),
                                          ties2 == "last three 200" & altruist_num == 2 & altruist_1 == 0 & altruist_4 == 0 ~ sample(c("BBAA", "BABA"),1), 
                                          # with three altruists excluding player 1, original hints suffice
                                          # three altruists including player 1, randomly pick one of the altruists to get A
                                          ties2 == "last three 200" & altruist_num == 3 & altruist_2 == 0 ~ sample(c("BABA", "BAAB"),1),
                                          ties2 == "last three 200" & altruist_num == 3 & altruist_3 == 0 ~ sample(c("BBAA", "BAAB"),1),
                                          ties2 == "last three 200" & altruist_num == 3 & altruist_4 == 0 ~ sample(c("BBAA", "BABA"),1),
                                          # with 4 altruists original hints suffice
                                          # 400 200 200 0, first always gets B, last always gets A, altruists among middle two get B
                                          # with one altruist it only matters if it is number 2 or 3
                                          # STILL ISSUES
                                          ties2 == "middle two 200" & altruist_num == 1 & altruist_2 == 1 ~ "BBAA",
                                          ties2 == "middle two 200" & altruist_num == 1 & altruist_3 == 1 ~ "BABA",
                                          # with two altruists, it only matters if one of them is number 2 or 3 and the other one is not (altruist gets B)
                                          ties2 == "middle two 200" & altruist_num == 2 & altruist_1 == 1 & altruist_2 == 1 ~ "BBAA", # this is group 3r9qg, and because we assume 
                                          ties2 == "middle two 200" & altruist_num == 2 & altruist_1 == 1 & altruist_3 == 1 ~ "BABA",
                                          # with one and four both altruists, original hints suffice (none of the middle two are altruists)
                                          # if two and three are both altruists, original hints suffice
                                          ties2 == "middle two 200" & altruist_num == 2 & altruist_2 == 1 & altruist_4 == 1 ~ "BBAA",
                                          ties2 == "middle two 200" & altruist_num == 2 & altruist_3 == 1 & altruist_4 == 1 ~ "BABA",
                                          # with three altruists only matters if number 2 and 3 have one altruist and one non-altruist (altruist gets B)
                                          ties2 == "middle two 200" & altruist_num == 3 & altruist_2 == 1 & altruist_3 == 0 ~ "BBAA",
                                          ties2 == "middle two 200" & altruist_num == 3 & altruist_2 == 0 & altruist_3 == 1 ~ "BABA",
                                          # with four altruists original hints suffice
                                          TRUE ~ hints2_original), # this defines that in all other cases, stick to original hints
          hints2_different = case_when(hints2_original == hints2_alternative ~ "equal", TRUE ~ "new"))
  
simulations <- simulations %>% 
  mutate(hint1_alt = case_when(num_in_groupb1 == 1 ~ substr(hints1_alternative, 1, 1),
                               num_in_groupb1 == 2 ~ substr(hints1_alternative, 2, 2),
                               num_in_groupb1 == 3 ~ substr(hints1_alternative, 3, 3),
                               num_in_groupb1 == 4 ~ substr(hints1_alternative, 4, 4),
                               TRUE ~ NA_character_),
         hint2_alt = case_when(num_in_group_b2 == 1 ~ substr(hints2_alternative, 1, 1),
                               num_in_group_b2 == 2 ~ substr(hints2_alternative, 2, 2),
                               num_in_group_b2 == 3 ~ substr(hints2_alternative, 3, 3),
                               num_in_group_b2 == 4 ~ substr(hints2_alternative, 4, 4),
                               TRUE ~ NA_character_),
         decision1_alt = case_when(follow_hint1 == 1 ~ hint1_alt,  # assume those who followed before will follow again
                                   TRUE ~ decision1),              # others will stick to their original decision
         decision2_alt = case_when(follow_hint2 == 1 ~ hint2_alt,
                                   TRUE ~ decision2),
         decision1_alt_ = hint1_alt, # assume everyone follows the new hint
         decision2_alt_ = hint2_alt,
         decision1_original_follow = hint1, # compare with original hints, what if everyone followed those
         decision2_original_follow = hint2,
         hint1_diff = case_when(hint1_alt == hint1 ~ "equal", TRUE ~ "different"),
         hint2_diff = case_when(hint2_alt == hint2 ~ "equal", TRUE ~ "different")
         ) 

# how many groups get different hints?
pct_diff_hints1_g <- round(prop.table(table(simulations$hints1_different))[2]*100,2)
num_diff_hints1_g <- table(simulations$hints1_different)[2]/4
pct_diff_hints2_g <- round(prop.table(table(simulations$hints2_different))[2]*100,2)
num_diff_hints2_g <- table(simulations$hints2_different)[2]/4

# how many individuals get different hints? 
num_diff_hints1_i <- table(simulations$hint1_diff)[1]
num_diff_hints2_i <- table(simulations$hint2_diff)[1]

######################## calculate new allocations and new payoffs

choices1 <- simulations %>% 
  mutate(allocation_alt1 = case_when(rank1 < 3 ~ decision1_alt, TRUE ~ ""), # if you are a high rank you always get what you choose)
         allocation_alt2 = case_when(rank2 < 3 ~ decision2_alt, TRUE ~ ""),
         allocation_original_follow1 = case_when(rank1 < 3 ~ decision1_original_follow, TRUE ~ ""),
         allocation_original_follow2 = case_when(rank2 < 3 ~ decision2_original_follow, TRUE ~ ""),
         allocation_alt1_follow = case_when(rank1 < 3 ~ decision1_alt_, TRUE ~ ""),
         allocation_alt2_follow = case_when(rank2 < 3 ~ decision2_alt_, TRUE ~ "")
                  )

# 524 in all allocations still missing (all the individuals with rank 3 or 4)

choices2 <- choices1 %>% 
  group_by(group_id) %>% 
  summarise(string = paste0(paste0(allocation_alt1, collapse = '')),
            string02 = paste0(paste0(allocation_alt2, collapse = '')),
            string03 = paste0(paste0(allocation_original_follow1, collapse = '')),
            string04 = paste0(paste0(allocation_original_follow2, collapse = '')),
            string05 = paste0(paste0(allocation_alt1_follow, collapse = '')),
            string06 = paste0(paste0(allocation_alt2_follow, collapse = '')),
            
            ) %>% 
  ungroup() 

choices3 <- merge(choices1, choices2, by = "group_id") 
rm(choices1, choices2)
choices3 <- choices3 %>% 
  mutate(allocation_alt1 = case_when(rank1 > 2 & string == "AA" ~ "B",
                                         rank1 > 2 & string == "BB" ~ "A",
                                         rank1 == 3 & string == "AB" ~ decision1_alt,
                                         rank1 == 3 & string == "BA" ~ decision1_alt,
                                         TRUE ~ allocation_alt1),
         allocation_alt2 = case_when(rank2 > 2 & string02 == "AA" ~ "B",
                                         rank2 > 2 & string02 == "BB" ~ "A",
                                         rank2 == 3 & string02 == "AB" ~ decision2_alt,
                                         rank2 == 3 & string02 == "BA" ~ decision2_alt,
                                         TRUE ~ allocation_alt2),
         allocation_original_follow1 = case_when(rank1 > 2 & string03 == "AA" ~ "B",
                                     rank1 > 2 & string03 == "BB" ~ "A",
                                     rank1 == 3 & string03 == "AB" ~ decision1_original_follow,
                                     rank1 == 3 & string03 == "BA" ~ decision1_original_follow,
                                     TRUE ~ allocation_original_follow1),
         allocation_original_follow2 = case_when(rank2 > 2 & string04 == "AA" ~ "B",
                                     rank2 > 2 & string04 == "BB" ~ "A",
                                     rank2 == 3 & string04 == "AB" ~ decision2_original_follow,
                                     rank2 == 3 & string04 == "BA" ~ decision2_original_follow,
                                     TRUE ~ allocation_original_follow2),
         allocation_alt1_follow = case_when(rank1 > 2 & string05 == "AA" ~ "B",
                                         rank1 > 2 & string05 == "BB" ~ "A",
                                         rank1 == 3 & string05 == "AB" ~ decision1_alt_,
                                         rank1 == 3 & string05 == "BA" ~ decision1_alt_,
                                         TRUE ~ allocation_alt1_follow),
         allocation_alt2_follow = case_when(rank2 > 2 & string06 == "AA" ~ "B",
                                         rank2 > 2 & string06 == "BB" ~ "A",
                                         rank2 == 3 & string06 == "AB" ~ decision2_alt_,
                                         rank2 == 3 & string06== "BA" ~ decision2_alt_,
                                         TRUE ~ allocation_alt2_follow),
         )


# Note: There are cases where this simulation randomly picks among a set of hints. 
# This may artificially increase the number of hint changes, since the original hint may be among the hints considered.

choices4 <- choices3 %>% 
  group_by(group_id) %>% 
  summarise(string2 = paste0(paste0(allocation_alt1, collapse = '')),
            string22 = paste0(paste0(allocation_alt2, collapse = '')),
            string3 = paste0(paste0(allocation_original_follow1, collapse = '')),
            string4 = paste0(paste0(allocation_original_follow2, collapse = '')),
            string5 = paste0(paste0(allocation_alt1_follow, collapse = '')),
            string6 = paste0(paste0(allocation_alt2_follow, collapse = '')),
            id_otree = id_otree) %>% 
  ungroup() 

# only rank 4 has no allocation here, to check run:
# table(choices4$string6, useNA = "always")
  
choices5 <- merge(choices4, choices3, by = "id_otree")  %>% 
  mutate(group_id = group_id.x) %>% 
  select(-group_id.x, -group_id.y)  

choices6 <- choices5 %>% select(id_otree, group_id, rank1, string, string2, allocation_alt1, decision1_alt, valuationsb1,
                                rank2, string02, string22, allocation_alt2, decision2_alt, valuationsb2,
                                string03, string3, allocation_original_follow1, string04, string4, allocation_original_follow2,
                                string05, string5, allocation_alt1_follow, string06, string6, allocation_alt2_follow) 

choices6 <- choices6 %>% 
  mutate(allocation_alt1 = case_when(rank1 == 4 & string2 == "AAB"  ~ "B",
                                         rank1 == 4 & string2 == "ABA"  ~ "B",
                                         rank1 == 4 & string2 == "BAA"  ~ "B",
                                         rank1 == 4 & string2 == "BBA"  ~ "A",
                                         rank1 == 4 & string2 == "BAB"  ~ "A",
                                         rank1 == 4 & string2 == "ABB"  ~ "A",
                                         TRUE ~ allocation_alt1),
         payoff_alt1 = case_when(allocation_alt1 == "A" ~ 300,
                                     allocation_alt1 == "B" ~ as.numeric(valuationsb1)),
         allocation_alt2 = case_when(rank2 == 4 & string22 == "AAB"  ~ "B",
                                         rank2 == 4 & string22 == "ABA"  ~ "B",
                                         rank2 == 4 & string22 == "BAA"  ~ "B",
                                         rank2 == 4 & string22 == "BBA"  ~ "A",
                                         rank2 == 4 & string22 == "BAB"  ~ "A",
                                         rank2 == 4 & string22 == "ABB"  ~ "A",
                                         TRUE ~ allocation_alt2),
         payoff_alt2 = case_when(allocation_alt2 == "A" ~ 300,
                                     allocation_alt2 == "B" ~ as.numeric(valuationsb2)),
         allocation_original_follow1 = case_when(rank1 == 4 & string3 == "AAB"  ~ "B",
                                     rank1 == 4 & string3 == "ABA"  ~ "B",
                                     rank1 == 4 & string3 == "BAA"  ~ "B",
                                     rank1 == 4 & string3 == "BBA"  ~ "A",
                                     rank1 == 4 & string3 == "BAB"  ~ "A",
                                     rank1 == 4 & string3 == "ABB"  ~ "A",
                                     TRUE ~ allocation_original_follow1),
         payoff_ori1 = case_when(allocation_original_follow1 == "A" ~ 300,
                                 allocation_original_follow1 == "B" ~ as.numeric(valuationsb1)),
         allocation_original_follow2 = case_when(rank2 == 4 & string4 == "AAB"  ~ "B",
                                     rank2 == 4 & string4 == "ABA"  ~ "B",
                                     rank2 == 4 & string4 == "BAA"  ~ "B",
                                     rank2 == 4 & string4 == "BBA"  ~ "A",
                                     rank2 == 4 & string4 == "BAB"  ~ "A",
                                     rank2 == 4 & string4 == "ABB"  ~ "A",
                                     TRUE ~ allocation_original_follow2),
         payoff_ori2 = case_when(allocation_original_follow2 == "A" ~ 300,
                                 allocation_original_follow2 == "B" ~ as.numeric(valuationsb2)),
         allocation_alt1_follow = case_when(rank1 == 4 & string5 == "AAB"  ~ "B",
                                         rank1 == 4 & string5 == "ABA"  ~ "B",
                                         rank1 == 4 & string5 == "BAA"  ~ "B",
                                         rank1 == 4 & string5 == "BBA"  ~ "A",
                                         rank1 == 4 & string5 == "BAB"  ~ "A",
                                         rank1 == 4 & string5 == "ABB"  ~ "A",
                                         TRUE ~ allocation_alt1_follow),
         payoff_alt1_follow = case_when(allocation_alt1_follow == "A" ~ 300,
                                 allocation_alt1_follow == "B" ~ as.numeric(valuationsb1)),
         allocation_alt2_follow = case_when(rank2 == 4 & string6 == "AAB"  ~ "B",
                                         rank2 == 4 & string6 == "ABA"  ~ "B",
                                         rank2 == 4 & string6 == "BAA"  ~ "B",
                                         rank2 == 4 & string6 == "BBA"  ~ "A",
                                         rank2 == 4 & string6 == "BAB"  ~ "A",
                                         rank2 == 4 & string6 == "ABB"  ~ "A",
                                         TRUE ~ allocation_alt2_follow),
         payoff_alt2_follow = case_when(allocation_alt2_follow == "A" ~ 300,
                                 allocation_alt2_follow == "B" ~ as.numeric(valuationsb2)),
         ) %>% 
  group_by(group_id) %>% 
  mutate(grouppay_alt1 = sum(payoff_alt1),
         grouppay_alt2 = sum(payoff_alt2),
         grouppay_ori1 = sum(payoff_ori1),
         grouppay_ori2 = sum(payoff_ori2),
         grouppay_altfollow1 = sum(payoff_alt1_follow),
         grouppay_altfollow2 = sum(payoff_alt2_follow)) 


simulations <- merge(choices6, df, by = "id_otree") 
rm(choices3, choices4, choices5, choices6)

# table(simulations$allocation1, simulations$allocation_alt1)
# table(simulations$allocation2, simulations$allocation_alt2)
# table(simulations$payoff1, simulations$payoff_alt1, simulations$redistribution)

# table(simulations$payoff_alt1, simulations$payoff_ori1)


simulations <- simulations %>% 
  mutate(allocation1_diff = case_when(allocation1 == allocation_alt1 ~ "equal", 
                                      TRUE ~ "different"),
         allocation2_diff = case_when(allocation2 == allocation_alt2 ~ "equal", 
                                      TRUE ~ "different"),
         # individual payoff
         payoff1_diff = case_when(payoff1 == payoff_alt1 ~ "equal", 
                                  payoff_alt1 > payoff1 ~ "increased",
                                  payoff_alt1 < payoff1 ~ "decreased",
                                  TRUE ~ "different"),
         payoff2_diff = case_when(payoff2 == payoff_alt2 ~ "equal",
                                  payoff_alt2 > payoff2 ~ "increased",
                                  payoff_alt2 < payoff2 ~ "decreased", 
                                  TRUE ~ "different"),
         payoff1_diff_follow = case_when(payoff_alt1_follow == payoff_ori1 ~ "equal", 
                                         payoff_alt1_follow > payoff_ori1 ~ "increased",
                                         payoff_alt1_follow < payoff_ori1 ~ "decreased",
                                  TRUE ~ "different"),
         payoff2_diff_follow = case_when(payoff_alt2_follow == payoff_ori2 ~ "equal", 
                                         payoff_alt2_follow > payoff_ori2 ~ "increased",
                                         payoff_alt2_follow < payoff_ori2 ~ "decreased",
                                  TRUE ~ "different"),
         # group payoff
         grouppay1_diff =case_when(group_payoff1 == grouppay_alt1 ~ "equal",
                                   group_payoff1 < grouppay_alt1 ~ "increased",
                                   group_payoff1 > grouppay_alt1 ~ "decreased"),
         grouppay2_diff = case_when(group_payoff2 == grouppay_alt2 ~ "equal", 
                                    group_payoff2 < grouppay_alt2 ~ "increased",
                                    group_payoff2 > grouppay_alt2~ "decreased"),
         allocation1_follow_diff = case_when(allocation_alt1_follow == allocation_original_follow1 ~ "equal",
                                             TRUE ~ "different"),
         allocation2_follow_diff = case_when(allocation_alt2_follow == allocation_original_follow2 ~ "equal",
                                             TRUE ~ "different"),
         grouppay1_follow_diff = case_when(grouppay_altfollow1 == grouppay_ori1 ~ "equal",
                                           grouppay_altfollow1 < grouppay_ori1 ~ "increased",
                                           grouppay_altfollow1 > grouppay_ori1 ~ "decreased"),
         grouppay2_follow_diff = case_when(grouppay_altfollow2 == grouppay_ori2 ~ "equal",
                                           grouppay_altfollow2 < grouppay_ori2 ~ "increased",
                                           grouppay_altfollow2 > grouppay_ori2 ~ "decreased"),
         )
         

# how different are the allocations?
# original decisions
pct_different_allocation1 <- round(prop.table(table(simulations$allocation1_diff))[1]*100,2)
num_different_allocation1 <- table(simulations$allocation1_diff)[1]
different_allocation1 <- paste0(num_different_allocation1, " (", pct_different_allocation1, "%)")
pct_different_allocation2 <- round(prop.table(table(simulations$allocation2_diff))[1]*100,2)
num_different_allocation2 <- table(simulations$allocation2_diff)[1]
different_allocation2 <- paste0(num_different_allocation2, " (", pct_different_allocation2, "%)")
rm(pct_different_allocation1, num_different_allocation1, pct_different_allocation2, num_different_allocation2)

# always follow
pct_diff_follow_allocation1 <- round(prop.table(table(simulations$allocation1_follow_diff))[1]*100,2)
num_diff_follow_allocation1 <- table(simulations$allocation1_follow_diff)[1]
diff_follow_allocation1 <- paste0(num_diff_follow_allocation1, " (", pct_diff_follow_allocation1, "%)")
pct_diff_follow_allocation2 <- round(prop.table(table(simulations$allocation2_follow_diff))[1]*100,2)
num_diff_follow_allocation2 <- table(simulations$allocation2_follow_diff)[1]
diff_follow_allocation2 <- paste0(num_diff_follow_allocation2, " (", pct_diff_follow_allocation2, "%)")
rm(pct_diff_follow_allocation1, num_diff_follow_allocation1, pct_diff_follow_allocation2, num_diff_follow_allocation2)


# how different is the individual payoff?
# split by altruism
# round(prop.table(table(simulations$payoff1_diff, simulations$redistribution))["increased", 1]*100, 2)
by_altruists_ori <- as.data.frame(rbind(paste0(table(simulations$payoff1_diff, simulations$redistribution)["increased",1],
                                                " (", round(prop.table(table(simulations$payoff1_diff, simulations$redistribution))["increased", 1]*100, 2), "%)"),
                                         paste0(table(simulations$payoff1_diff, simulations$redistribution)["increased",2],
                                                " (", round(prop.table(table(simulations$payoff1_diff, simulations$redistribution))["increased", 2]*100, 2), "%)"),
                                         paste0(table(simulations$payoff1_diff, simulations$redistribution)["decreased",1],
                                                " (", round(prop.table(table(simulations$payoff1_diff, simulations$redistribution))["decreased", 1]*100, 2), "%)"),
                                         paste0(table(simulations$payoff1_diff, simulations$redistribution)["decreased",2],
                                                " (", round(prop.table(table(simulations$payoff1_diff, simulations$redistribution))["increased", 1]*100, 2), "%)")))

# table(simulations$payoff2_diff, simulations$redistribution)
by_altruists_ori2 <- as.data.frame(rbind(paste0(table(simulations$payoff2_diff, simulations$redistribution)["increased",1],
                                                " (", round(prop.table(table(simulations$payoff2_diff, simulations$redistribution))["increased", 1]*100, 2), "%)"),
                                         paste0(table(simulations$payoff2_diff, simulations$redistribution)["increased",2],
                                                " (", round(prop.table(table(simulations$payoff2_diff, simulations$redistribution))["increased", 2]*100, 2), "%)"),
                                         paste0(table(simulations$payoff2_diff, simulations$redistribution)["decreased",1],
                                                " (", round(prop.table(table(simulations$payoff2_diff, simulations$redistribution))["decreased", 1]*100, 2), "%)"),
                                         paste0(table(simulations$payoff2_diff, simulations$redistribution)["decreased",2],
                                                " (", round(prop.table(table(simulations$payoff2_diff, simulations$redistribution))["increased", 1]*100, 2), "%)")))

# split by altruism
# table(simulations$payoff1_diff_follow, simulations$redistribution)
by_altruists_f <- as.data.frame(rbind(paste0(table(simulations$payoff1_diff_follow, simulations$redistribution)["increased",1],
                                             " (", round(prop.table(table(simulations$payoff1_diff_follow, simulations$redistribution))["increased", 1]*100, 2), "%)"),
                                      paste0(table(simulations$payoff1_diff_follow, simulations$redistribution)["increased",2],
                                             " (", round(prop.table(table(simulations$payoff1_diff_follow, simulations$redistribution))["increased", 2]*100, 2), "%)"),
                                      paste0(table(simulations$payoff1_diff_follow, simulations$redistribution)["decreased",1],
                                             " (", round(prop.table(table(simulations$payoff1_diff_follow, simulations$redistribution))["decreased", 1]*100, 2), "%)"),
                                      paste0(table(simulations$payoff1_diff_follow, simulations$redistribution)["decreased",2],
                                             " (", round(prop.table(table(simulations$payoff1_diff_follow, simulations$redistribution))["increased", 1]*100, 2), "%)")))
# table(simulations$payoff2_diff_follow, simulations$redistribution)
by_altruists_f2 <- as.data.frame(rbind(paste0(table(simulations$payoff2_diff_follow, simulations$redistribution)["increased",1],
                                              " (", round(prop.table(table(simulations$payoff2_diff_follow, simulations$redistribution))["increased", 1]*100, 2), "%)"),
                                       paste0(table(simulations$payoff2_diff_follow, simulations$redistribution)["increased",2],
                                              " (", round(prop.table(table(simulations$payoff2_diff_follow, simulations$redistribution))["increased", 2]*100, 2), "%)"),
                                       paste0(table(simulations$payoff2_diff_follow, simulations$redistribution)["decreased",1],
                                              " (", round(prop.table(table(simulations$payoff2_diff_follow, simulations$redistribution))["decreased", 1]*100, 2), "%)"),
                                       paste0(table(simulations$payoff2_diff_follow, simulations$redistribution)["decreased",2],
                                              " (", round(prop.table(table(simulations$payoff2_diff_follow, simulations$redistribution))["increased", 1]*100, 2), "%)")))
titles <- rbind("   increased welfare non-altruists",
                "   increased welfare altruists",
                "   decreased welfare non-altruists",
                "   decreased welfare altruists")
by_altruists <- cbind(titles, by_altruists_ori, by_altruists_f)
by_altruists2 <- cbind(titles, by_altruists_ori2, by_altruists_f2)
rm(titles, by_altruists_ori, by_altruists_ori2, by_altruists_f, by_altruists_f2)

# how different is the group payoff?
# table(simulations$grouppay1_diff) # all are equal
diff_grouppay1 <- paste0("0 (0%)")

pct_lower_grouppay2 <- round(prop.table(table(simulations$grouppay2_diff))[1]*100,1)
pct_same_grouppay2 <- round(prop.table(table(simulations$grouppay2_diff))[2]*100,1)
pct_higher_grouppay2 <- round(prop.table(table(simulations$grouppay2_diff))[3]*100,1)
pct_diff_grouppay2 <- 100 - pct_same_grouppay2
# 1 group (4 subjects) decreased, 3 groups (12 subjects) increased
diff_grouppay2 <- paste0(table(simulations$grouppay2_diff)[1] + table(simulations$grouppay2_diff)[3], " (", pct_diff_grouppay2, "%)")
rm(pct_lower_grouppay2, pct_same_grouppay2, pct_higher_grouppay2, pct_diff_grouppay2)

# always follow

# table(simulations$grouppay1_follow_diff) # all are equal
grouppay1_follow_diff <- paste0("0 (0%)")

# table(simulations$grouppay2_follow_diff) # all are equal
grouppay2_follow_diff <- paste0("0 (0%)")


# checking the cases with different group payoffs:

# simulations %>% 
#   filter(grouppay2_diff != "equal") %>% 
#   select(id_otree, group_id.x, rank2.x,  hint2, decision2, allocation2,payoff2, allocation_alt2, valuationsb2.x, payoff_alt2, redistribution, group_payoff2, grouppay_alt2) %>% 
#   arrange(group_id.x, -valuationsb2.x)

# first case
# what happens with this group is that originally all players choose A, therefore the ones with rank 1 and 2 get it, and the other two get B, one of which has 0 
# if they would have followed the hints it would have worked out
# in the new case, two hints are swapped, and the player with rank 2 now follows another hint (namely B) so there is a A spot left for the player with rank 3

# simulations %>% 
#   filter(group_id.x == "3r9qg") %>% 
#   select(hint2, decision2, payoff2, follow_hint2, allocation_alt2, decision2_alt, rank2.x, payoff_alt2, valuationsb2.x) %>% 
#   arrange(rank2.x)

# second case
# in this group originally the first three ranked players choose A, therefore player 3 does not get it, who has 0
# if they would have followed the hints it would have worked out
# in the new case, sometimes the two hints are swapped, and the player with rank 1 now follows the hint (to choose B) so there is a spot left for the player with rank 3

# simulations %>% 
#   filter(group_id.x == "5adeg") %>% 
#   select(hint2, decision2, payoff2, follow_hint2, allocation_alt2, decision2_alt, rank2.x, payoff_alt2, valuationsb2.x)  %>% 
#   arrange(rank2.x)

# third case
# in this group the player with rank 4 originally got B (which was 0)
# with the new hints the hint-following players ranked 1 and 2 now both select B, so there is an A spot left

# simulations %>% 
#   filter(group_id.x == "i4ei5") %>% 
#   select(hint2, decision2, payoff2, follow_hint2, allocation_alt2, decision2_alt, rank2.x, payoff_alt2, valuationsb2.x)  %>% 
#   arrange(rank2.x)

# last case 
# in this group it works the other way around: player with rank 1 did not follow the original hint so selects A as he did originally
# but the player with rank 2 followed the hint before, and now again, which now also leads to A,
# so the player with rank 3 now gets B, which is 0 for this person

# simulations %>% 
#   filter(group_id.x == "r4v7q") %>% 
#   select(hint2, decision2, payoff2, follow_hint2, allocation_alt2, decision2_alt, rank2.x, payoff_alt2, valuationsb2.x)  %>% 
#   arrange(rank2.x)

Under these assumptions, using our experimental data on the subjects’ social preferences, we find that 6 recommendations across 3 groups (1.15%) would change in Scenario 1 and 164 recommendations across 80 groups (30.53%) would change in Scenario 2.

Table F1: Welfare changes under alternative recommendations

# Analysis already done in in-text statistics code chunk. Here only stitching together in a table. 

r2 <- as.data.frame(cbind("groups with different welfare", diff_grouppay1, grouppay1_follow_diff))                  
r3 <- as.data.frame(cbind("individuals with different welfare", different_allocation1, diff_follow_allocation1)) 
r5 <- as.data.frame(cbind("groups with different welfare", diff_grouppay2, grouppay2_follow_diff))                  
r6 <- as.data.frame(cbind("individuals with different welfare", different_allocation2, diff_follow_allocation2)) 
names(r2) <-  c("V1", "V2", "V3")
names(r3) <-  c("V1", "V2", "V3")
names(by_altruists) <- c("V1", "V2", "V3")
names(r5) <-  c("V1", "V2", "V3")
names(r6) <-  c("V1", "V2", "V3")
names(by_altruists2) <- c("V1", "V2", "V3")

table <- rbind(r2, r3, by_altruists, r5, r6, by_altruists2)  
rm(r2, r3, r5, r6, by_altruists, by_altruists2)
names(table) <- NULL
  
kable(table, align=c("lcc"), format = "html") %>% # use "format = 'latex', booktabs = TRUE" if you want to output latex 
  add_header_above(header = c(" " = 1, "Original decisions" = 1, "Always follow" = 1)) %>% 
  pack_rows("Scenario 1", start_row = 1, end_row = 6) %>% 
  add_indent(c(3:6)) %>% 
  pack_rows("Scenario 2", start_row = 7, end_row = 12) %>% 
  add_indent(c(9:12)) %>% 
  footnote(general = "Table reports the number (%) of groups and individuals for which welfare changes under new recommendations. First column reports changes under original decisions (follow if originally followed, else original decision); second column compares welfare under the new vs original recommendations if these were always followed.")
Original decisions
Always follow
Scenario 1
groups with different welfare 0 (0%) 0 (0%)
individuals with different welfare 4 (0.38%) 6 (0.57%)
increased welfare non-altruists 2 (0.19%) 3 (0.29%)
increased welfare altruists 0 (0%) 0 (0%)
decreased welfare non-altruists 0 (0%) 0 (0%)
decreased welfare altruists 2 (0.19%) 3 (0.29%)
Scenario 2
groups with different welfare NA (98.9%) 0 (0%)
individuals with different welfare 90 (8.59%) 164 (15.65%)
increased welfare non-altruists 41 (3.91%) 77 (7.35%)
increased welfare altruists 4 (0.38%) 5 (0.48%)
decreased welfare non-altruists 7 (0.67%) 12 (1.15%)
decreased welfare altruists 38 (3.91%) 70 (7.35%)
Note:
Table reports the number (%) of groups and individuals for which welfare changes under new recommendations. First column reports changes under original decisions (follow if originally followed, else original decision); second column compares welfare under the new vs original recommendations if these were always followed.
rm(diff_grouppay1, grouppay1_follow_diff, different_allocation1, diff_follow_allocation1, diff_grouppay2, grouppay2_follow_diff, different_allocation2, diff_follow_allocation2)

Table G1: Balance table

vars <- c("clean_time", "riskaverse", "exp_demand", "bayesian_abs", "female", "work_student", "instr_check", "instr_fail", "attention_fail", "redistribution")
table_out <- CreateTableOne(
  data = subset(df, bot == 0), 
  strata = "treatment", 
  vars = vars,
  factorVars = c("female", "work_student", "instr_check", "instr_fail", "attention_fail", "redistribution"),
  test = TRUE)

table_out_stats <- print(table_out, varLabels = TRUE, format = "p", catDigits = 2, printToggle = FALSE)
table_out_stats <- gsub(" \\(.*\\)", "", table_out_stats)
table_out_stats <- as.data.frame(table_out_stats)[2:11,]

table_out_stats <- cbind(table_out_stats[,1], table_out_stats[,3], table_out_stats[,2], table_out_stats[,4])
rownames(table_out_stats) <- c("Completion time (min)", "Risk aversion", "Experimenter demand", "Bayesian deviation",
                               "% female", "% student status", "% instructions check", "% instructions failure", "% attention failure", "% altruist")

kable(table_out_stats, align = ("cccc")) %>% 
  pack_rows(" ", start_row = 1, end_row = 10) %>% 
  add_header_above(c(" " = 1, "Full Info" = 1, "Partial Info" = 1, "No Info" = 1, "p-value" = 1)) 
Full Info
Partial Info
No Info
p-value
Completion time (min) 13.26 14.20 14.14 0.023
Risk aversion 58.14 59.22 60.70 0.141
Experimenter demand 11.73 11.65 13.52 0.224
Bayesian deviation 13.00 13.50 14.82 0.165
% female 42.28 43.12 39.61 0.603
% student status 6.13 9.10 7.06 0.053
% instructions check 10.41 8.13 10.20 0.219
% instructions failure 60.57 65.76 63.92 0.068
% attention failure 47.75 49.73 50.59 0.601
% altruist 28.06 29.90 28.78 0.679

Table G2: Optimal choices across Scenarios (sequential logit)

df.logit <- df %>% 
  mutate(optimal_choice_logit = case_when(
           follow_hint1 == 1 & decision2 == "A" ~ 3,
           follow_hint1 == 0 & decision2 == "B" ~ 4,
           follow_hint1 == 0 & decision2 == "A" ~ 2,
           follow_hint1 == 1 & decision2 == "B" ~ 1,
           bot == 1 ~ NA_real_),
         hint1A = case_when(hint1 == "A"~ 1, 
                            bot == 1 ~ NA_real_,
                            TRUE ~ 0),
         hint2A = case_when(hint2 == "A"~ 1,
                            bot == 1 ~ NA_real_,
                            TRUE ~ 0)) %>% 
  filter(treatment == "Partial Info")

col1 <- df.logit %>% 
  filter(bot == 0 & !is.na(bayesian)) %>% 
  select(decision1) %>% 
  table() %>% 
  as.data.frame() %>% 
  rename(rowname = "decision1")

col1_1 <- df.logit %>% 
  filter(bot == 0 & !is.na(bayesian) & follow_hint1 == 1) %>% 
  select(follow_hint1) %>% 
  mutate(follow_hint1 = case_when(follow_hint1 == 1 ~ "chose optimally")) %>% 
  table() %>% 
  as.data.frame() %>% 
  rename(rowname = "follow_hint1")

col1 <- rbind(col1, col1_1) %>% 
  column_to_rownames()

rownames(col1) <- c("N (chose A)", "N (chose B)", "N (chose optimally)")

col2 <- df.logit %>% 
  filter(bot == 0 & !is.na(bayesian)) %>% 
  select(optimal_choice_logit) %>% 
  table() %>% 
  as.data.frame()

col23 <- data.frame(col2 = c(col2[3,2], col2[1,2], col2[3,2]),
                   col3=c(col2[2,2], col2[4,2], col2[2,2]))
rownames(col23) <- c("N (chose A)", "N (chose B)", "N (chose optimally)")
cols <- cbind(col1, col23)
rm(col1_1, col2, col1, col23)

# Filter the data
filtered_data <- subset(df.logit, bot == 0) %>% 
    mutate(recommendedA = hint1A)

# Fit the logistic regression model using `glm`
model1 <- glm(follow_hint1 ~ recommendedA + riskaverse * bayesian_abs + attention_fail + instr_fail + redistribution + exp_demand, data = filtered_data, family = binomial())

# Compute cluster-robust standard errors using `coeftest` and `vcovCL` functions
cluster_se <- sqrt(diag(vcovCL(model1, cluster = ~ id_otree, type = "HC1")))
coef_tests1 <- coeftest(model1, vcov = sandwich)[,4]
# compute odds ratios by exponentiating
odds_ratios <- exp(coef(model1))

# model before model1

model1b <- glm(follow_hint1 ~ recommendedA, data = filtered_data, family = binomial())

# Compute cluster-robust standard errors using `coeftest` and `vcovCL` functions
cluster_se1b <- sqrt(diag(vcovCL(model1b, cluster = ~ id_otree, type = "HC1")))
coef_tests1b <- coeftest(model1b, vcov = sandwich)[,4]
# compute odds ratios by exponentiating
odds_ratios1b <- exp(coef(model1b))

# sequence to get coefficients conditional on choices in scenario 2
filtered_data2 <- filtered_data %>% 
  mutate(recommendedA = hint2A) %>% 
  filter(follow_hint1 == 1)

model2 <- glm(choseA2 ~ recommendedA + riskaverse * bayesian_abs + attention_fail + instr_fail + redistribution + exp_demand, data = filtered_data2, family = binomial())

odds_ratios2 <- exp(coef(model2))
cluster_se2 <- sqrt(diag(vcovCL(model2, cluster = ~ id_otree, type = "HC1")))
coef_tests2 <- coeftest(model2, vcov = sandwich)[,4]

model2a <- glm(choseA2 ~ recommendedA + riskaverse * bayesian_abs + attention_fail + redistribution + exp_demand +check1 +check2 + check3 + check4, data = filtered_data2, family = binomial())

odds_ratios2a <- exp(coef(model2a))
cluster_se2a <- sqrt(diag(vcovCL(model2a, cluster = ~ id_otree, type = "HC1")))
coef_tests2a <- coeftest(model2a, vcov = sandwich)[,4]



filtered_data3 <- filtered_data %>% 
  mutate(recommendedA = hint2A) %>% 
  filter(follow_hint1 == 0)

model3 <- glm(choseA2 ~ recommendedA + riskaverse * bayesian_abs + attention_fail + instr_fail + redistribution + exp_demand, data = filtered_data3, family = binomial())

odds_ratios3 <- exp(coef(model3))
cluster_se3 <- sqrt(diag(vcovCL(model3, cluster = ~ id_otree, type = "HC1")))
coef_tests3 <- coeftest(model3, vcov = sandwich)[,4]

model3a <- glm(choseA2 ~ recommendedA + riskaverse * bayesian_abs + attention_fail + redistribution + exp_demand +check1 + check2 + check3 + check4, data = filtered_data3, family = binomial())

odds_ratios3a <- exp(coef(model3a))
cluster_se3a <- sqrt(diag(vcovCL(model3a, cluster = ~ id_otree, type = "HC1")))
coef_tests3a <- coeftest(model3a, vcov = sandwich)[,4]


# in the next code chunk I'm running this, but not including it (otherwise it gets knitted twice)

# stargazer_output <- stargazer(model1, model2, model3,  
#           coef =list(odds_ratios, odds_ratios2, odds_ratios3), 
#           p = list(coef_tests1, coef_tests2, coef_tests3), 
#           se = list(cluster_se, cluster_se2, cluster_se3), 
#           type = 'text',
#           omit.table.layout = 'sn',
#           order = c(1:3,8,4,5,6,7), # to get the interaction term to the right place
#           star.cutoffs = c(0.05, 0.01, 0.001),
#           header = FALSE
#           ) 
output_df <- data.frame(Variable = character(), Coefficient_1 = character(), Coefficient_2 = character(), 
                        Coefficient_3 = character(), Coefficient_4 =character(), 
                        Coefficient_5 = character(), Coefficient_6 =character(),
                        stringsAsFactors = FALSE)

# Extract variable names and coefficients using regular expressions
for (i in 1:length(stargazer_output)) {
  line <- stargazer_output[i]
  
  # Skip the top rows with "Dependent variable"
  if (grepl("Dependent variable:", line)) {
    next
  }
  
  # Check if the line matches the pattern for variable names and coefficients
  if (grepl("^[^=]+\\s+[^\\s]+\\s+[^\\s]+\\s+[^\\s]+", line)) {
    # Extract variable name and coefficients using regex
    line_parts <- trimws(strsplit(line, "\\s+")[[1]])
    var_name <- line_parts[1]
    coefficient_1 <- line_parts[2]
    coefficient_2 <- line_parts[3]
    coefficient_3 <- line_parts[4]
    coefficient_4 <- line_parts[5]
    coefficient_5 <- line_parts[6]
    coefficient_6 <- line_parts[7]
    
    
    # Append variable name and coefficients to the dataframe
    output_df <- rbind(output_df, data.frame(Variable = var_name, Coefficient_1 = coefficient_1, Coefficient_2 = coefficient_2, 
                                             Coefficient_3 = coefficient_3, Coefficient_4 = coefficient_4, 
                                             Coefficient_5 = coefficient_5, Coefficient_6 = coefficient_6,
                                             stringsAsFactors = FALSE))
  }
}

output_df <- output_df %>% 
  filter(!is.na(Coefficient_2)) 

output_df <- output_df[2:26,] # to remove the constant
output_df[4:17,3:7] <- output_df[4:17,2:6] # make column 1 almost empty
output_df[17:25,5] <- output_df[17:25,2]
output_df[17:25,7] <- output_df[17:25,3]
output_df[4:25,2] <- c(" ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", " ", "")
output_df[18:25,3] <- c("", "", "", "", "", "", "", "")
output_df[18:25,4] <- c("", "", "", "", "", "", "", "")
output_df[18:25,6] <- c("", "", "", "", "", "", "", "")

output_df[12:13,6] <- output_df[12:13,5]
output_df[12:13,7] <- c("", "")
output_df[12:13,5] <- c("", "")

names(output_df) <- c("var", "col1", "col2", "col3", "col4", "col5", "col6") 

output_df <- output_df %>% 
  mutate(var = case_when(var == "recommendedA" ~ "Recommended A", 
                         var == "riskaverse" ~ "Risk averse",
                         var == "riskaverse:bayesian_abs" ~ "Risk averse X non-Bayesian",
                         var == "bayesian_abs" ~ "Non Bayesian",
                         var == "attention_fail" ~ "Inattention",
                         var == "instr_fail" ~ "# attempts comprehension Qs",
                         var == "redistribution" ~ "Altruist",
                         var == "exp_demand" ~ "Experimenter demand",
                         var == "check1" ~ "Failed comprehension Q1",
                         var == "check2" ~ "Failed comprehension Q2",
                         var == "check3" ~ "Failed comprehension Q3",
                         var == "check4" ~ "Failed comprehension Q4",
                         TRUE ~ var))
names(output_df) <- c(" ", "col1", "col2", "col3", "col4", "col5", "col6") 
# 
# total <- data.frame(" " = "N (total)", 
#                     col1 = 869,
#                     col2 = 869,
#                     col3 = 869,
#                     col4 = 869,
#                     col5 = 869, 
#                     col6 = 869)
# names(total) <- c(" ", "col1", "col2", "col3", "col4", "col5", "col6") 

tabled2 <- output_df
# tabled2 <- tabled2[2:26,]
rownames(tabled2) <- NULL
names(tabled2) <- NULL

kable(tabled2, align = c("lcccccc")) %>%  # use "format = 'latex', booktabs = TRUE" if you want to output latex 
  add_header_above(c(" " = 1,  "Scenario 1" = 2, "Scenario 2" = 4)) %>% 
  # pack_rows(" ", start_row = 25, end_row = 25) %>% 
  footnote(general = "The table shows the odds ratios for making the theoretically optimal choices in the Partial Info treatment. Optimal strategies follow the Recommendation in the first Scenario, and select A in the second Scenario. Only subjects who completed all ancillary tasks are included. Simple logit estimation. Robust standard errors are in parentheses. Clustering on individual level. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001")
Scenario 1
Scenario 2
Recommended A 1.795*** 1.697*** 31.075*** 33.454*** 3.131* 2.909*
(0.151) (0.159) (0.347) (0.362) (0.460) (0.465)
Risk averse 1.005 1.021* 1.020 1.007 1.005
(0.006) (0.010) (0.010) (0.014) (0.014)
Non Bayesian 1.025 1.067* 1.067 0.975 0.969
(0.016) (0.033) (0.034) (0.058) (0.058)
Risk averse X non-Bayesian 1.000 0.999 0.999 1.000 1.000
(0.0003) (0.001) (0.001) (0.001) (0.001)
Inattention 0.949 0.941 0.987 0.680 0.692
(0.159) (0.228) (0.233) (0.422) (0.405)
# attempts comprehension Qs 1.328 0.542* 0.632
(0.163) (0.256) (0.466)
Altruist 1.013 0.729 0.718 0.406* 0.423
(0.173) (0.243) (0.245) (0.440) (0.481)
Experimenter demand 1.006 1.007 1.006 0.978 0.979
(0.005) (0.007) (0.005) (0.017) (0.005)
Failed comprehension Q1 1.129 0.840
(0.190) (0.258)
Failed comprehension Q2 0.846 1.355
(0.222) (0.525)
Failed comprehension Q3 0.650** 0.701
(0.162) (0.342)
Failed comprehension Q4 0.642* 0.698
(0.226) (0.513)
Note:
The table shows the odds ratios for making the theoretically optimal choices in the Partial Info treatment. Optimal strategies follow the Recommendation in the first Scenario, and select A in the second Scenario. Only subjects who completed all ancillary tasks are included. Simple logit estimation. Robust standard errors are in parentheses. Clustering on individual level. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001
# kable(tabled2, align = c("lcccccc"), format = 'latex', booktabs = TRUE) %>%  
#   add_header_above(c(" " = 1,  "Scenario 1" = 2, "Scenario 2" = 4)) %>% 
#   # pack_rows(" ", start_row = 25, end_row = 25) %>% 
#   footnote(general = "The table shows the odds ratios for making the theoretically optimal choices in the Partial Info treatment. Optimal strategies follow the Recommendation in the first Scenario, and select A in the second Scenario. Only subjects who completed all ancillary tasks are included. Simple logit estimation. Robust standard errors are in parentheses. Clustering on individual level. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001")
# 
# kable(newtable, align = c("lccc"), format = 'latex', booktabs = TRUE) %>%
#   add_header_above(c(" "= 1, "At least wrong once" =1, "Chose A in Scenario 1" = 1, "Chose A in Scenario 2" = 1)) %>% 
#   pack_rows("overall in Partial Info", start_row = 5, end_row = 5)

Table G3: Treatment effects excluding bots

# dropping all groups that contain any bots
dropped_groups <- df %>% 
  group_by(group_id) %>% 
  summarise(bots = paste0(paste0(bot, collapse = ''))) %>% 
  filter(bots == "0000")

dropped_groups <- merge(dropped_groups, maingroups, by = "group_id", all.x = TRUE)

library(DescTools)
dropped_groups <- dropped_groups %>% 
  mutate(treatment_trend = ordered(treatment, levels = c("Partial Info", "Full Info", "No Info")),
         treatment_trend2 = ordered(treatment, levels = c("Full Info", "Partial Info", "No Info"))) %>% 
  filter(bot == 0) 

set.seed(1234)
jt1 <- JonckheereTerpstraTest(dropped_groups$group_payoff1, dropped_groups$treatment_trend, 
                       alternative = "decreasing", nperm = 100000, exact = TRUE)$p.value
set.seed(1234)
jt2 <- JonckheereTerpstraTest(dropped_groups$group_payoff2, dropped_groups$treatment_trend2, 
                       alternative = "decreasing", nperm = 100000, exact = TRUE)$p.value

mww1.1 <- wilcox.test(group_payoff1 ~ treatment, data = subset(dropped_groups, treatment %in% c("Full Info", "Partial Info")))$p.value
mww1.2 <- wilcox.test(group_payoff1 ~ treatment, data = subset(dropped_groups, treatment %in% c("Partial Info", "No Info")))$p.value
mww1.3 <- wilcox.test(group_payoff1 ~ treatment, data = subset(dropped_groups, treatment %in% c("Full Info", "No Info")))$p.value
mww2.1 <- wilcox.test(group_payoff2 ~ treatment, data = subset(dropped_groups, treatment %in% c("Full Info", "Partial Info")))$p.value
mww2.2 <- wilcox.test(group_payoff2 ~ treatment, data = subset(dropped_groups, treatment %in% c("Partial Info", "No Info")))$p.value
mww2.3 <- wilcox.test(group_payoff2 ~ treatment, data = subset(dropped_groups, treatment %in% c("Full Info", "No Info")))$p.value


p_values <- c(mww1.1, mww1.2, mww1.3, mww2.1, mww2.2, mww2.3)

bky_adjusted <- round(fdr_sharpened_qvalues(p_values),3)
table3 <- data.frame(jt = c(jt1, jt2),
                     mww1 = c(mww1.1, mww2.1),
                     mww2 = c(mww1.2, mww2.2),
                     mww3 = c(mww1.3, mww2.3)) %>% 
  mutate_all(add_stars)

table3top <- table3[1,]
table3bottom <- table3[2,]

# add sharpened q-values
# paste square brackets
brackets <- function(x){
  ifelse(x == "", paste0(x), paste0("[", x, "]"))
}

q1 <- as.data.frame(cbind("", t(bky_adjusted[1:3]))) %>% 
  mutate_all(brackets)
q2 <- as.data.frame(cbind("", t(bky_adjusted[4:6])))%>% 
  mutate_all(brackets)

names(q1) <- names(table3)
names(q2) <- names(table3)

table3 <- rbind(table3top, q1, table3bottom, q2, c(395, 348, 229, 213))
names(table3) <- NULL
rownames(table3) <- NULL

table3[1,1] <- "0.040*" # coded by hand because this is very tedious to fix computerized...

kable(table3, align=c("cccc")) %>% 
  add_header_above(header = c("All comparisons" = 1, "Full vs. Partial" = 1, "Partial vs. No" = 1, "Full vs. No" = 1)) %>% 
  add_header_above(header = c("(Jonckheere-Terpstra)" = 1, "(Mann-Whitney-U)" = 3)) %>% 
  pack_rows("H1: Partial > Full > No", start_row = 1, end_row = 3) %>% 
  pack_rows("H2: Full > Partial > No", start_row = 3, end_row = 5) %>% 
  pack_rows("N(groups)", start_row = 5, end_row = 5) %>% 
  footnote(general = "The first column lists p-values from the Jonckheere-Terpstra trend test for the ordered aggregate social welfare levels, and columns 2-4 list p-values for two-sided pairwise comparisons using the Mann-Whitney-U test. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the six pairwise tests (Anderson, 2008) are in brackets.")
(Jonckheere-Terpstra)
(Mann-Whitney-U)
All comparisons
Full vs. Partial
Partial vs. No
Full vs. No
H1: Partial > Full > No
0.040* 0.499 0.001*** 0.000***
[0.091] [0.002] [0.002]
H2: Full > Partial > No
0.000*** 0.029* 0.059 0.002**
[0.023] [0.037] [0.003]
N(groups)
395 348 229 213
Note:
The first column lists p-values from the Jonckheere-Terpstra trend test for the ordered aggregate social welfare levels, and columns 2-4 list p-values for two-sided pairwise comparisons using the Mann-Whitney-U test. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the six pairwise tests (Anderson, 2008) are in brackets.

Table G4: First best excluding bots

# H1: W* = Partial
firstbest1 <- dropped_groups %>% 
  filter(treatment == "Partial Info") 
firstbest1_p <- round(wilcox.test(firstbest1$deviationwelfare1, mu = 0)$p.value,3)
firstbest1_n <- nrow(firstbest1)

# H2: W* = Full
firstbest2 <- dropped_groups %>% 
  filter(treatment == "Full Info") 
firstbest2_n <- nrow(firstbest2)
firstbest2_p <- round(wilcox.test(firstbest2$deviationwelfare2, mu = 0)$p.value,3)

# add stars
tabled4 <- data.frame(scenario1 = c(firstbest1_p, firstbest1_n),
                     scenario2 = c(firstbest2_p, firstbest2_n)) %>% 
  mutate_all(add_stars)

# add sharpened q-values
bky_adjusted <- round(fdr_sharpened_qvalues(c(firstbest1_p, firstbest2_p)),3)
q1 <- as.data.frame(cbind(t(bky_adjusted[1:2]))) %>% 
  mutate_all(brackets)

names(q1) <- names(tabled4)
tabled4top <- tabled4[1,]
tabled4bottom <- tabled4[2,]

tabled4 <- rbind(tabled4top, q1, tabled4bottom)
names(tabled4) <- NULL
rownames(tabled4) <- NULL

kable(tabled4, align=c("cccc")) %>% 
  add_header_above(header = c("H1: W* = Partial" = 1, "H2: W* > Full" = 1)) %>% 
  add_header_above(header = c("Scenario 1" = 1, "Scenario 2" = 1)) %>% 
  pack_rows("N(groups)", start_row = 3, end_row = 3) %>% 
  footnote(general = "The table lists the p-values for the one-sample two-sided Wilcoxon sign rank test, comparing the treatment that was
hypothesized to equal (Scenario 1) or fail to reach (Scenario 2) the aggregate social welfare optimum to the theoretical first best. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the two tests (Anderson, 2008) are in brackets.") 
Scenario 1
Scenario 2
H1: W* = Partial
H2: W* > Full
0.000*** 0.000***
[0.001] [0.001]
N(groups)
182 166
Note:
The table lists the p-values for the one-sample two-sided Wilcoxon sign rank test, comparing the treatment that was
hypothesized to equal (Scenario 1) or fail to reach (Scenario 2) the aggregate social welfare optimum to the theoretical first best. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the two tests (Anderson, 2008) are in brackets.

Table G5: Treatment effects excluding groups with multiple comprehension mistakes

# dropping all groups that contain any people who failed attention check at least twice
dropped_groups <- df %>% 
  mutate(instr_fail_total = check1 + check2 + check3 + check4,
         instr_fail2 = case_when(instr_fail_total < 2 ~ 0, 
                                TRUE ~ 1)) %>%  # failed check at least twice
  group_by(group_id) %>% 
  summarise(instr_fail2 = paste0(paste0(instr_fail2, collapse = ''))) %>% 
  filter(instr_fail2 == "0000")

dropped_groups <- merge(dropped_groups, maingroups, by = "group_id", all.x = TRUE)

library(DescTools)
dropped_groups <- dropped_groups %>% 
  mutate(treatment_trend = ordered(treatment, levels = c("Partial Info", "Full Info", "No Info")),
         treatment_trend2 = ordered(treatment, levels = c("Full Info", "Partial Info", "No Info"))) 

set.seed(1234)
jt1 <- JonckheereTerpstraTest(dropped_groups$group_payoff1, dropped_groups$treatment_trend, 
                       alternative = "decreasing", nperm = 100000, exact = TRUE)$p.value
set.seed(1234)
jt2 <- JonckheereTerpstraTest(dropped_groups$group_payoff2, dropped_groups$treatment_trend2, 
                       alternative = "decreasing", nperm = 100000, exact = TRUE)$p.value

mww1.1 <- wilcox.test(group_payoff1 ~ treatment, data = subset(dropped_groups, treatment %in% c("Full Info", "Partial Info")))$p.value
mww1.2 <- wilcox.test(group_payoff1 ~ treatment, data = subset(dropped_groups, treatment %in% c("Partial Info", "No Info")))$p.value
mww1.3 <- wilcox.test(group_payoff1 ~ treatment, data = subset(dropped_groups, treatment %in% c("Full Info", "No Info")))$p.value
mww2.1 <- wilcox.test(group_payoff2 ~ treatment, data = subset(dropped_groups, treatment %in% c("Full Info", "Partial Info")))$p.value
mww2.2 <- wilcox.test(group_payoff2 ~ treatment, data = subset(dropped_groups, treatment %in% c("Partial Info", "No Info")))$p.value
mww2.3 <- wilcox.test(group_payoff2 ~ treatment, data = subset(dropped_groups, treatment %in% c("Full Info", "No Info")))$p.value


p_values <- c(mww1.1, mww1.2, mww1.3, mww2.1, mww2.2, mww2.3)

bky_adjusted <- round(fdr_sharpened_qvalues(p_values),3)
table3 <- data.frame(jt = c(jt1, jt2),
                     mww1 = c(mww1.1, mww2.1),
                     mww2 = c(mww1.2, mww2.2),
                     mww3 = c(mww1.3, mww2.3)) %>% 
  mutate_all(add_stars)

table3top <- table3[1,]
table3bottom <- table3[2,]

# add sharpened q-values
# paste square brackets
brackets <- function(x){
  ifelse(x == "", paste0(x), paste0("[", x, "]"))
}

q1 <- as.data.frame(cbind("", t(bky_adjusted[1:3]))) %>% 
  mutate_all(brackets)
q2 <- as.data.frame(cbind("", t(bky_adjusted[4:6])))%>% 
  mutate_all(brackets)

names(q1) <- names(table3)
names(q2) <- names(table3)

table3 <- rbind(table3top, q1, table3bottom, q2, c(116, 106, 59, 67))
names(table3) <- NULL
rownames(table3) <- NULL

table3[3,1] <- "0.020*" # coded by hand because this is very tedious to fix computerized...

kable(table3, align=c("cccc")) %>% 
  add_header_above(header = c("All comparisons" = 1, "Full vs. Partial" = 1, "Partial vs. No" = 1, "Full vs. No" = 1)) %>% 
  add_header_above(header = c("(Jonckheere-Terpstra)" = 1, "(Mann-Whitney-U)" = 3)) %>% 
  pack_rows("H1: Partial > Full > No", start_row = 1, end_row = 3) %>% 
  pack_rows("H2: Full > Partial > No", start_row = 3, end_row = 5) %>% 
  pack_rows("N(groups)", start_row = 5, end_row = 5) %>% 
  footnote(general = "The first column lists p-values from the Jonckheere-Terpstra trend test for the ordered aggregate social welfare levels, and columns 2-4 list p-values for two-sided pairwise comparisons using the Mann-Whitney-U test. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the six pairwise tests (Anderson, 2008) are in brackets.")
(Jonckheere-Terpstra)
(Mann-Whitney-U)
All comparisons
Full vs. Partial
Partial vs. No
Full vs. No
H1: Partial > Full > No
0.038* 0.353 0.041* 0.132
[0.308] [0.215] [0.215]
H2: Full > Partial > No
0.020* 0.178 0.205 0.059
[0.215] [0.215] [0.215]
N(groups)
116 106 59 67
Note:
The first column lists p-values from the Jonckheere-Terpstra trend test for the ordered aggregate social welfare levels, and columns 2-4 list p-values for two-sided pairwise comparisons using the Mann-Whitney-U test. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the six pairwise tests (Anderson, 2008) are in brackets.

Table G6: First best excluding groups with multiple comprehension mistakes

# H1: W* = Partial
firstbest1 <- dropped_groups %>% 
  filter(treatment == "Partial Info") 
firstbest1_p <- round(wilcox.test(firstbest1$deviationwelfare1, mu = 0)$p.value,3)
firstbest1_n <- nrow(firstbest1)

# H2: W* = Full
firstbest2 <- dropped_groups %>% 
  filter(treatment == "Full Info") 
firstbest2_n <- nrow(firstbest2)
firstbest2_p <- round(wilcox.test(firstbest2$deviationwelfare2, mu = 0)$p.value,3)

# add stars
tabled4 <- data.frame(scenario1 = c(firstbest1_p, firstbest1_n),
                     scenario2 = c(firstbest2_p, firstbest2_n)) %>% 
  mutate_all(add_stars)

# add sharpened q-values
bky_adjusted <- round(fdr_sharpened_qvalues(c(firstbest1_p, firstbest2_p)),3)
q1 <- as.data.frame(cbind(t(bky_adjusted[1:2]))) %>% 
  mutate_all(brackets)

names(q1) <- names(tabled4)
tabled4top <- tabled4[1,]
tabled4bottom <- tabled4[2,]

tabled4 <- rbind(tabled4top, q1, tabled4bottom)
names(tabled4) <- NULL
rownames(tabled4) <- NULL

kable(tabled4, align=c("cccc")) %>% 
  add_header_above(header = c("H1: W* = Partial" = 1, "H2: W* > Full" = 1)) %>% 
  add_header_above(header = c("Scenario 1" = 1, "Scenario 2" = 1)) %>% 
  pack_rows("N(groups)", start_row = 3, end_row = 3) %>% 
  footnote(general = "The table lists the p-values for the one-sample two-sided Wilcoxon sign rank test, comparing the treatment that was
hypothesized to equal (Scenario 1) or fail to reach (Scenario 2) the aggregate social welfare optimum to the theoretical first best. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the two tests (Anderson, 2008) are in brackets.") 
Scenario 1
Scenario 2
H1: W* = Partial
H2: W* > Full
0.000*** 0.000***
[0.001] [0.001]
N(groups)
49 57
Note:
The table lists the p-values for the one-sample two-sided Wilcoxon sign rank test, comparing the treatment that was
hypothesized to equal (Scenario 1) or fail to reach (Scenario 2) the aggregate social welfare optimum to the theoretical first best. * p-val< 0.05, ** p-val < 0.01, *** p-val < 0.001 Sharpened false discovery rate q-values for the two tests (Anderson, 2008) are in brackets.

Figure G1: Treatment effects excluding groups with multiple comprehension mistakes

# dropping all groups that contain any people who failed attention check at least twice
dropped_groups <- df %>% 
    mutate(instr_fail_total = check1 + check2 + check3 + check4,
         instr_fail2 = case_when(instr_fail_total < 2 ~ 0, 
                                TRUE ~ 1)) %>%  # failed check at least twice
  group_by(group_id) %>% 
  summarise(instr_fail2 = paste0(paste0(instr_fail2, collapse = ''))) %>% 
  filter(instr_fail2 == "0000") 

df_summary <- merge(dropped_groups, df, by = "group_id", all.x = TRUE) %>%  
  group_by(treatment) %>% 
  summarize(mean_payoff1 = mean(group_payoff1), 
            sd_payoff1 = sd(group_payoff1), 
            n = n()) %>% 
  mutate(high = mean_payoff1 + qt(0.975, n-1)*(sd_payoff1/sqrt(n)), 
         low = mean_payoff1 - qt(0.975, n-1)*(sd_payoff1/sqrt(n)))

left <- ggplot(df_summary, aes(x = treatment, y = mean_payoff1, fill = treatment)) +
  geom_col(position = position_dodge()) +
  geom_errorbar(aes(ymin = low, ymax = high), 
                width = 0.2, position = position_dodge(0.9)) +
  geom_hline(yintercept=1348.48, linetype = 2, size = 1) +
  xlab("Treatment") + ylab("Mean Payoff") +
  scale_fill_manual(values = c("#A0A0A0", "#D9D9D9", "#505050"), guide = FALSE) +
  theme_bw() +
  scale_x_discrete(limits = c("No Info", "Partial Info", "Full Info")) +
  scale_y_continuous(breaks = seq(1000, 1350, by = 50)) +
  coord_cartesian(ylim = c(1000, 1380)) + labs(y = "", x = "\n (a) Scenario 1", title = "Aggregate Social Welfare")

df_summary <- merge(dropped_groups, df, by = "group_id", all.x = TRUE) %>%  
  group_by(treatment) %>% 
  summarize(mean_payoff2 = mean(group_payoff2), 
            sd_payoff2 = sd(group_payoff2), 
            n = n()) %>% 
  mutate(high = mean_payoff2 + qt(0.975, n-1)*(sd_payoff2/sqrt(n)), 
         low = mean_payoff2 - qt(0.975, n-1)*(sd_payoff2/sqrt(n)))

right <- ggplot(df_summary, aes(x = treatment, y = mean_payoff2, fill = treatment)) +
  geom_col(position = position_dodge()) +
  geom_errorbar(aes(ymin = low, ymax = high), 
                width = 0.2, position = position_dodge(0.9)) +
  geom_hline(yintercept = 1148.48, linetype = 2, size = 1) +
  xlab("Treatment") + ylab("Mean Payoff") +
  scale_fill_manual(values = c("#A0A0A0", "#D9D9D9", "#505050"), guide = FALSE) +
  theme_bw() +
  scale_x_discrete(limits = c("No Info", "Partial Info", "Full Info")) +
  scale_y_continuous(breaks = seq(800, 1150, by = 50)) +
  coord_cartesian(ylim = c(800, 1180)) + labs(y = "", x = "\n (b) Scenario 2", title = "Aggregate Social Welfare")


g <- grid.arrange(left, right, nrow=1)
The figure plots the aggregate social welfare reached in Scenarios 1 and 2 with 95% confidence intervals. The dashed line indicates the theoretical maximum aggregate welfare (first best) that can be achieved in expectation.

The figure plots the aggregate social welfare reached in Scenarios 1 and 2 with 95% confidence intervals. The dashed line indicates the theoretical maximum aggregate welfare (first best) that can be achieved in expectation.

ggsave(here("output", "figures", "treatment_effects_robustness.pdf"), g, device = 'pdf',
       width = 20, height = 10, units = "cm")
rm(df_summary, left, right)