Supplementary Material
Notes of homogeneity conditions:
For a non-IV estimator, the non-parametric bias formulas presented by Vanderweele and Arah (2011) demonstrate that covariate balance only appears as a bias component under homogeneity assumptions regarding the effect of treatment on the outcome across levels of the unmeasured confounder. In the main text, we followed derivations by Brookhart and Schneeweiss (2007) that rely on a constant treatment effect assumption or alternative homogeneity assumptions implied by the error terms. These homogeneity assumptions lead to covariate balance being a bias component for both the IV and non-IV estimator. Without these homogeneity assumptions, covariate balance would not necessarily be a bias component for the IV estimator.
Brookhart, M. A., & Schneeweiss, S. (2007). Preference-based instrumental variable methods for the estimation of treatment effects: assessing validity and interpreting results.The International Journal of Biostatistics,3(1).
VanderWeele, T. J., & Arah, O. A. (2011). Unmeasured Confounding for General Outcomes, Treatments, and Confounders: Bias Formulas for Sensitivity Analysis.Epidemiology,22(1), 42.
R code:
library(ggplot2)
library(gridExtra)
#This code relies on the ggplot2 and gridExtra packages. To install type:
#install.packages(c("ggplot2","gridExtra"),dependencies=TRUE)
path<-"E:\\"
#######################################################
## function to graph unscaled and scaled graphs
## for single instrumental variable
##
## inputs:
## - dataset (see below for organization example)
## - location of file to save
## - name for files
## - limits on the axis (default=[-1,1])
#######################################################
graphfun_singleIV <- function(dat,location,graphname,lb=-1,ub=1) {
col1 <- dat[,1]
col2 <- dat[,2]
col3 <- as.numeric(dat[,3])
col4 <- as.numeric(dat[,4])
datun <- data.frame(col1,col2,1,col4,"(a) Unscaled")
colnames(datun) <- c("covariate","contrast","scale","difference","plot")
datsc <- data.frame(col1,col2,col3,col4,"(b) Scaled")
colnames(datsc) <- c("covariate","contrast","scale","difference","plot")
datstacked <- rbind(datun,datsc)
index <- with(datstacked,order(contrast,difference))
datstacked <- datstacked[index,]
datstacked$covariate <- factor(datstacked$covariate,
as.character(unique(datstacked$covariate)))
p <- ggplot(data=datstacked,
aes(x=covariate,y=difference*scale,shape=factor(contrast))) + geom_point(size=2.5)
p <- p + scale_shape_manual(name=" Contrast",values=c(1,17)) + ylab("Bias Component")
p <- p + xlab("Covariate") + coord_flip() + ylim(lb,ub)
p <- p + theme(legend.position="bottom",panel.margin = unit(1.25, "lines"),aspect.ratio=3/2)
plot <- p + facet_wrap(~plot)
ggsave(filename=paste(location,graphname,".png"),height=6,width=8,units="in",dpi=300,plot=plot)
list(plot=plot,plotdata=datstacked)
}
####################
##graph mcclellan 1994
####################
tx <- ' By Treatment'
iv <- 'By Instrument'
## scaling factor: 1/IV denominator
c <- 1/0.067
mcclellan <- matrix(c(
'Female',tx,1,.397-.535,
'Black',tx,1,.043-.060,
#'Rural',tx,1,.738-.696,
'Cancer',tx,1,.008-.022,
'Pulmonary disease',tx,1,.093-.111,
'Dementia',tx,1,.001-.012,
'Diabetes',tx,1,.171-.183,
'Renal disease',tx,1,.007-.023,
'Cerebrovascular disease',tx,1,.028-.054,
'Female',iv,c,.513-.495,
'Black',iv,c,.071-.043,
#'Rural',iv,c,.065-.524,
'Cancer',iv,c,.019-.019,
'Pulmonary disease',iv,c,.104-.109,
'Dementia',iv,c,.0099-.0094,
'Diabetes',iv,c,.181-.180,
'Renal disease',iv,c,.020-.019,
'Cerebrovascular disease',iv,c,.048-.048
),nrow=16,ncol=4,byrow=TRUE)
f <- graphfun_singleIV(mcclellan,location=path,"mcclellan",-.5,.5)
##bias ratios
mcclellan2 <- matrix(c(
'Female',.397-.535,.513-.495,
'Black',.043-.060,.071-.043,
'Rural',.738-.696,.065-.524,
'Cancer',.008-.022,.019-.019,
'Pulmonary disease',.093-.111,.104-.109,
'Dementia',.001-.012,.0099-.0094,
'Diabetes',.171-.183,.181-.180,
'Renal disease',.007-.023,.020-.019,
'Cerebrovascular disease',.028-.054,.048-.048),
nrow=9,ncol=3,byrow=TRUE)
biasratio <- (as.numeric(mcclellan2[,3])*c)/as.numeric(mcclellan2[,2])
biasratio
#######################################################
## function to graph unscaled and scaled graphs
## for multiple instrumental variables
##
## inputs:
## - dataset (see below for organization example)
## - location of file(s)
## - root name for files
## - limits on the axis (default=[-1,1])
#######################################################
graphfun_multIV <- function (dat,location,graphname,lb=-1,ub=1) {
col1 <- dat[,1]
col2 <- dat[,2]
col3 <- as.numeric(dat[,3])
col4 <- as.numeric(dat[,4])
datraw <- data.frame(col1,col2,col3,col4)
N <- length(unique(datraw[,2]))
col5 <- sort(rep(1:N,nrow(datraw)/N))
datraw <- data.frame(datraw,col5)
colnames(datraw) <- c("covariate","instrument","scale","difference","table")
index <- which(datraw$table!=1)
datwide <- data.frame(datraw[index,],rep(datraw[-index,"difference"],N-1))
colnames(datwide) <- c("covariate","instrument","scale","iv.difference","table","tx.difference")
datwide$table <- NULL
datwidetx <- data.frame(datwide)
datwidetx$plot.difference <- datwide$tx.difference
datwidetx$contrast <- "By Treatment"
datwideun <- data.frame(datwide)
datwideun$plot.difference <- datwide$iv.difference
datwideun$contrast <- "By Unscaled IV"
datwidesc <- data.frame(datwide)
datwidesc$plot.difference <- datwide$iv.difference*datwide$scale
datwidesc$contrast <- "By Scaled IV"
datstacked <- rbind(datwidetx,datwideun,datwidesc)
index <- with(datstacked,order(contrast,instrument,tx.difference))
datstacked <- datstacked[index,]
datstacked$covariate <- factor(datstacked$covariate,
as.character(unique(datstacked$covariate)))
datstacked$contrast <- factor(datstacked$contrast,levels=c("By Treatment","By Unscaled IV","By Scaled IV"))
label_wrap <- function(variable, value) {
lapply(strwrap(as.character(value), width=15, simplify=FALSE),
paste, collapse="\n")
}
#UNSCALED
p <- ggplot(data=datstacked[which(datstacked$contrast!="By Scaled IV"),], aes(x=covariate,y=plot.difference,shape=factor(contrast))) + ylim(lb,ub)
p <- p + geom_point(size=2) + coord_flip()
p <- p + facet_grid(.~instrument,labeller=label_wrap)
p <- p + scale_shape_manual(name="Contrast",values=c(1,17,2))
plot1 <- p + ylab("Bias Component") + xlab("Covariate") + theme(legend.position="bottom",panel.margin = unit(1.25, "lines"),aspect.ratio=3/2) + facet_grid(.~instrument,labeller=label_wrap)
#SCALED
p <- ggplot(data=datstacked[which(datstacked$contrast!="By Unscaled IV"),], aes(x=covariate,y=plot.difference,shape=factor(contrast))) + ylim(lb,ub)
p <- p + geom_point(size=2) + coord_flip()
p <- p + facet_grid(.~instrument,labeller=label_wrap)
p <- p + scale_shape_manual(name="Contrast",values=c(1,17,2))
plot2 <- p + ylab("Bias Component") + xlab("Covariate") + theme(legend.position="bottom",panel.margin = unit(1.25, "lines"),aspect.ratio=3/2) + facet_grid(.~instrument,labeller=label_wrap)
#ALL
p <- ggplot(data=datstacked, aes(x=covariate,y=plot.difference,shape=factor(contrast))) + ylim(lb,ub)
p <- p + geom_point(size=2) + coord_flip()
p <- p + facet_grid(.~instrument,labeller=label_wrap)
p <- p + scale_shape_manual(name="Contrast",values=c(1,17,2))
plot3 <- p + ylab("Bias Component") + xlab("Covariate") + theme(legend.position="bottom",panel.margin = unit(1.25, "lines"),aspect.ratio=3/2) + facet_grid(.~instrument,labeller=label_wrap)
ggsave(filename=paste(location,graphname,"unscaled.png"),height=6,width=9,units="in",dpi=300,plot=plot1)
ggsave(filename=paste(location,graphname,"scaled.png"),height=6,width=9,units="in",dpi=300,plot=plot2)
ggsave(filename=paste(location,graphname,"combined.png"),height=6,width=9,units="in",dpi=300,plot=plot3)
list(unscaled=plot1,scaled=plot2,combined=plot3,plotdata=datstacked)
}
####################
##graph fang 2010
####################
#note: post-exposure covariates are commented out
#note: all but one reference level for cat covariates are commented out
tx <- 'Thiazide Diuretic Treatment'
iv1 <- 'IV #1: Area Practice Style (DACC)'
iv2 <- 'IV #2: Area Practice Style (PCSA)'
iv3 <- 'IV #3: Physician Preference'
c1 <- 1/.188
c2 <- 1/.16
c3 <- 1/.045
fang <- matrix(c(
#treatment
#'Age 55 and below ',tx,1,-0.016,
#'Age 56-65 ',tx,1,-0.017,
#'Age 66-75 ',tx,1,0.01,
'Age 76+ ',tx,1,0.023,
'Female ',tx,1,0.084,
'Poverty 0%-50% ',tx,1,0.013,
#'Poverty 51%-100% ',tx,1,0.011,
#'Poverty 101%-150% ',tx,1,-0.023,
#'Poverty 151% and above ',tx,1,-0.002,
'Chronic kidney disease ',tx,1,-0.029,
'Diabetes ',tx,1,-0.06,
'Hyperlipidemia ',tx,1,-0.019,
'Myocardial infarction ',tx,1,-0.011,
'Coronary revascularization ',tx,1,-0.005,
'Ischemic heart disease ',tx,1,-0.09,
'Heart failure ',tx,1,-0.053,
'Atrium fibrillation ',tx,1,-0.016,
'Stroke ',tx,1,-0.004,
'Transient ischemic attack ',tx,1,-0.01,
'Statins use',tx,2,-0.023,
#'0 antihypertensive drugs',tx,3,-0.004,
#'1-2 antihypertensive drugs ',tx,4,-0.208,
'3+ antihypertensive drugs ',tx,5,0.212,
#'ACE inhibitors/ARBs ',tx,6,0.176,
#'Beta blockers ',tx,7,-0.033,
#'Calcium channel blockers',tx,8,-0.182,
#'Potassium sparing diuretics ',tx,9,-0.002,
#'Other hypertension drugs',tx,10,-0.021,
#dacc
#'Age 55 and below ',iv1,c1,-0.004,
#'Age 56-65 ',iv1,c1,-0.007,
#'Age 66-75 ',iv1,c1,-0.011,
'Age 76+ ',iv1,c1,0.023,
'Female ',iv1,c1,-0.007,
'Poverty 0%-50% ',iv1,c1,0.008,
#'Poverty 51%-100% ',iv1,c1,-0.003,
#'Poverty 101%-150% ',iv1,c1,-0.002,
#'Poverty 151% and above ',iv1,c1,-0.002,
'Chronic kidney disease ',iv1,c1,-0.001,
'Diabetes ',iv1,c1,0.012,
'Hyperlipidemia ',iv1,c1,-0.008,
'Myocardial infarction ',iv1,c1,0.003,
'Coronary revascularization ',iv1,c1,0.003,
'Ischemic heart disease ',iv1,c1,-0.005,
'Heart failure ',iv1,c1,0,
'Atrium fibrillation ',iv1,c1,0.002,
'Stroke ',iv1,c1,-0.002,
'Transient ischemic attack ',iv1,c1,-0.009,
'Statins use',iv1,c1,0,
#'0 antihypertensive drugs',iv1,c1,0,
#'1-2 antihypertensive drugs ',iv1,c1,0,
'3+ antihypertensive drugs ',iv1,c1,-0.001,
#'ACE inhibitors/ARBs ',iv1,c1,0.003,
#'Beta blockers ',iv1,c1,0.003,
#'Calcium channel blockers',iv1,c1,-0.045,
#'Potassium sparing diuretics ',iv1,c1,0.006,
#'Other hypertension drugs',iv1,c1,-0.008,
#pcsa
#'Age 55 and below ',iv2,c2,0.006,
#'Age 56-65 ',iv2,c2,0,
#'Age 66-75 ',iv2,c2,-0.008,
'Age 76+ ',iv2,c2,-0.002,
'Female ',iv2,c2,0.005,
'Poverty 0%-50% ',iv2,c2,0.009,
#'Poverty 51%-100% ',iv2,c2,-0.004,
#'Poverty 101%-150% ',iv2,c2,-0.003,
#'Poverty 151% and above ',iv2,c2,-0.002,
'Chronic kidney disease ',iv2,c2,-0.002,
'Diabetes ',iv2,c2,0.025,
'Hyperlipidemia ',iv2,c2,-0.005,
'Myocardial infarction ',iv2,c2,0.002,
'Coronary revascularization ',iv2,c2,0.003,
'Ischemic heart disease ',iv2,c2,-0.012,
'Heart failure ',iv2,c2,0.004,
'Atrium fibrillation ',iv2,c2,0.003,
'Stroke ',iv2,c2,0,
'Transient ischemic attack ',iv2,c2,-0.008,
'Statins use',iv2,c2,-0.001,
#'0 antihypertensive drugs',iv2,c2,-0.001,
#'1-2 antihypertensive drugs ',iv2,c2,-0.01,
'3+ antihypertensive drugs ',iv2,c2,0.01,
#'ACE inhibitors/ARBs ',iv2,c2,0.023,
#'Beta blockers ',iv2,c2,0.016,
#'Calcium channel blockers',iv2,c2,-0.052,
#'Potassium sparing diuretics ',iv2,c2,-0.004,
#'Other hypertension drugs',iv2,c2,-0.004,
#preference
#'Age 55 and below ',iv3,c3,0.002,
#'Age 56-65 ',iv3,c3,0.004,
#'Age 66-75 ',iv3,c3,0.003,
'Age 76+ ',iv3,c3,-0.01,
'Female ',iv3,c3,-0.007,
'Poverty 0%-50% ',iv3,c3,-0.012,
#'Poverty 51%-100% ',iv3,c3,0.008,
#'Poverty 101%-150% ',iv3,c3,0.001,
#'Poverty 151% and above ',iv3,c3,0.002,
'Chronic kidney disease ',iv3,c3,-0.006,
'Diabetes ',iv3,c3,-0.003,
'Hyperlipidemia ',iv3,c3,-0.013,
'Myocardial infarction ',iv3,c3,-0.004,
'Coronary revascularization ',iv3,c3,-0.003,
'Ischemic heart disease ',iv3,c3,-0.011,
'Heart failure ',iv3,c3,-0.009,
'Atrium fibrillation ',iv3,c3,-0.001,
'Stroke ',iv3,c3,0.001,
'Transient ischemic attack ',iv3,c3,-0.003,
'Statins use',iv3,c3,0.002,
#'0 antihypertensive drugs',iv3,c3,-0.001,
#'1-2 antihypertensive drugs ',iv3,c3,-0.01,
'3+ antihypertensive drugs ',iv3,c3,0.01#,
#'ACE inhibitors/ARBs ',iv3,c3,0.022,
#'Beta blockers ',iv3,c3,-0.004,
#'Calcium channel blockers',iv3,c3,-0.003,
#'Potassium sparing diuretics ',iv3,c3,0.002,
#'Other hypertension drugs',iv3,c3,-0.005,
),nrow=60,ncol=4,byrow=TRUE) #nrow=92
f <- graphfun_multIV(dat=fang,location=path,graphname="fang",lb=-.4,ub=.4)