#######################################################
## Replication data for Bell and Quek, Authoritarian 
## Public Opinion and the Democratic Peace
#######################################################

# Set working directory to location of data file
# e.g., setwd("/Users/markbell/Documents/My Documents/Minnesota/Work/Research/Democratic Peace/IO_final_version")


#Load data
data<-read.csv("Bell_Quek_data_final.csv")

#Load library (install gplots library using install.packages() if necessary)
library(gplots)

#Get basic treatment effects
reg<-lm(approve~democracy,data=data); summary(reg)
reg<-lm(approve~trade,data=data); summary(reg)
reg<-lm(approve~alliance,data=data); summary(reg)

##########FIGURES

###### Main Paper Figure  1 (Treatment effect compared to Tomz and Weeks, Johns and Davies)
one<-mean(data[data$democracy==1,]$approve, na.rm=T)   ##Mean support when dem=1
one2<-mean(data[data$democracy==0,]$approve, na.rm=T) ##Mean support when dem=0


matrix1<-matrix(data=c(one,one2,0.347,0.41,0.447,0.494,0.209,0.342,0.419,0.533), ncol=5, byrow=F) ##create matrix for barplot

lower1<-matrix(c(0,0,0,0,0,0,0,0,0,0), ncol=5, byrow=T) ##Set CIs to zero since not displaying them
upper1<-matrix(c(0,0,0,0,0,0,0,0,0,0), ncol=5, byrow=T)

##Create Figure 1
pdf("fig1.pdf", height=7, width=11)
barplot2(matrix1, beside=T, col=c("GRAY", "WHITE"), space=c(0,0,1.5,0,0.5,0,1.5,0,0.5,0), ylim=c(0,.6), axes=F, ylab="Percentage Support", names.arg=c("China", "United Kingdom ", "United States", "United Kingdom", "United States"), cex.names=0.8,plot.ci=T, ci.l=lower1, ci.u=upper1, ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.5, grid.lty=2)

legend("topleft", fill=c("WHITE","GRAY"), legend=c("Target non-democracy","Target democracy"), cex=1, bty="n")
text(x=5.75, y=0.51, labels=c("Johns and Davies (2012)"), cex=0.8, pos=3)
lines(x=c(3.4,8.1), y=c(0.51,0.51))
lines(x=c(3.4,3.4), y=c(0.51,0.49))
lines(x=c(8.1,8.1), y=c(0.51,0.49))

text(x=11.75, y=0.549, labels=c("Tomz and Weeks (2013)"), cex=0.8, pos=3)
lines(x=c(9.4,14.1), y=c(0.549,0.549))
lines(x=c(9.4,9.4), y=c(0.549,0.529))
lines(x=c(14.1,14.1), y=c(0.549,0.529))
axis(side=2, labels=seq(0,60,10), at=seq(0,.6,0.1))
dev.off()


###### Appendix Figure 1 (All treatment effects)

one<-mean(data[data$democracy==1,]$approve, na.rm=T)   
one2<-mean(data[data$democracy==0,]$approve, na.rm=T) 
two<-mean(data[data$alliance==1,]$approve, na.rm=T)
three<-mean(data[data$trade==1,]$approve, na.rm=T)


one2<-mean(data[data$democracy==0,]$approve, na.rm=T)
two2<-mean(data[data$alliance==0,]$approve, na.rm=T)
three2<-mean(data[data$trade==0,]$approve, na.rm=T)

mat<-matrix(data=c(one,one2,two,two2,three,three2), ncol=3, byrow=F)

