.Rhistory

dorm.clusters <- unique(d1$dormid)
dorm.clusters <- unique(d$dormid)
estimate_ate <- function(df) {
tmt_clusters <- sample(c(0,1), length(dorm.clusters), replace=TRUE)
return (mean( (d1 %>% filter(dormid %in% (dorm.clusters %>% filter(tmt_clusters==1)) ))$turnout) -
mean( (d1 %>% filter(dormid %in% (dorm.clusters %>% filter(tmt_clusters==0)) ))$turnout) )
}
distribution.under.sharp.null <- replicate(5000, estimate_ate(d))
estimate_ate <- function(df) {
tmt_clusters <- sample(c(0,1), length(dorm.clusters), replace=TRUE)
return (mean( (df %>% filter(dormid %in% (dorm.clusters %>% filter(tmt_clusters==1)) ))$turnout) -
mean( (df %>% filter(dormid %in% (dorm.clusters %>% filter(tmt_clusters==0)) ))$turnout) )
}
distribution.under.sharp.null <- replicate(5000, estimate_ate(d))
dorm.clusters <- unique(d$dormid)
tmt_clusters <- sample(c(0,1), length(dorm.clusters), replace=TRUE)
tmt_clus <- dorm.clusters %>% filter(tmt_clusters==1)
head(dorm.clusters)
tmt_clus <- dorm.clusters[tmt_clusters==1,]
tmt_clus <- dorm.clusters[tmt_clusters==1]
head(tmt_clus)
dorm.clusters <- unique(d$dormid)
#tmt_clusters <- sample(c(0,1), length(dorm.clusters), replace=TRUE)
#tmt_clus <- dorm.clusters[tmt_clusters==1]
#con_clus <- dorm.clusters[tmt_clusters==0]
estimate_ate <- function(df) {
tmt_clusters <- sample(c(0,1), length(dorm.clusters), replace=TRUE)
tmt_clus <- dorm.clusters[tmt_clusters==1]
con_clus <- dorm.clusters[tmt_clusters==0]
return (mean( (df %>% filter(dormid %in% (dorm.clusters[tmt_clusters==1]) ))$turnout) -
mean( (df %>% filter(dormid %in% (dorm.clusters[tmt_clusters==0]) ))$turnout) )
}
distribution.under.sharp.null <- replicate(5000, estimate_ate(d))
hist(distribution.under.sharp.null,
main="Histogram of ITT",
freq=FALSE)
freq=FALSE)
hist(distribution.under.sharp.null, main = "Histogram of ITT", freq=FALSE)
mean(itt < distribution.under.sharp.null)
mean(ITT < distribution.under.sharp.null)
ITT
head(distribution.under.sharp.null)
CACE_1 <- ITT/ITT_d
paste0("CACE calculated using ITT and ITT_d is: ", round(CACE_1,4))
m1 <- lm(turnout ~ treat2, data=d)
coeftest(m1, vcovHC(m1))
m2 <- lm(contact ~ treat2, data=d)
summary(m2)
cace_2 <- coef(m1) / coef(m2)
cace_2
m3 <- lm(contact ~ treat2, data=d)
m4 <- lm(turnout ~ predict(m3), data=d)
summary(m4)
confint(mod4)
confint(m4)
d4 <- data.frame(treatment = c(rep(0,2572), rep(1, 2572), rep(0, 2579)),
contacted = c(rep(0,2572), rep(1, 486), rep(0, 2086), rep(1, 470), rep(0, 2109)),
placebo = c(rep(0,2572), rep(0,2572), rep(1, 2579)),
voted = c(rep(1, 803), rep(0, 1769),
rep(1, 190), rep(0, 296), rep(1, 683), rep(0, 1403),
rep(1, 140), rep(0, 330), rep(1, 678), rep(0, 1431))
)
table(d4$treatment, d4$placebo)
table(d4$contacted, d4$placebo)
treat_compliers <- nrow(d4 %>% filter(contacted==1 & treated==1))
treat_compliers <- nrow(d4 %>% filter(contacted==1 & treatment==1))
paste0("Number of compliers in treatment group: ", treat_compliers)
placebo_compliers <- nrow(d4 %>% filter(contacted==1 & placebo==1))
paste0("Number of compliers in placebo group: ", placebo_compliers)
tmt <- d4 %>% filter(contacted==1 & placebo==0)
se_tmt <- sd(tmt$contacted) / sqrt(length(tmt$contacted))
se_tmt
tmt <- d4 %>% filter(contacted==1 & treatment==1)
se_tmt <- sd(tmt$contacted) / sqrt(length(tmt$contacted))
se_tmt
treat_compliers <- nrow(d4 %>% filter(contacted==1 & treatment==1))
paste0("Number of compliers in treatment group: ", treat_compliers)
placebo_compliers <- nrow(d4 %>% filter(contacted==1 & placebo==1))
paste0("Number of compliers in placebo group: ", placebo_compliers)
ITT_D_Treatment <- treat_compliers / length(d4 %>% filter(treatment==1))
ITT_D_Treatment
tmt <- d4 %>% filter(treatment==1)
486/2572
ITT_D_Treatment <- treat_compliers / nrow(d4 %>% filter(treatment==1))
ITT_D_Treatment
placebo_compliers <- nrow(d4 %>% filter(contacted==1 & placebo==1))
paste0("Number of compliers in placebo group: ", placebo_compliers)
ITT_D_Placebo <- placebo_compliers / nrow(d4 %>% filter(placebo==1))
paste0("Proportion of compliers in treatment group: ", ITT_D_Placebo)
d4_placebo <- d4 %>% filter(placebo==1)
d4_baseline <- d4 %>% filter(treatment==0 & placebo==0)
d4_treatment <- d4 %>% filter(treatment==1)
cace_placebo <- lm(voted ~ contacted, data=d4_placebo)
summary(cace_placebo)
d4_ITT <- mean(d4_placebo$voted) - mean(d4_baseline$voted)
d4_ITT
mean((d4_placebo %>% filter(contacted==1))$voted) -
mean((d4_placebo %>% filter(contacted==0))$voted)
mean((d4_placebo %>% filter(contacted==1))$voted) -
mean((d4_baseline %>% filter(contacted==0))$voted)
mean(d4_placebo$contacted)
d4_ITT <- mean((d4_placebo %>% filter(contacted==1))$voted) -
mean((d4_placebo %>% filter(contacted==0))$voted)
d4_ITT
d4_ITT / 0.188224
0.004968802/0.18822
tmt_ITT <- mean((d4_treatment %>% filter(assigned==1))$voted) - mean(d4_baseline$voted)
tmt_ITT <- mean((d4_treatment %>% filter(contacted==1))$voted) - mean(d4_baseline$voted)
tmt_ITT_d <- mean(d4_treatment$contacted)
tmt_CACE <- tmt_ITT / tmt_ITT_d
paste0("CACE for treatment group: ", tmt_CACE)
placebo_ITT <- mean((d4_placebo %>% filter(contacted==1))$voted) - mean(d4_baseline$voted)
placebo_ITT_d <- mean(d4_placebo$contacted)
placebo_CACE <- placebo_ITT / placebo_ITT_d
paste0("CACE for placebo group: ", placebo_CACE)
tmt_turnout_rate <- .3909
placebo_turnout_rate <- .2979
CACE_2 <- tmt_turnout_rate - placebo_turnout_rate
print(CACE_2)
mean((d4_treatment %>% filter(contacted==1))$voted)
mean((d4_placebo %>% filter(contacted==1))$voted)
tmt_turnout_rate <- .3909
placebo_turnout_rate <- .2979
CACE_2 <- tmt_turnout_rate - placebo_turnout_rate
paste0("THe value of CACE calculated by comparing turnout rates among compliers in treatment and placebo groups: ", CACE_2)
paste0("CACE by comparing turnout rates among compliers in treatment and placebo groups: ", CACE_2)
# load packages
library(foreign)
library(dplyr)
library(sandwich)
library(lmtest)
library(AER)
library(multiwayvcov)
data1 = read.csv("Compliers_NeverTakers_Data.csv")
data1
ATE <- mean(data1$Treatment_Yi) - mean(data1$Control_Yi)
CACE <- mean((data1 %>% filter(IsComplier==1))$Treatment_Yi) -
mean((data1 %>% filter(IsComplier==1))$Control_Yi)
paste0("The value of ATE is: ", round(ATE, 2))
paste0("The value of CACE is: ", round(CACE, 2))
CACE_500 <- ((400/1000) - (700/2000)) / (500/1000)
paste0("If 500 subjects were actually contacted, CACE will be: ", CACE_500)
CACE_250 <- ((400/1000) - (700/2000)) / (250/1000)
paste0("If 250 subjects were actually contacted, CACE will be: ", CACE_250)
library(foreign)
d <- read.dta("./data/Guan_Green_CPS_2006.dta")
d <- na.omit(d)
head(d)
ITT <- mean((d %>% filter(treat2==1))$turnout) - mean((d %>% filter(treat2==0))$turnout)
paste0("The value of ITT calculated as differnece in 2-group means is: ", ITT)
ITT_d <- mean((d %>% filter(treat2==1))$contact)
itt_fit <- lm(turnout ~ treat2, data=d)
summary(itt_fit)
# rubust SEs
coeftest(itt_fit, vcovHC(itt_fit))
# calculating CACE using means calculation
CACE <- ITT/ITT_d
paste0("CACE calculated using ITT and ITT_d is: ", round(CACE,4))
# calculating CACE using reg
cace_fit <- ivreg(turnout ~ contact, ~treat2, data=d)
paste0("CACE calculated using regression is: ", round(cace_fit$coefficients[2],4))
# perform randomization inference
dorm.clusters <- unique(d$dormid)
estimate_ate <- function(df) {
tmt_clusters <- sample(c(0,1), length(dorm.clusters), replace=TRUE)
tmt_clus <- dorm.clusters[tmt_clusters==1]
con_clus <- dorm.clusters[tmt_clusters==0]
return (mean( (df %>% filter(dormid %in% (dorm.clusters[tmt_clusters==1]) ))$turnout) -
mean( (df %>% filter(dormid %in% (dorm.clusters[tmt_clusters==0]) ))$turnout) )
}
distribution.under.sharp.null <- replicate(5000, estimate_ate(d))
hist(distribution.under.sharp.null,
main="Histogram of ITT",
freq=FALSE)
CACE_1 <- ITT/ITT_d
paste0("CACE calculated using ITT and ITT_d is: ", round(CACE_1,4))
m1 <- lm(turnout ~ treat2, data=d)
coeftest(m1, vcovHC(m1))
m2 <- lm(contact ~ treat2, data=d)
summary(m2)
cace_2 <- coef(m1) / coef(m2)
cace_2
m3 <- lm(contact ~ treat2, data=d)
m4 <- lm(turnout ~ predict(m3), data=d)
summary(m4)
d4 <- data.frame(treatment = c(rep(0,2572), rep(1, 2572), rep(0, 2579)),
contacted = c(rep(0,2572), rep(1, 486), rep(0, 2086), rep(1, 470), rep(0, 2109)),
placebo = c(rep(0,2572), rep(0,2572), rep(1, 2579)),
voted = c(rep(1, 803), rep(0, 1769),
rep(1, 190), rep(0, 296), rep(1, 683), rep(0, 1403),
rep(1, 140), rep(0, 330), rep(1, 678), rep(0, 1431))
)
treat_compliers <- nrow(d4 %>% filter(contacted==1 & treatment==1))
paste0("Number of compliers in treatment group: ", treat_compliers)
ITT_D_Treatment <- treat_compliers / nrow(d4 %>% filter(treatment==1))
paste0("Proportion of compliers in treatment group: ", ITT_D_Treatment)
placebo_compliers <- nrow(d4 %>% filter(contacted==1 & placebo==1))
paste0("Number of compliers in placebo group: ", placebo_compliers)
ITT_D_Placebo <- placebo_compliers / nrow(d4 %>% filter(placebo==1))
paste0("Proportion of compliers in treatment group: ", ITT_D_Placebo)
# calculating CACE using reg
d4_placebo <- d4 %>% filter(placebo==1)
d4_baseline <- d4 %>% filter(treatment==0 & placebo==0)
d4_treatment <- d4 %>% filter(treatment==1)
cace_placebo <- lm(voted ~ contacted, data=d4_placebo)
summary(cace_placebo)
d4_ITT <- mean(d4_placebo$voted) - mean(d4_baseline$voted)
# calculating CACE for placebo group
placebo_ITT <- mean(d4_placebo$voted) - mean(d4_baseline$voted)
placebo_ITT_d <- mean(d4_placebo$contacted)
placebo_CACE <- placebo_ITT / placebo_ITT_d
paste0("CACE for placebo group: ", placebo_CACE)
# calculating CACE for treatment group
tmt_ITT <- mean(d4_treatment$voted) - mean(d4_baseline$voted)
tmt_ITT_d <- mean(d4_treatment$contacted)
tmt_CACE <- tmt_ITT / tmt_ITT_d
paste0("CACE for treatment group: ", tmt_CACE)
tmt_turnout_rate <- .3909
placebo_turnout_rate <- .2979
CACE_2 <- tmt_turnout_rate - placebo_turnout_rate
paste0("CACE by comparing turnout rates among compliers in treatment and placebo groups: ", CACE_2)
d5 <- read.dta("./data/Hough_WorkingPaper_2010.dta")
head(d5)
d5 <- na.omit(d5)
d5_ate <- mean((d5 %>% filter(run==1))$tetris) - mean((d5 %>% filter(run==0))$tetris)
paste0("ATE for effect of running on Tetris score: ", d5_ate)
d5_mod <- lm(tetris ~ run, data=d5)
summary(d5_mod)
paste0("ATE calculated using regression: ", d5_mod$coefficients[2], " with p-value of 0.01035")
d5$yesterday_tetris <- c(NA, d5$tetris[1:23])
ate_yesterday <- lm(yesterday_tetris ~ run, data=d5)
summary(ate_yesterday)
d5_ate <- mean((d5 %>% filter(run==1))$tetris) - mean((d5 %>% filter(run==0))$tetris)
paste0("ATE for effect of running on Tetris score: ", d5_ate)
d5_mod <- lm(tetris ~ run, data=d5)
summary(d5_mod)
paste0("ATE calculated using regression: ", d5_mod$coefficients[2], " with p-value of 0.01035")
d5_gre <- lm(gre ~ run, data=d5)
summary(d5_gre)
d5_energy <- lm(energy ~ run, data=d5)
summary(d5_energy)
mod_both <- lm(tetris ~ run + yesterday_tetris, data=d5)
summary(mod_both)
d5$yesterday_run <- c(NA, d5$run[1:23])
View(d5)
mod_both <- lm(tetris ~ run + yesterday_run, data=d5)
summary(mod_both)
mean(d4_placebo$voted) - mean(d4_baseline$voted)
mean(d4_placebo$voted)
mean(d4_treatment$voted)
mean(d4_baseline$voted)
mean(d4_treatment$voted) - mean(d4_placebo$voted)
mean(d4_treatment$voted)
mean(d4_placebo$voted)
setwd("~/GitHub/ExperimentsForTrust")
# Load the expriment's data
dmt <- read.csv("./cleaneddata/mturk.csv")
head(dmt)
summary(dmt)
str(dmt)
plot(dmt$Q6)
plot(dmt$Q7)
table(dmt$Q7)
# Load the expriment's data
dmt_lu <- read.csv("./cleaneddata/lucid.csv")
head(dmt_lu)
table(dmt_lu$Q7)
View(dmt_lu)
table(dmt_lu$Q7, dmt_lu$FL_3_DO)
table(dmt$Q7, dmt$FL_3_DO)
ggplot(dmt_lu, aes(x=FL_3_DO, fill=Q6, y=Q7)) +
geom_dotplot(binaxis="y", method = "histodot", binwidth = 5, stackdir = "center")
library(ggplot2)
geom_dotplot(binaxis="y", method = "histodot", binwidth = 5, stackdir = "center")
ggplot(dmt_lu, aes(x=FL_3_DO, fill=Q6, y=Q7))
ggplot(dmt_lu, aes(x=FL_3_DO, y=Q7))
ggplot(dmt_lu, aes(x=FL_3_DO))
ggplot(dmt_lu, aes(x=FL_3_DO)) + geom_bar
geom_bar(dmt_lu, aes(x=FL_3_DO))
p2 <- ggplot(data=dmt_lu, aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=Q6), stat = "identity") +
labs(title="Counts", y="", x = "Tmt vs. Control") +
#ylim(0, 300) +
#scale_y_continuous(labels=function(n) {format(n, scientific = FALSE)} )
#expand_limits(y=0) + #axis(labels=format(scientific=FALSE)) +
theme_classic() + theme(legend.position = "none")
p2
p2 <- ggplot(data=dmt_lu, aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=Q6), stat = "identity") +
labs(title="Counts", y="", x = "Tmt vs. Control") +
#ylim(0, 300) +
#scale_y_continuous(labels=function(n) {format(n, scientific = FALSE)} )
#expand_limits(y=0) + #axis(labels=format(scientific=FALSE)) +
theme_classic()
p2
p2 <- ggplot(data=dmt_lu, aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=Q7), stat = "identity") +
labs(title="Counts", y="", x = "Tmt vs. Control") +
#ylim(0, 300) +
#scale_y_continuous(labels=function(n) {format(n, scientific = FALSE)} )
#expand_limits(y=0) + #axis(labels=format(scientific=FALSE)) +
theme_classic()
p2
str(dmt_lu)
dmt_lu$new_Q7 <- ifelse(dmt_lu$Q7=="I do not trust what the speaker was saying.", "Trust", ifelse(dmt_lu$Q7=="I trust what the speaker was saying.", "No Trust", "NA"))
View(dmt_lu)
p2 <- ggplot(data=dmt_lu, aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=Q7), stat = "identity") +
labs(title="Counts", y="", x = "Tmt vs. Control") +
#ylim(0, 300) +
#scale_y_continuous(labels=function(n) {format(n, scientific = FALSE)} )
#expand_limits(y=0) + #axis(labels=format(scientific=FALSE)) +
theme_classic()
p2
p2 <- ggplot(data=dmt_lu, aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y="", x = "Tmt vs. Control") +
#ylim(0, 300) +
#scale_y_continuous(labels=function(n) {format(n, scientific = FALSE)} )
#expand_limits(y=0) + #axis(labels=format(scientific=FALSE)) +
theme_classic() #+ theme(legend.position = "none")
p2
p2 <- ggplot(data=dmt_lu %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y="", x = "Tmt vs. Control") +
#ylim(0, 300) +
#scale_y_continuous(labels=function(n) {format(n, scientific = FALSE)} )
#expand_limits(y=0) + #axis(labels=format(scientific=FALSE)) +
theme_classic() #+ theme(legend.position = "none")
library(dplyr)
p2 <- ggplot(data=dmt_lu %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y="", x = "Tmt vs. Control") +
#ylim(0, 300) +
#scale_y_continuous(labels=function(n) {format(n, scientific = FALSE)} )
#expand_limits(y=0) + #axis(labels=format(scientific=FALSE)) +
theme_classic() #+ theme(legend.position = "none")
p2
p2 <- ggplot(data=dmt_lu %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y=" ", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p2
dmt$new_Q7 <- ifelse(dmt$Q7=="I do not trust what the speaker was saying.", "Trust", ifelse(dmt$Q7=="I trust what the speaker was saying.", "No Trust", "NA"))
summary(dmt)
p1 <- ggplot(data=dmt_lu %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y=" ", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p1 <- ggplot(data=dmt %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y=" ", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p2
p1
library(ggthemr)
ggthemr('earth', type="outer", layout='scientific', spacing=2)
p1
p1 <- ggplot(data=dmt %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y="Total Counts", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p1
p1 <- ggplot(data=dmt %>% filter(new_Q7 != "NA")) + #, aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y="Total Counts", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p1
p1 <- ggplot(data=dmt %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
#geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y="Total Counts", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p1
p1 <- ggplot(data=dmt %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO), stat = "identity") +
labs(title="Counts", y="Total Counts", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p1
p1 <- ggplot(data=dmt %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(x=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y="Total Counts", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p1
View(dmt_lu)
p1 <- ggplot(data=dmt %>% filter(new_Q7 != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=new_Q7), stat = "identity") +
labs(title="Counts", y="Total Counts", x = "Tmt vs. Control") +
ylab("") +
theme_classic() #+ theme(legend.position = "none")
p1
dmt$Outcome <- ifelse(dmt$Q7=="I do not trust what the speaker was saying.", "Trust", ifelse(dmt$Q7=="I trust what the speaker was saying.", "No Trust", "NA"))
p1 <- ggplot(data=dmt %>% filter(Outcome != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=Outcome), stat = "identity") +
labs(title="Counts", y="Total Counts", x = "Tmt vs. Control") +
ylab("") +
theme_classic() #+ theme(legend.position = "none")
p1
dmt_lu$Outcome <- ifelse(dmt_lu$Q7=="I do not trust what the speaker was saying.", "Trust", ifelse(dmt_lu$Q7=="I trust what the speaker was saying.", "No Trust", "NA"))
p2 <- ggplot(data=dmt_lu %>% filter(Outcome != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=Outcome), stat = "identity") +
labs(title="Counts", y=" ", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p2
# Load the expriment's data
dff <- read.csv("./cleaneddata/ff.csv")
head(dff)
dff$Outcome <- ifelse(dmt_lu$Q7=="I do not trust what the speaker was saying.", "Trust", ifelse(dmt_lu$Q7=="I trust what the speaker was saying.", "No Trust", "NA"))
dff$Outcome <- ifelse(dff$Q7=="I do not trust what the speaker was saying.", "Trust", ifelse(dff$Q7=="I trust what the speaker was saying.", "No Trust", "NA"))
summary(dff)
p3 <- ggplot(data=dff %>% filter(Outcome != "NA"), aes(x=FL_3_DO)) +
geom_bar(aes(y=FL_3_DO, fill=Outcome), stat = "identity") +
labs(title="Counts", y=" ", x = "Tmt vs. Control") +
theme_classic() #+ theme(legend.position = "none")
p2
p3
# Load all data
dff <- read.csv("./cleaneddata/ff.csv")
dl <- read.csv("./cleaneddata/lucid.csv")
dmt <- read.csv("./cleaneddata/mturk.csv")
# Adding columns for control, treatment and clusters (M, F)
t = c("FTreatment", "MTreatment")
s = c("FTreatment", "FControl")
dmt$treat = ifelse(dmt$FL_3_DO %in% t, 1, 0 )
dmt$female = ifelse(dmt$FL_3_DO %in% s, 1, 0)
dmt$trust = ifelse(dmt$Q7 == "I trust what the speaker was saying.", 1,
ifelse (dmt$Q6 != "Pritikin" | dmt$Q7 == "", NA, 0))
dl$treat = ifelse(dl$FL_3_DO %in% t, 1, 0 )
dl$female = ifelse(dl$FL_3_DO %in% s, 1, 0)
dl$trust = ifelse(dl$Q7 == "I trust what the speaker was saying.", 1,
ifelse (dl$Q6 != "Pritikin" | dl$Q7 == "", NA, 0))
dff$treat = ifelse(dff$FL_3_DO %in% t, 1, 0 )
dff$female = ifelse(dff$FL_3_DO %in% s, 1, 0)
dff$trust = ifelse(dff$Q7 == "I trust what the speaker was saying.", 1,
ifelse (dff$Q6 != "Pritikin" | dff$Q7 == "", NA, 0))
dall <- rbind(dff, dl, dmt)
# Models with female covariate
lmff = lm(trust ~ treat + female * treat, data = dff)
lml = lm(trust ~ treat + female * treat, data = dl)
lmmt = lm(trust ~ treat + female * treat, data = dmt)
lmall = lm(trust ~ treat + female * treat, data = dall)
# Results from the regression with all subjects
summary(lmall)
ate = lmall$coefficients['treat']
paste0('The ATE of the dissonace video is: ', ate)
# Stargazer
library(stargazer)
stargazer(lmff, lml, lmmt, lmall, title="Experimental Results", column.labels=c("FriendsFamily","Lucid","MTurk","AllSources"))
# Randomization inference among all the subjects
# Not worrying about M & F distributions in this one
dall$treat_rnd = sample(dall$treat)
head(dall)
summary(dall)
# Generic function to randomly(pseudo) pick the treatment assignments
coeff_treat_rnd <- function(d) {
d$treat_rnd = sample(d$treat)
lmr = lm(trust ~ treat_rnd + female * treat, data=d, na.action = na.omit)
return (lmr$coefficients['treat_rnd'])
}
# Trying another randomization
coeff_treat_rnd(dall)
# 10K Randomizations
h.distribution.under.sharp.null <- replicate(10000, coeff_treat_rnd(dall))
h.ate.mean <- mean(h.distribution.under.sharp.null)
paste0("Mean fromm 1k randomizations: ", h.ate.mean)
# Graph for the estimates
plot(density(h.distribution.under.sharp.null),
main = "Density under Sharp Null")
# Adding our original ATE to the plot
abline(v = ate, col = "blue")
# Adding the mean to the plot
abline(v = mean(h.distribution.under.sharp.null), col = 'orange')
# Values close to what we get
pn = sum(abs(h.distribution.under.sharp.null) >= abs(ate))
paste0('p value: ', pn)
# 2 tailed p value
pv = mean(abs(h.distribution.under.sharp.null) >= abs(ate))
paste0('p value: ', pv)
# RI by randomizing within male and female groups
dallf = dall[which(dall$female == 1), ]
dallf$treat_rnd = sample(dallf$treat)
head(dallf)
summary(dallf)
# RI by randomizing within male and female groups
dallm = dall[which(dall$female == 0), ]
dallm$treat_rnd = sample(dallm$treat)
head(dallm)
summary(dallm)
dall2 = rbind(dallf, dallm)
head(dall2)
summary(dall2)
# 100K Randomizations
h.distribution.under.sharp.null.2 <- replicate(100000, coeff_treat_rnd(dall2))
h.ate.mean.2 <- mean(h.distribution.under.sharp.null.2)
paste0("Mean fromm 100k randomizations: ", h.ate.mean.2)
# Graph for the estimates
plot(density(h.distribution.under.sharp.null.2),
main = "Density under Sharp Null")
# Adding our original ATE to the plot
abline(v = ate, col = "blue")
# Adding the mean to the plot
abline(v = mean(h.distribution.under.sharp.null.2), col = 'orange')
# Values close to what we get
pn2 = sum(abs(h.distribution.under.sharp.null.2) >= abs(ate))
paste0('Number of observations bigger than the experiment: ', pn2)
# 2 tailed p value
pv2 = mean(abs(h.distribution.under.sharp.null.2) >= abs(ate))
paste0('p value: ', pv2)
library(ggplot2)
library(ggthemr)
install.packages("ggthemr", lib="/Library/Frameworks/R.framework/Versions/3.4/Resources/library")
install.packages("ggthemr", lib="/Library/Frameworks/R.framework/Versions/3.4/Resources/library")
install.packages("httpuv")