#code for Chapter 3 "Propensity Score Weighting" of book: #Leite, W. L. (2017). Practical propensity score methods using R. #Thousand Oaks, CA: Sage. # #this is the code that was used to generate the example results in the book #As the R software and the R packages used in this example are updated frequently #some incompatibilities between the current code and new R versions or package versions #may appear #Any updates to the code will be posted at: # http://www.practicalpropensityscore.com #Estimate the effect of high school student participation in #career academies on future income #using data from the Education Longitudinal Study (ELS) #===================================================================== #Doubly-robust (DR) estimation with regression estimation of #propensity score weighted means using the propensity score as a covariate #load required package #load required survey library require(survey) options(survey.lonely.psu = "adjust") #load data load(file="Chapter_3_treatment_effect_estimates.Rdata") #define survey design for treated observations surveyDesign2006T <- subset(surveyDesign2006Boot, treat==1) #define survey design for untreated observations surveyDesign2006U <- subset(surveyDesign2006Boot, treat==0) #fit model for treated and untreated observations modelT <- svyglm(F2ERN5P2~pScores+I(pScores^2)+I(pScores^3),surveyDesign2006T) modelU <- svyglm(F2ERN5P2~pScores+I(pScores^2)+I(pScores^3),surveyDesign2006U) #obtain predicted outcomes with treated observations Yt1 <- predict(modelT, newdata=data.frame(pScores=with(ELS.data.imputed,pScores[treat==1])), vcov=TRUE, type="response") Yt0 <- predict(modelU, newdata=data.frame(pScores=with(ELS.data.imputed,pScores[treat==1])), vcov=TRUE, type="response") #obtain the difference between predicted outcomes diff <- Yt1 - Yt0 ATT2006DR <- svycontrast(diff, contrasts=rep(1/length(diff),length(diff))) #save all objects save(list=ls(), file="Chapter_3_treatment_effect_estimates.Rdata",compress=T)