treatSens: Performing Dual-Parameter Sensitivity Analysis

Description

treatSens performs dual-parameter sensitivity analysis for unobserved confounding for both continuous and binary treatment variables.

Usage

treatSens(formula, sensParam = "coef", resp.family = gaussian, trt.family = gaussian, 
theta = 0.5, grid.dim = c(8,4), standardize = TRUE, 
nsim = 20, zero.loc = 1/3, verbose = FALSE, 
buffer = 0.1, weights = NULL, data = NULL, seed = 1234, 
iter.j = 10, offset = TRUE, core = NULL, spy.range = NULL, 
spz.range = NULL, jitter = FALSE, trim.wt = 10)

Arguments

formula

an object of class formula (or one that can be coerced to that class): a symbolic description of the model to be fitted.

sensParam

a string giving the choice of the form of the sensitivity parameters. Supports "coef" for model coefficient and "cor" for partial correlations (with continuous treatment only)

resp.family

an assumed family for GLM of the response model. Currently only supports gaussian.

trt.family

an assumed family for GLM of the treatment model. The default is gaussian. binomial(link="probit") can be also specified.

theta

this option specifies the marginal probability that a binary unobserved confounder takes one (i.e Pr(U=1)). The default is 0.5.

grid.dim

the final dimensions of output grid. GLM.sens draws nsim of unobserved confounders and corresponding treatment effect estimates for each grid cell. The first argument specifies the horizontal (treatment) dimension, and the second

standardize

logical. If TRUE all variables except binary ones are standardized to have mean 0 and standard deviation 1. The default is TRUE.

nsim

number of simulated Us to generate (and hence to average over) per cell in grid. The default is 20.

zero.loc

A parameter used to automatically determine the range of each axis by specifying at what point along the y=x line the "treatment effect=0" contour will cross it. The parameter specifies the fraction of the length of the y=x line that this intersection oc

verbose

logical. If TRUE, GLM.sens will show detailed messages.

buffer

restriction to range of coefficients on U to ensure stability around the edges. The default is 0.1. This option is rarely used.

weights

the user can supply a vector of weights or specify one of three causal estimads -- "ATE", "ATT", or "ATC" -- in which case the algorithm will automatically create appropriate weights. The default is NULL

data

an optional data frame. If not found in data, the variables are taken from environment(formula), typically the environment from which glm is called. The default is NULL.

seed

random seed. The default is 1234.

iter.j

