loess.psa: Graphic for data and loess-based estimate of effect size after propensity score adjustment

Description

Plots data points using propesity scores vs. the response, separately for treatment and control groups; points are distinguished by both type and color for the two groups. Also shows (non-linear, loess-based) regression curves for both groups. The loess regresion curves are then used to derive an overall estimate of effect size (based on number and/or location of strata as set by the user). Several other statistics are also provided, for both description and inference. Graphic motivated by a suggestion of R. L. Obenchain.

Usage

loess.psa(
  response,
  treatment = NULL,
  propensity = NULL,
  family = "gaussian",
  span = 0.7,
  degree = 1,
  minsize = 5,
  xlim = c(0, 1),
  colors = c("dark blue", "dark green", "blue", "dark green"),
  legend.xy = "topleft",
  legend = NULL,
  int = 10,
  lines = TRUE,
  strata.lines = TRUE,
  rg = TRUE,
  xlab = "Estimated Propensity Scores",
  ylab = "Response",
  pch = c(16, 1),
  ...
)

Value

In addition to the plot, the function returns a list with the following components:

ATE: Estimated effect size based upon (number of) strata defined by int; that is, this is the Average Treatment Effect, after propensity-based adjustment.
se.wtd: Weighted standard error based on pooling of within-strata variance estimates.
CI.95: Approximate 95% confidence interval for the overall effect size (conditional on the specification of int).
summary.strata: A table with rows corresponding to strata; first two columns show counts (by statum) for both control and treatment; followed by mean differences for all strata. for control and treatment, followed by mean differences for all strata. The weighted average difference yields the effect size noted above.

Arguments

response: Either a numeric vector containing the response of interest in a propensity score analysis, or a three column array containing response, treatment and strata.
treatment: Binary variable of same length as response; 0 for 'control,' 1 for 'treatment.'
propensity: Numeric vector of estimated propensity scores.
family: Passed to loess. Either "gaussian" (default) or "symmetric".
span: Parameter passed to loess governing degree of smoothing. Default = 0.7.
degree: Parameter passed to loess governing degree of polynomials used. Default = 1
minsize: Integer. Determines the minimum number of observations in each stratum treatment group allowed. If one of the treatment groups in a given statum does not meet this minsize, then all observations in this stratum are ignored as far as the effect size calculation is concerned.
xlim: Binary vector (min, max) providing the horizontal axis minimum and maximum. Default is c(0, 1).
colors: List of four colors used for control points, treatment points, control loess line, treatment loess line respectively. Default = c("seagreen3", "goldenrod1", "seagreen4", "goldenrod3").
legend.xy: Coordinates for legend box, see legend. Default = "topleft".
legend: Binary character vector containing the text of the legend. Default is taken from treatment.
int: Integer or ordered vector. If an integer is used, it represents the maximum number of equally sized strata. Alternatively, it may be a vector of cuts of the unit interval. Lower and upper ends need not be included. See examples. Default = 10.
lines: Logical; fitted loess values are plotted by default as points. If true, values are plotted as two lines.
strata.lines: Logical; default = TRUE. Creates light vertical lines that delineate strata.
rg: Logical; if TRUE (default) then rug plots are given for treatment and control propensity score and response distributions.
xlab: X axis label, default = "Estimated Propensity Scores".
ylab: Y axis label, default = "Response".
pch: Character types for plotted points, default = c(16, 1). Note: must be of length 2 to allow different plotting points for each treatment.
...: Optional parameters passed to points command.

Author

James E. Helmreich James.Helmreich@Marist.edu

Robert M. Pruzek RMPruzek@yahoo.com

Examples

Run this code


#Artificial example where ATE should be 1 over all of (0,1).
response1 <- c(rep(1, 100), rep(2, 100), rep(3, 100)) + rnorm(300, 0, .5)
response0 <- c(rep(0, 100), rep(1, 100), rep(2, 100)) + rnorm(300, 0, .5)
response <- c(response1, response0)
treatment <- c(rep(1, 300), rep(0, 300))
propensity <- rep(seq(.01, .99, (.98/299)), 2)
a <- data.frame(response, treatment, propensity)
loess.psa(a, span = .15, degree = 1, int = c(0, .33, .67, 1))


#Artificial example where estimates are unstable with varying
#numbers of strata. Note: sometimes get empty treatment/strata error.
rr <- c(rnorm(150, 3, .75), rnorm(700, 0, .75), rnorm(150, 3, .75),
     rnorm(150, -3, .75), rnorm(700, 0, .75), rnorm(150, -3, .75))
tt <- c(rep(1, 1000),rep(0, 1000))
pp <- NULL
for(i in 1:1000){pp <- c(pp, rnorm(1, 0, .05) + .00045*i + .25)}
for(i in 1:1000){pp <- c(pp, rnorm(1, 0, .05) + .00045*i + .4)}
a <- data.frame(rr, tt, pp)
loess.psa(a, span=.5, cex = .6)

#Using strata of possible interest as determined by loess lines.
data(lindner)
attach(lindner)
lindner.ps <- glm(abcix ~ stent + height + female +
      diabetic + acutemi + ejecfrac + ves1proc,
      data = lindner, family = binomial)
loess.psa(log(cardbill), abcix, lindner.ps$fitted,
         int = c(.37, .56, .87, 1), lines = TRUE)
         abline(v=c(.37, 56, .87))