###################################################################################
# Replication file for "Experimental evidence for a pure collaboration effect" 	#
# Figure 2													#
# author: Mary McGrath (mary.mcgrath@northwestern.edu) ORCID: 0000-0002-8437-4039 	#
# created: 10/25/2018	R version 3.3.2								#
 ###################################################################################


library(foreign)
library(sandwich)

##Experiments 1-3##
dat<-read.dta("Collab_3-5.dta")

pdf(file="Fig2.pdf", width=7.09, height=4, pointsize=7)


given<-dat$given
Collab<-dat$Collab
effdiff<-dat$effdiff*10
effdiff2<-(dat$effdiff*10)^2
effdiff3<-(dat$effdiff*10)^3
effdiff4<-(dat$effdiff*10)^4
rd<-dat$rd
leveldiff<-dat$effdiff
leveldiff[dat$eff_level=="(+) large"]<-(-25.03)
leveldiff[dat$eff_level=="(+) small"]<-(-.03)
leveldiff[dat$eff_level=="(-) small"]<-.03
leveldiff[dat$eff_level=="(-) large"]<-25.03
leveldiff[dat$rd!=3]<-NA

data<-cbind(given,Collab,effdiff, effdiff2,effdiff3,effdiff4,rd, leveldiff)

trt<-as.data.frame(data[Collab==1,])
ctl<-as.data.frame(data[Collab==0,])

z<-qnorm(1-(.05/2))

plot(jitter(ctl$given[ctl$rd!=3], factor=30)~ctl$effdiff[ctl$rd!=3], pch=183, col="purple", xlab=NA, ylab=NA)
points(jitter(trt$given[trt$rd!=3], factor=30)~trt$effdiff[trt$rd!=3], pch=183, col="red", xlab=NA, ylab=NA)

points(jitter(ctl$given[ctl$rd==3], factor=30)~jitter(ctl$leveldiff[ctl$rd==3], amount=4), pch=183, col="purple", xlab=NA, ylab=NA)
points(jitter(trt$given[trt$rd==3], factor=30)~jitter(trt$leveldiff[trt$rd==3], amount=4), pch=183, col="red", xlab=NA, ylab=NA)




Creg<-lm(given ~ effdiff + effdiff2 + as.factor(rd), data=ctl)
Crbse<-sqrt(diag(vcovHC(Creg, type="HC0")))
cCEF<-function(x){
	Creg$coeff[1] + Creg$coeff[2]*x + Creg$coeff[3]*(x^2)
	}


plot(cCEF,from=-60,to=60,col="purple", ylim=c(0,.25), lwd=2, xlab=NA, ylab=NA, add=T)



#Confidence BANDS: yhat +/- z*residual std err*sqrt[(1/n) + (x-xbar)^2 / (n-1)*sd_x ]
xbar<-mean(ctl$effdiff) ##
xSD<-sd(ctl$effdiff)

cl95<-function(x){
	Creg$coeff[1] + Creg$coeff[2]*x + Creg$coeff[3]*(x^2) - #yhat
	(z*summary(Creg)$sigma * #residual std error
		sqrt( (1/nobs(Creg)) + ( ((x-xbar)^2) / ((nobs(Creg)-1)*xSD) ) ) )
}

cu95<-function(x){
	Creg$coeff[1] + Creg$coeff[2]*x + Creg$coeff[3]*(x^2) + #yhat
	(z*summary(Creg)$sigma * #residual std error
		sqrt( (1/nobs(Creg)) + ( ((x-xbar)^2) / ((nobs(Creg)-1)*xSD) ) ) )
}


plot(cl95,from=-60,to=60,col="purple",lty=3,add=T)
plot(cu95,from=-60,to=60,col="purple",lty=3,add=T)


#COLLABORATIVELY WORKING

Treg<-lm(given ~ effdiff + effdiff2 + as.factor(rd), data=trt)
Trbse<-sqrt(diag(vcovHC(Treg, type="HC0")))
tCEF<-function(x){
	Treg$coeff[1] + Treg$coeff[2]*x + Treg$coeff[3]*(x^2)
	}



#Confidence BANDS: yhat +/- z*residual std err*sqrt[(1/n) + (x-xbar)^2 / (n-1)*sd_x ]
xbar<-mean(trt$effdiff) ##
xSD<-sd(trt$effdiff)

tl95<-function(x){
	Treg$coeff[1] + Treg$coeff[2]*x + Treg$coeff[3]*(x^2) - #yhat
	(z*summary(Treg)$sigma * #residual std error
		sqrt( (1/nobs(Treg)) + ( ((x-xbar)^2) / ((nobs(Treg)-1)*xSD) ) ) )
}

tu95<-function(x){
	Treg$coeff[1] + Treg$coeff[2]*x + Treg$coeff[3]*(x^2) + #yhat
	(z*summary(Treg)$sigma * #residual std error
		sqrt( (1/nobs(Treg)) + ( ((x-xbar)^2) / ((nobs(Treg)-1)*xSD) ) ) )
}


plot(tCEF,from=-60,to=60, lwd=2, col="red",  add=T)
plot(tl95,from=-60,to=60, col="red", lty=3,add=T)
plot(tu95,from=-60,to=60, col="red", lty=3,add=T)


legend(x=30, y=.42, legend=c("Collaboratively Working","Separately Working"), fill=c("red","purple"), cex=.75)


mtext(side=1,"Partner's Effort Level Relative to Participant's Effort Level\n (Partner's Effort Score - Your Effort Score)",  line=4)
mtext(side=2,"Amount Given to Partner (cents)", line=2.5)


dev.off()