number of iterations used to draw simulated U from its conditional distribution in the iterated algorithm required for the binary treatment variable model. The default is 10. This option is used only when trt.family=binomial(link="probi

offset

Models are fit by setting zeta.z*U or zeta.y*U explicitly as an offset (rather than simply including generated U in the model). Only consistent with sensParam = "coef". The default is TRUE.

core

number of CPU cores used for parallel processing (Max=8). The default is NULL.

spy.range

custom range for the sensitiviy parameter associated with Y (the coefficient on U in the outcome model or partial correlation of U with Y given X), e.g. c(0,2). When this option is specified, zetaz.range must be also specified an

spz.range

custom range for the sensitivity parameter associated with Z (the coefficient on U in the treatment model or partial correlation of U with Z given X), e.g. c(-2,2). When this option is specified, zetay.range must be also specifie

jitter

logical. If TRUE, this option adds jitter to grids near the axis. The default is FALSE.

trim.wt

the maximum size of weight as a percentage of the sample size of the inferential group for the causal estimand. This option is used only when weights option is specified as "ATE", "ATT" or "ATC". The default is 10. For example,

Value

treatSens returns an object of class "sensitivity". An object of class "sensitivity" is a list containing the following components:
model.typeestimation methods used (i.e. "GLM""). Currently only GLM methods are available, but in future other types of models may be added.
sensParamtype of sensitivity parameter used (model coefficient or partial correlation).
tauindividual estimated treatment effects corresponding to each grid point (combination of sensitivity parameters) and each drawing of simulated U.
se.tauindividual standard error estimates of treatment effects for each grid point and each drawing of simulated U.
zeta.zindividual estimated coefficients on U in the treatment model corresponding to each grid point and each drawing of simulated U.
zeta.yindividual estimated coefficients on U in the outcome model corresponding to each grid point and each drawing of simulated U.
se.zzindividual standard error estimates of the coefficient on U in the treatment model corresponding to each grid point and each drawing of simulated U.
se.zyindividual standard error estimates of the coefficient on U in the outcome model corresponding to each grid point and each drawing of simulated U.
Youtcome variable
Ztreatment variable
Xcovariates
sig2.trtindividual residual variances for the treatment model corresponding to each grid point and each drawing of simulated U.
sig2.respindividual residual variances for the outcome model corresponding to each grid point and each drawing of simulated U.
tau0averaged estimated treatment effect assuming no unobserved confounding (i.e. naive model).
se.tau0averaged standard error estimate of the treatment effect without unobserved confounding (i.e. naive model)
Xcoefthe coefficients on the covariates in the treatment model (1st column) and the outcome model (2nd column) without unobserved confounding (i.e. naive model).
Xcoef.plotcoefficients as for Xcoef, but after transforming covariates to have positive coefficients in the response model (X = -X, used for plotting purposes).

Details

This function performs a dual-parameter sensitivity analysis for unobserved confounding by drawing simulated potential confounders U from the conditional distribution of U given observed response, treatment and covariates. Weights are calculated for any of three different estimands. The weights for the ATE are calculated as $1/\hat{e}(x)$ for the treatment group observations and $1/(1-\hat{e}(x))$ for the control group observations, where $\hat{e}(x)$ denotes the estimated modification score (that is, the estimate of the probability of being treated conditional on observed covariates, $x$. The weights for the ATT are 1 for the treatment group observations and $\hat{e}(x)/(1-\hat{e}(x))$ for the control group observations. The weights for the ATC are $(1-\hat{e}(x))/\hat{e}(x)$ for the treatment group observations and 1 for the control group observations. These weights are discussed in Carnegie et al (2014) and similar versions using the propensity score are discussed in a variety of places including Gelman and Hill, 2007.

References

Carnegie NB, Hill JH and Harada M. Assessing sensitivity to unmeasured confounding using simulated potential confounders (under review)

Examples

Run this code

N = 500          #number of observations

zetay = .5             #coefficient on U in the outcome model
zetaz = .5             #coefficient on U in the treatment model
betaz = c(.75,-.5,.25) #coefficients of X in the treatment model
betay = c(.5,1,-1.5)   #coefficients of X in the outcome model
tau = .25              #treatment effect

X = matrix(rnorm(3*N),N,3)           #covariates

###################################
# Example with continuous treatment

U = rnorm(N, 0, 1)         #unmeasured confounder
Z = rnorm(N,X %*% betaz + zetaz*U,1)         #treatment variable
Y = rnorm(N,X %*% betay + zetay*U + tau*Z,2) #outcome variable

# sensitivity analysis
out.cont <- treatSens(Y~Z+X, grid.dim = c(10,5), nsim = 5,
                     standardize = FALSE, verbose = TRUE, zero.loc = "full")

summary(out.cont) #Prints average value of tau (trt effect) for each cell in grid.

sensPlot(out.cont) # draw contour plot
# You can see when zetay=zetaz=.5, tau is about .25.

###############################
# Example with binary treatment

U = rbinom(N,1,.5)                   #unmeasured confounder
ps = pnorm(X%*%betaz + zetaz*(U-.5)) #propensity score
Z = rbinom(N,1,ps)                   #treatment variable
epsilon = rnorm(N,0,2)               #error term
Y0 = X%*%betay + zetay*(U-.5) + epsilon       #potential outcome(Z=0)
Y1 = X%*%betay + zetay*(U-.5) + tau + epsilon #potential outcome(Z=1)
Y = Y0*(1-Z) + Y1*Z                  #realization of potential outcome

# sensitivity analysis
out.bin <- treatSens(Y~Z+X, trt.family = binomial(link="probit"), nsim = 3, 
                      spy.range = c(0,4), spz.range = c(-2,2),grid.dim = c(5,3),
                      standardize = FALSE, verbose = TRUE)

summary(out.bin) #Prints average value of tau (trt effect) for each cell in grid.

sensPlot(out.bin) # draw contour plot

Run the code above in your browser using DataLab