lower<-matrix(c(t.test(data[data$democracy==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$alliance==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[1], t.test(data[data$democracy==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[1]), ncol=3, byrow=T)

upper<-matrix(c(t.test(data[data$democracy==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$alliance==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[2], t.test(data[data$democracy==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[2]), ncol=3, byrow=T)

pdf("fig1-app.pdf", height=7, width=11)

barplot2(mat, beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.5), axes=F, ylab="Percentage Support", names.arg=c("Democracy", "Alliance", "Trade"), cex.names=1.2,plot.ci=T, ci.l=lower, ci.u=upper, ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
legend(x=1, y=0.5, fill=c("WHITE","GRAY"), legend=c("Non-democracy","Democracy"), cex=1, bty="n")
legend(x=4,y=0.5, fill=c("WHITE","GRAY"), legend=c("Non-ally","Ally"), cex=1, bty="n")
legend(x=7,y=0.5, fill=c("WHITE","GRAY"), legend=c("Low trade","High trade"), cex=1, bty="n")
text(x=c(1.5,2.5), y=0.25, labels=c("32.9 %","38.8 %"), cex=1.1)
text(x=c(4.5,5.5), y=0.25, labels=c("33.3 %","38.5 %"), cex=1.1)
text(x=c(7.5,8.5), y=0.25, labels=c("35.6 %","36.1 %"), cex=1.1)
axis(side=2, labels=seq(0,50,10), at=seq(0,.5,0.1))

dev.off()

###### Appendix Table 1 (Comparison with national census)

#Gender
data$sex[which(data$sex==-99)]<-NA
data$male<-ifelse(data$sex==1,1,0)
1-(sum(data$male,na.rm=T)/length(data$male))

#Area
data$area<-ifelse((data$region==3 | data$region==30 | data$region==11 | data$region==28 | data$region==17),1,0)
data$area<-ifelse((data$region==21 | data$region==20 | data$region==12), 2, data$area)
data$area<-ifelse((data$region==27 | data$region==18 | data$region==34 | data$region==2 | data$region==5 | data$region==19 | data$region==26),3, data$area)
data$area<-ifelse((data$region==13 | data$region==15 | data$region==16 | data$region==7 | data$region==8 | data$region==10),4, data$area)
data$area<-ifelse((data$region==4 | data$region==29 | data$region==9 | data$region==33 | data$region==31), 5, data$area)
data$area<-ifelse((data$region==25 | data$region==6 | data$region==24 | data$region==23 | data$region==32), 6, data$area)

length(na.omit(data$area[data$area==1]))/length((data$area))
length(na.omit(data$area[data$area==2]))/length((data$area))
length(na.omit(data$area[data$area==3]))/length((data$area))
length(na.omit(data$area[data$area==4]))/length((data$area))
length(na.omit(data$area[data$area==5]))/length((data$area))
length(na.omit(data$area[data$area==6]))/length((data$area))

#HH Income
data$hhold[which(data$hhold==-99)]<-NA
data$pcincome<-as.numeric(data$income2)/data$hhold
mean(data$pcincome, na.rm=T)

#Ethnicity
data$han<-ifelse(data$ethnic==1,1,0)
(sum(data$han,na.rm=T)/length(na.omit(data$han)))


####Appendix Figure 2 (Balance Checks)

pdf("fig2-app.pdf", height=11, width=11)

par(mfrow=c(3,3))

data$yob2<-2003-data$yob
plot(density(data[data$democracy==1,]$yob2, na.rm=T), xlab="Year of birth", ylim=c(0,0.04), xlim=c(1930,2010), main="Year of birth", lty=1, cex=1.2)
par(new=TRUE)
plot(density(data[data$democracy==0,]$yob2, na.rm=T), xlab="", ylim=c(0,0.04), xlim=c(1930,2010), main="", axes=F, ylab="", lty=2)
legend("topleft", legend=c("Target democracy", "Target non-democracy"), lty=c(1,2), cex=1, bty="n")

data$religion1[which(data$religion1==-99)]<-NA
data$religion_bin<-ifelse(data$religion1>0,1,0)
plot(density(data[data$democracy==1,]$religion_bin, na.rm=T), xlab="Religious (0 = No, 1 = Yes)", ylim=c(0,4), main="Religion", lty=1, cex=1.2)
par(new=TRUE)
plot(density(data[data$democracy==0,]$religion_bin, na.rm=T), xlab="", ylim=c(0,4), main="", axes=F, ylab="", lty=2)
legend("topright", legend=c("Target democracy", "Target non-democracy"), lty=c(1,2), cex=1, bty="n")

data$income3<-data$income2/1000
plot(density(data[data$democracy==1,]$income3, na.rm=T), xlab="Household income (Thousand Yuan)", xlim=c(-50,500), main="Household income", lty=1, cex=1.2)
par(new=TRUE)
plot(density(data[data$democracy==0,]$income3, na.rm=T), axes=F, xlab="", xlim=c(-50,500), main="", lty=2, ylab="", cex=1.2)
legend("topright", legend=c("Target democracy", "Target non-democracy"), lty=c(1,2), cex=1, bty="n")

plot(density(data[data$democracy==1,]$male, na.rm=T), xlab="Gender (1 = male)", ylim=c(0,2.5), main="Gender", lty=1, cex=1.2)
par(new=TRUE)
plot(density(data[data$democracy==0,]$male, na.rm=T), xlab="", ylim=c(0,2.5), main="", axes=F, ylab="", lty=2)
legend("topleft", legend=c("Target democracy", "Target non-democracy"), lty=c(1,2), cex=1, bty="n")

plot(density(data[data$democracy==1,]$area, na.rm=T), xlab="Region (1 = North, 2 = Northeast, etc)", ylim=c(0,0.5), main="Region", lty=1, cex=1.2)
par(new=TRUE)
plot(density(data[data$democracy==0,]$area, na.rm=T), xlab="", ylim=c(0,0.5), main="", axes=F, ylab="", lty=2)
legend("topleft", legend=c("Target democracy", "Target non-democracy"), lty=c(1,2), cex=1, bty="n")

data$party[which(data$party==-99)]<-NA
data$party_bin<-ifelse(data$party>1,1,0)
plot(density(data[data$democracy==1,]$party_bin, na.rm=T), xlab="Party member (1 = member)", ylim=c(0,4), main="Communist Party membership", lty=1, cex=1.2)
par(new=TRUE)
plot(density(data[data$democracy==0,]$party_bin, na.rm=T), xlab="", ylim=c(0,4), main="", axes=F, ylab="", lty=2)
legend("topleft", legend=c("Target democracy", "Target non-democracy"), lty=c(1,2), cex=1, bty="n")

plot(density(data[data$democracy==1,]$han, na.rm=T), xlab="Ethnicity (1 = Han, 0 = Other)", ylim=c(0,8), main="Han ethnicity", lty=1, cex=1.2)
par(new=TRUE)
plot(density(data[data$democracy==0,]$han, na.rm=T), xlab="", ylim=c(0,8), main="", axes=F, ylab="", lty=2)
legend("topleft", legend=c("Target democracy", "Target non-democracy"), lty=c(1,2), cex=1, bty="n")

data$marital[which(data$marital==-99)]<-NA
data$marital_bin<-ifelse(data$marital==2,1,0)
plot(density(data[data$democracy==1,]$marital_bin, na.rm=T), xlab="Marital Status (1 = Married) ", ylim=c(0,2.6), main="Marital status", lty=1, cex=1.2)
par(new=TRUE)
plot(density(data[data$democracy==0,]$marital_bin, na.rm=T), xlab="", ylim=c(0,2.6), main="", axes=F, ylab="", lty=2)
legend("topleft", legend=c("Target democracy", "Target non-democracy"), lty=c(1,2), cex=1, bty="n")

plot(density(data[data$democracy==1,]$hhold, na.rm=T), xlab="Household size", ylim=c(0,1.3), xlim=c(0,10), main="Household size", lty=1, cex=1.2)
par(new=TRUE)
plot(density(data[data$democracy==0,]$hhold, na.rm=T), xlab="", ylim=c(0,1.3), xlim=c(0,10), main="", axes=F, ylab="", lty=2)
legend("topright", legend=c("Target democracy", "Target non-democracy"), lty=c(1,2), cex=1, bty="n")
dev.off()


####Appendix Figure 3 (Treatment effects for those more likely to have a western conception of democracy)

data2<-subset(data, data$global1==1)
reg<-lm(approve~democracy,data=data2); summary(reg) #-5.8

data2<-subset(data, data$global2==1)
reg<-lm(approve~democracy,data=data2); summary(reg) #Not significant

data2<-subset(data, data$global3==1)
reg<-lm(approve~democracy,data=data2); summary(reg) #-4

data2<-subset(data, data$edu1>4)
reg<-lm(approve~democracy,data=data2); summary(reg) #-5.7


one<-mean(data[data$global1==1 & data$democracy==1,]$approve, na.rm=T)
two<-mean(data[data$global2==1 & data$democracy==1,]$approve, na.rm=T)
three<-mean(data[data$global3==1 & data$democracy==1,]$approve, na.rm=T)
four<-mean(data[(data$worldnew==4 | data$worldnew==3) & data$democracy==1,]$approve, na.rm=T)

one2<-mean(data[data$global1==1 & data$democracy==0,]$approve, na.rm=T)
two2<-mean(data[data$global2==1 & data$democracy==0,]$approve, na.rm=T)
three2<-mean(data[data$global3==1 & data$democracy==0,]$approve, na.rm=T)
four2<-mean(data[(data$worldnew==4 | data$worldnew==3) & data$democracy==0,]$approve, na.rm=T)


mat<-matrix(data=c(one,one2,two,two2,three,three2,four,four2), ncol=4, byrow=F)
lower<-matrix(c(t.test(data[data$global1==1 & data$democracy==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$global2==1 & data$democracy==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$global3==1 & data$democracy==1,]$approve, na.rm=T)$conf.int[1], t.test(data[(data$worldnew==4 | data$worldnew==3) & data$democracy==1,]$approve, na.rm=T)$conf.int[1], t.test(data[data$global1==1 & data$democracy==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$global2==1 & data$democracy==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$global3==1 & data$democracy==0,]$approve, na.rm=T)$conf.int[1], t.test(data[(data$worldnew==4 | data$worldnew==3) & data$democracy==0,]$approve, na.rm=T)$conf.int[1]), ncol=4, byrow=T)
upper<-matrix(c(t.test(data[data$global1==1 & data$democracy==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$global2==1 & data$democracy==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$global3==1 & data$democracy==1,]$approve, na.rm=T)$conf.int[2], t.test(data[(data$worldnew==4 | data$worldnew==3) & data$democracy==1,]$approve, na.rm=T)$conf.int[2], t.test(data[data$global1==1 & data$democracy==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$global2==1 & data$democracy==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$global3==1 & data$democracy==0,]$approve, na.rm=T)$conf.int[2], t.test(data[(data$worldnew==4 | data$worldnew==3) & data$democracy==0,]$approve, na.rm=T)$conf.int[2]), ncol=4, byrow=T)

pdf("fig3-app.pdf", height=7, width=14)

barplot2(mat, beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.5), axes=F, ylab="Percentage Support", names.arg=c("Speak foreign language", "Been abroad in last year", "Friends who have been abroad", "Interested in international news"), cex.names=1.2,plot.ci=T, ci.l=lower, ci.u=upper, ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
legend(x=1, y=0.5, fill=c("WHITE","GRAY"), legend=c("Non-democracy","Democracy"), cex=1, bty="n")
text(x=c(1.5,2.5), y=0.25, labels=c("33.4 %","39.2 %"), cex=1.1)
text(x=c(4.5,5.5), y=0.25, labels=c("41.5 %","42.7 %"), cex=1.1)
text(x=c(7.5,8.5), y=0.25, labels=c("34.3 %","38.3 %"), cex=1.1)
text(x=c(10.5,11.5), y=0.25, labels=c("35.4 %","42.3 %"), cex=1.1)
axis(side=2, labels=seq(0,50,10), at=seq(0,.5,0.1))

dev.off()



####Appendix Figure 4 (Treatment effects for those who understood the scenario)

understanders1<-subset(data,((data$tcheck1==1 & data$democracy==1) | (data$tcheck1==2 & data$democracy==0)))
understanders2<-subset(data,data$tcheck2==3)

one<-mean(understanders1[understanders1$democracy==1,]$approve, na.rm=T)
two<-mean(understanders2[understanders2$democracy==1,]$approve, na.rm=T)

one2<-mean(understanders1[understanders1$democracy==0,]$approve, na.rm=T)
two2<-mean(understanders2[understanders2$democracy==0,]$approve, na.rm=T)

mat<-matrix(data=c(one,one2,two,two2), ncol=2, byrow=F)
lower<-matrix(c(t.test(understanders1[understanders1$democracy==1,]$approve, na.rm=T)$conf.int[1],t.test(understanders2[understanders2$democracy==1,]$approve, na.rm=T)$conf.int[1], t.test(understanders1[understanders1$democracy==0,]$approve, na.rm=T)$conf.int[1],t.test(understanders2[understanders2$democracy==0,]$approve, na.rm=T)$conf.int[1]), ncol=2, byrow=T)
upper<-matrix(c(t.test(understanders1[understanders1$democracy==1,]$approve, na.rm=T)$conf.int[2],t.test(understanders2[understanders2$democracy==1,]$approve, na.rm=T)$conf.int[2], t.test(understanders1[understanders1$democracy==0,]$approve, na.rm=T)$conf.int[2],t.test(understanders2[understanders2$democracy==0,]$approve, na.rm=T)$conf.int[2]), ncol=2, byrow=T)

pdf("fig4-app.pdf", height=7, width=8)
barplot2(mat, beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.5), axes=F, ylab="Percentage Support", names.arg=c("Passed test 1", "Passed test 2"), cex.names=1.2,plot.ci=T, ci.l=lower, ci.u=upper, ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
legend(x=1, y=0.5, fill=c("WHITE","GRAY"), legend=c("Non-democracy","Democracy"), cex=1, bty="n")
text(x=c(1.5,2.5), y=0.22, labels=c("32.2 %","38.8 %"), cex=1.1)
text(x=c(4.5,5.5), y=0.22, labels=c("29.9 %","36.3 %"), cex=1.1)
axis(side=2, labels=seq(0,50,10), at=seq(0,.5,0.1))
dev.off()

####Appendix Figure 5 (Treatment effects by education)

data$edu1[which(data$edu1==-99)]<-NA
one<-mean(data[data$edu1<5 & data$democracy==1,]$approve, na.rm=T)
two<-mean(data[data$edu1>=5 & data$democracy==1,]$approve, na.rm=T)
three<-mean(data[data$edu1>=7 & data$democracy==1,]$approve, na.rm=T)

one2<-mean(data[data$edu1<5 & data$democracy==0,]$approve, na.rm=T)
two2<-mean(data[data$edu1>=5 & data$democracy==0,]$approve, na.rm=T)
three2<-mean(data[data$edu1>=7 & data$democracy==0,]$approve, na.rm=T)


mat<-matrix(data=c(one,one2,two,two2,three,three2), ncol=3, byrow=F)
lower<-matrix(c(t.test(data[data$edu1<5 & data$democracy==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$edu1>=5 & data$democracy==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$edu1>=7 & data$democracy==1,]$approve, na.rm=T)$conf.int[1], t.test(data[data$edu1<5 & data$democracy==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$edu1>=5 & data$democracy==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$edu1>=7 & data$democracy==0,]$approve, na.rm=T)$conf.int[1]), ncol=3, byrow=T)
upper<-matrix(c(t.test(data[data$edu1<5 & data$democracy==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$edu1>=5 & data$democracy==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$edu1>=7 & data$democracy==1,]$approve, na.rm=T)$conf.int[2], t.test(data[data$edu1<5 & data$democracy==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$edu1>=5 & data$democracy==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$edu1>=7 & data$democracy==0,]$approve, na.rm=T)$conf.int[2]), ncol=3, byrow=T)

pdf("fig5-app.pdf", height=7, width=11)
barplot2(mat, beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.6), axes=F, ylab="Percentage Support", names.arg=c("High school or less", "Completed degree", "Completed advanced degree"), cex.names=1.2,plot.ci=T, ci.l=lower, ci.u=upper, ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
legend(x=1, y=0.6, fill=c("WHITE","GRAY"), legend=c("Non-democracy","Democracy"), cex=1, bty="n")
text(x=c(1.5,2.5), y=0.15, labels=c("30.8 %","37.0 %"), cex=1.1)
text(x=c(4.5,5.5), y=0.15, labels=c("33.5 %","39.3 %"), cex=1.1)
text(x=c(7.5,8.5), y=0.15, labels=c("26.0 %","43.1 %"), cex=1.1)
axis(side=2, labels=seq(0,60,10), at=seq(0,.6,0.1))
dev.off()

####Appendix Figure 6 (Treatment effects by gender)

men<-subset(data,data$male==1)
women<-subset(data,data$male==0)

one<-mean(men[men$democracy==1,]$approve, na.rm=T)
two<-mean(women[women$democracy==1,]$approve, na.rm=T)

one2<-mean(men[men$democracy==0,]$approve, na.rm=T)
two2<-mean(women[women$democracy==0,]$approve, na.rm=T)

mat<-matrix(data=c(one,one2,two,two2), ncol=2, byrow=F)
lower<-matrix(c(t.test(men[men$democracy==1,]$approve, na.rm=T)$conf.int[1],t.test(women[women$democracy==1,]$approve, na.rm=T)$conf.int[1], t.test(men[men$democracy==0,]$approve, na.rm=T)$conf.int[1],t.test(women[women$democracy==0,]$approve, na.rm=T)$conf.int[1]), ncol=2, byrow=T)
upper<-matrix(c(t.test(men[men$democracy==1,]$approve, na.rm=T)$conf.int[2],t.test(women[women$democracy==1,]$approve, na.rm=T)$conf.int[2], t.test(men[men$democracy==0,]$approve, na.rm=T)$conf.int[2],t.test(women[women$democracy==0,]$approve, na.rm=T)$conf.int[2]), ncol=2, byrow=T)

pdf("fig6-app.pdf", height=7, width=8)
barplot2(mat, beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.5), axes=F, ylab="Percentage Support", names.arg=c("Men", "Women"), cex.names=1.2,plot.ci=T, ci.l=lower, ci.u=upper, ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
legend(x=1, y=0.5, fill=c("WHITE","GRAY"), legend=c("Non-democracy","Democracy"), cex=1, bty="n")
text(x=c(1.5,2.5), y=0.22, labels=c("33.9 %","38.7 %"), cex=1.1)
text(x=c(4.5,5.5), y=0.22, labels=c("31.4 %","39.0 %"), cex=1.1)
axis(side=2, labels=seq(0,50,10), at=seq(0,.5,0.1))
dev.off()


####Appendix Figure 7 (Treatment effects by age)

one<-mean(data[data$yob2>1989 & data$democracy==1,]$approve, na.rm=T)
two<-mean(data[data$yob2>1969 & data$yob2<1990 & data$democracy==1,]$approve, na.rm=T)
three<-mean(data[data$yob2<1970 & data$democracy==1,]$approve, na.rm=T)

one2<-mean(data[data$yob2>1989 & data$democracy==0,]$approve, na.rm=T)
two2<-mean(data[data$yob2>1969 & data$yob2<1990 & data$democracy==0,]$approve, na.rm=T)
three2<-mean(data[data$yob2<1970 & data$democracy==0,]$approve, na.rm=T)

mat<-matrix(data=c(one,one2,two,two2,three,three2), ncol=3, byrow=F)
lower<-matrix(c(t.test(data[data$yob2>1989 & data$democracy==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$yob2>1969 & data$yob2<1990 & data$democracy==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$yob2<1970 & data$democracy==1,]$approve, na.rm=T)$conf.int[1], t.test(data[data$yob2>1989 & data$democracy==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$yob2>1969 & data$yob2<1990 & data$democracy==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$yob2<1970 & data$democracy==0,]$approve, na.rm=T)$conf.int[1]), ncol=3, byrow=T)
upper<-matrix(c(t.test(data[data$yob2>1989 & data$democracy==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$yob2>1969 & data$yob2<1990 & data$democracy==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$yob2<1970 & data$democracy==1,]$approve, na.rm=T)$conf.int[2], t.test(data[data$yob2>1989 & data$democracy==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$yob2>1969 & data$yob2<1990 & data$democracy==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$yob2<1970 & data$democracy==0,]$approve, na.rm=T)$conf.int[2]), ncol=3, byrow=T)

pdf("fig7-app.pdf", height=7, width=11)
barplot2(mat, beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.5), axes=F, ylab="Percentage Support", names.arg=c("Under 25", "25 - 45", "Over 45"), cex.names=1.2,plot.ci=T, ci.l=lower, ci.u=upper, ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
legend(x=1, y=0.5, fill=c("WHITE","GRAY"), legend=c("Non-democracy","Democracy"), cex=1, bty="n")
text(x=c(1.5,2.5), y=0.15, labels=c("32.5 %","38.0 %"), cex=1.1)
text(x=c(4.5,5.5), y=0.15, labels=c("34.8 %","39.3 %"), cex=1.1)
text(x=c(7.5,8.5), y=0.15, labels=c("25.9 %","38.2 %"), cex=1.1)
axis(side=2, labels=seq(0,50,10), at=seq(0,.5,0.1))
dev.off()



####Appendix Figure 8 (Treatment effects by alliance treatment)

one<-mean(data[data$democracy==1 & data$alliance==0,]$approve, na.rm=T)
two<-mean(data[data$alliance==1,]$approve, na.rm=T)
three<-mean(data[data$trade==1,]$approve, na.rm=T)

one2<-mean(data[data$democracy==0 & data$alliance==0,]$approve, na.rm=T)
two2<-mean(data[data$alliance==0,]$approve, na.rm=T)
three2<-mean(data[data$trade==0,]$approve, na.rm=T)

mat<-matrix(data=c(one,one2,two,two2,three,three2), ncol=3, byrow=F)
lower<-matrix(c(t.test(data[data$democracy==1 & data$alliance==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$alliance==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[1], t.test(data[data$democracy==0 & data$alliance==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[1]), ncol=3, byrow=T)
upper<-matrix(c(t.test(data[data$democracy==1 & data$alliance==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$alliance==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[2], t.test(data[data$democracy==0 & data$alliance==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[2]), ncol=3, byrow=T)

pdf("fig8-app.pdf", height=7, width=9)

par(mfrow=c(1,2))
barplot2(mat[,1], beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.5), axes=F, ylab="Percentage support for military action", names.arg=c("Democracy", "Non-democracy"), main="Ally = 0", cex.names=1,plot.ci=T, ci.l=lower[,1], ci.u=upper[,1], ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
# legend(x=1, y=0.5, fill=c("WHITE","GRAY"), legend=c("Non-democracy","Democracy"), cex=1, bty="n")
# legend(x=4,y=0.5, fill=c("WHITE","GRAY"), legend=c("Non-ally","Ally"), cex=1, bty="n")
# legend(x=7,y=0.5, fill=c("WHITE","GRAY"), legend=c("Low trade","High trade"), cex=1, bty="n")
text(x=c(1.5,2.5), y=0.25, labels=c("35.1 %","41.9 %"), cex=1.1)
# text(x=c(4.5,5.5), y=0.25, labels=c("33.3 %","38.5 %"), cex=1.1)
# text(x=c(7.5,8.5), y=0.25, labels=c("35.6 %","36.1 %"), cex=1.1)
axis(side=2, labels=seq(0,50,10), at=seq(0,.5,0.1))

one<-mean(data[data$democracy==1 & data$alliance==1,]$approve, na.rm=T)
two<-mean(data[data$alliance==1,]$approve, na.rm=T)
three<-mean(data[data$trade==1,]$approve, na.rm=T)
one2<-mean(data[data$democracy==0 & data$alliance==1,]$approve, na.rm=T)
two2<-mean(data[data$alliance==0,]$approve, na.rm=T)
three2<-mean(data[data$trade==0,]$approve, na.rm=T)

mat<-matrix(data=c(one,one2,two,two2,three,three2), ncol=3, byrow=F)
lower<-matrix(c(t.test(data[data$democracy==1 & data$alliance==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$alliance==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[1], t.test(data[data$democracy==0 & data$alliance==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[1]), ncol=3, byrow=T)
upper<-matrix(c(t.test(data[data$democracy==1 & data$alliance==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$alliance==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[2], t.test(data[data$democracy==0 & data$alliance==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[2]), ncol=3, byrow=T)

barplot2(mat[,1], beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.5), axes=F, ylab="Percentage support for military action", main="Ally = 1", names.arg=c("Democracy", "Non-democracy"), cex.names=1,plot.ci=T, ci.l=lower[,1], ci.u=upper[,1], ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
text(x=c(1.5,2.5), y=0.25, labels=c("30.7 %","35.8 %"), cex=1.1)
axis(side=2, labels=seq(0,50,10), at=seq(0,.5,0.1))

dev.off()

####Appendix Figure 9 (Treatment effects by trade treatment)

one<-mean(data[data$democracy==1 & data$trade==0,]$approve, na.rm=T)
two<-mean(data[data$alliance==1,]$approve, na.rm=T)
three<-mean(data[data$trade==1,]$approve, na.rm=T)

one2<-mean(data[data$democracy==0 & data$trade==0,]$approve, na.rm=T)
two2<-mean(data[data$alliance==0,]$approve, na.rm=T)
three2<-mean(data[data$trade==0,]$approve, na.rm=T)

mat<-matrix(data=c(one,one2,two,two2,three,three2), ncol=3, byrow=F)
lower<-matrix(c(t.test(data[data$democracy==1 & data$trade==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[1], t.test(data[data$democracy==0 & data$trade==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[1]), ncol=3, byrow=T)
upper<-matrix(c(t.test(data[data$democracy==1 & data$trade==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[2], t.test(data[data$democracy==0 & data$trade==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[2]), ncol=3, byrow=T)

pdf("fig9-app.pdf", height=7, width=9)

par(mfrow=c(1,2))

barplot2(mat[,1], beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.5), main="Trade = 0",axes=F, ylab="Percentage support for military action", names.arg=c("Democracy", "Non-democracy"), cex.names=1,plot.ci=T, ci.l=lower[,1], ci.u=upper[,1], ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
text(x=c(1.5,2.5), y=0.25, labels=c("32.7 %","39.5 %"), cex=1.1)
axis(side=2, labels=seq(0,50,10), at=seq(0,.5,0.1))

one<-mean(data[data$democracy==1 & data$trade==1,]$approve, na.rm=T)
two<-mean(data[data$alliance==1,]$approve, na.rm=T)
three<-mean(data[data$trade==1,]$approve, na.rm=T)

one2<-mean(data[data$democracy==0 & data$trade==1,]$approve, na.rm=T)
two2<-mean(data[data$alliance==0,]$approve, na.rm=T)
three2<-mean(data[data$trade==0,]$approve, na.rm=T)

mat<-matrix(data=c(one,one2,two,two2,three,three2), ncol=3, byrow=F)
lower<-matrix(c(t.test(data[data$democracy==1 & data$trade==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[1], t.test(data[data$democracy==0 & data$trade==1,]$approve, na.rm=T)$conf.int[1],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[1],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[1]), ncol=3, byrow=T)
upper<-matrix(c(t.test(data[data$democracy==1 & data$trade==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==1,]$approve, na.rm=T)$conf.int[2], t.test(data[data$democracy==0 & data$trade==1,]$approve, na.rm=T)$conf.int[2],t.test(data[data$alliance==0,]$approve, na.rm=T)$conf.int[2],t.test(data[data$trade==0,]$approve, na.rm=T)$conf.int[2]), ncol=3, byrow=T)

barplot2(mat[,1], beside=T, col=c("GRAY", "WHITE"), space=c(0,1), ylim=c(0,.5), axes=F, main="Trade = 1", ylab="Percentage support for military action", names.arg=c("Democracy", "Non-democracy"), cex.names=1,plot.ci=T, ci.l=lower[,1], ci.u=upper[,1], ci.lwd=0.5, lwd=0.5, plot.grid=T, grid.lwd=0.25, grid.lty=2)
text(x=c(1.5,2.5), y=0.25, labels=c("33.1 %","38.0 %"), cex=1.1)
axis(side=2, labels=seq(0,50,10), at=seq(0,.5,0.1))

dev.off()