## Not run:
#
# ##Load Blackwell data
#
# data(Blackwell)
#
# form1<-"d.gone.neg ~ d.gone.neg.l1 + d.gone.neg.l2 + d.neg.frac.l3 + camp.length + camp.length +
# deminc + base.poll + year.2002 + year.2004 + year.2006 + base.und + office"
#
#
# ##Fitting the models in Imai and Ratkovic (2014)
# ##Warning: may take a few mintues; setting time.vary to FALSE
# ##Results in a quicker fit but with poorer balance
#
# fit1<-CBMSM(formula = form1, time=Blackwell$time,id=Blackwell$demName,data=Blackwell, type="MSM",
# iterations = NULL, twostep = TRUE, msm.variance = "full", time.vary = TRUE)
#
# fit2<-CBMSM(formula = form1, time=Blackwell$time,id=Blackwell$demName,data=Blackwell, type="MSM",
# iterations = NULL, twostep = TRUE, msm.variance = "approx", time.vary = TRUE)
#
#
# ##Assessing balance
#
# bal1<-balance.CBMSM(fit1)
# bal2<-balance.CBMSM(fit2)
#
# ##Effect estimation: Replicating Effect Estimates in
# ##Table 3 of Imai and Ratkovic (2014)
#
# lm1<-lm(demprcnt[time==1]~fit1$treat.hist,data=Blackwell,weights=fit1$glm.weights)
# lm2<-lm(demprcnt[time==1]~fit1$treat.hist,data=Blackwell,weights=fit1$weights)
# lm3<-lm(demprcnt[time==1]~fit1$treat.hist,data=Blackwell,weights=fit2$weights)
#
# lm4<-lm(demprcnt[time==1]~fit1$treat.cum,data=Blackwell,weights=fit1$glm.weights)
# lm5<-lm(demprcnt[time==1]~fit1$treat.cum,data=Blackwell,weights=fit1$weights)
# lm6<-lm(demprcnt[time==1]~fit1$treat.cum,data=Blackwell,weights=fit2$weights)
#
#
#
# ### Example: Multiple Binary Treatments Administered at the Same Time
# n<-200
# k<-4
# set.seed(1040)
# X1<-cbind(1,matrix(rnorm(n*k),ncol=k))
#
# betas.1<-betas.2<-betas.3<-c(2,4,4,-4,3)/5
# probs.1<-probs.2<-probs.3<-(1+exp(-X1 %*% betas.1))^-1
#
# treat.1<-rbinom(n=length(probs.1),size=1,probs.1)
# treat.2<-rbinom(n=length(probs.2),size=1,probs.2)
# treat.3<-rbinom(n=length(probs.3),size=1,probs.3)
# treat<-c(treat.1,treat.2,treat.3)
# X<-rbind(X1,X1,X1)
# time<-c(rep(1,nrow(X1)),rep(2,nrow(X1)),rep(3,nrow(X1)))
# id<-c(rep(1:nrow(X1),3))
# y<-cbind(treat.1,treat.2,treat.3) %*% c(2,2,2)+X1 %*% c(-2,8,7,6,2) + rnorm(n,sd=5)
#
# multibin1<-CBMSM(treat~X,id=id,time=time,type="MultiBin",twostep=TRUE)
# summary(lm(y~-1+treat.1+treat.2+treat.3+X1, weights=multibin1$w))
# ## End(Not run)
Run the code above in your browser using DataLab