Skip to content

Commit

Permalink
Added propensity score weighting.
Browse files Browse the repository at this point in the history
  • Loading branch information
jbryer committed May 16, 2023
1 parent 496083d commit 3f4d036
Show file tree
Hide file tree
Showing 11 changed files with 177 additions and 22 deletions.
6 changes: 4 additions & 2 deletions .github/workflows/render-readme.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -11,8 +11,10 @@ jobs:
runs-on: macOS-latest
steps:
- uses: actions/checkout@v2
- uses: r-lib/actions/setup-r@v1
- uses: r-lib/actions/setup-pandoc@v1
- uses: r-lib/actions/setup-pandoc@v2
- uses: r-lib/actions/setup-r@v2
with:
use-public-rspm: true
- name: Install rmarkdown
run: Rscript -e 'install.packages(c("rmarkdown", "badger"))'
- name: Render README
Expand Down
4 changes: 4 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@ export(boot.matching)
export(boot.matchit)
export(boot.rpart)
export(boot.strata)
export(boot.weighting)
export(getPSAbootMethods)
export(matrixplot)
export(psa.strata)
Expand All @@ -40,14 +41,17 @@ importFrom(party,where)
importFrom(psych,describe)
importFrom(psych,describeBy)
importFrom(stats,as.formula)
importFrom(stats,binomial)
importFrom(stats,cor)
importFrom(stats,density)
importFrom(stats,fitted)
importFrom(stats,glm)
importFrom(stats,lm)
importFrom(stats,predict)
importFrom(stats,qnorm)
importFrom(stats,qt)
importFrom(stats,quantile)
importFrom(stats,quasibinomial)
importFrom(stats,sd)
importFrom(stats,t.test)
importFrom(stats,update.formula)
Expand Down
1 change: 1 addition & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
# PSAboot 1.3.7

* Added propensity score weighting.
* Fixed a bug where the balance table wasn't combined correctly if the covariates were not specified in the correct order.

# PSAboot 1.3.6
Expand Down
2 changes: 1 addition & 1 deletion Poster/Poster.R
Original file line number Diff line number Diff line change
Expand Up @@ -74,7 +74,7 @@ matrixplot(tutoring.boot)
dev.off()

tutoring.balance <- balance(tutoring.boot)
summary(tutoring.balance)
print(tutoring.balance)

plot(tutoring.balance) + theme(legend.position = 'bottom')
ggsave('Poster/balance.pdf', width = 12, height = 8, units = 'in', dpi = 100)
Expand Down
3 changes: 2 additions & 1 deletion R/PSAboot.R
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,8 @@ getPSAbootMethods <- function() {
'ctree' = boot.ctree,
'rpart' = boot.rpart,
'Matching' = boot.matching,
'MatchIt' = boot.matchit))
'MatchIt' = boot.matchit,
'Weighting' = boot.weighting))
invisible(methods)
}

Expand Down
72 changes: 72 additions & 0 deletions R/boot.weighting.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,72 @@
#' Calculates propensity score weights.
#'
#' @param treatment a logical vector for treatment status.
#' @param ps numeric vector of propensity scores
#' @param estimand character string indicating which estimand to be used. Possible
#' values are
#' ATE (average treatment effect),
#' ATT (average treatment effect for the treated),
#' ATC (average treatement effect for the controls),
#' ATM (Average Treatment Effect Among the Evenly Matchable),
#' ATO (Average Treatment Effect Among the Overlap Populatio)
calculate_ps_weights <- function(treatment, ps, estimand = 'ATE') {
# TODO: this is a copy of a function from the psa package. Use that once it is on CRAN
weights <- NA
if(estimand == 'ATE') {
weights <- (treatment / ps) + ((1 - treatment) / (1 - ps))
} else if(estimand == 'ATT') {
weights <- ((ps * treatment) / ps) + ((ps * (1 - treatment)) / (1 - ps))
} else if(estimand == 'ATC') {
weights <- (((1 - ps) * treatment) / ps) + (((1 - ps) * (1 - treatment)) / (1 - ps))
} else if(estimand == 'ATM') {
weights <- pmin(ps, 1 - ps) / (treatment * ps + (1 - treatment) * (1 - ps))
} else if(estimand == 'ATO') {
weights <- (1 - ps) * treatment + ps * (1 - treatment)
} else {
stop(paste0('Invalid estimand specified: ', estimand))
}
return(weights)
}

#' Propensity score weighting implementation for bootstrapping.
#'
#' @inherit boot.strata return params
#' @param estimand which treatment effect to estimate. Values can be ATE, ATT,
#' ATC, or ATM.
#' @importFrom stats quasibinomial binomial lm
#' @export
boot.weighting <- function(Tr, Y, X, X.trans, formu, estimand = 'ATE', ...) {
formu.treat <- update.formula(formu, 'treat ~ .')
df <- cbind(treat = Tr, X)
ps <- fitted(glm(formu.treat,
data = df,
family = binomial(link = 'logit')))

weights <- calculate_ps_weights(treatment = Tr,
ps = ps,
estimand = estimand)

X.trans.scaled <- lapply(X.trans, scale) |> as.data.frame()
X.trans.scaled$treat <- Tr
lm_balance <- glm(treat ~ .,
data = X.trans.scaled,
family = quasibinomial(link = 'logit'),
weights = weights) |> summary()

te_lm <- lm(formula = Y ~ treat,
data = data.frame(Y = Y, treat = Tr),
weights = weights) |> summary()

se <- te_lm$coefficients[2,]['Std. Error'] |> unname()
te <- te_lm$coefficients[2,]['Estimate'] |> unname()

return(list(
summary = c(estimate = te,
ci.min = te - 1.96 * se,
ci.max = te + 1.96 * se,
se.wtd = se,
approx.t = unname(te_lm$coefficients[2,]['t value'])),
details = te_lm,
balance = lm_balance$coefficients[,'Estimate'][-1]
))
}
5 changes: 4 additions & 1 deletion README.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@ editor_options:
---


# <img src="man/figures/PSAboot.png" align="right" width="120" align="right" /> An R Package for Bootstrapping Propensity Score Analysis
# <img src="man/figures/PSAboot.png" align="right" width="120" align="right" /> Bootstrapping Propensity Score Analysis

<!-- badges: start -->
[![R-CMD-check](https://github.com/jbryer/PSAboot/actions/workflows/R-CMD-check.yaml/badge.svg)](https://github.com/jbryer/PSAboot/actions/workflows/R-CMD-check.yaml)
Expand All @@ -14,6 +14,9 @@ editor_options:
[![CRAN Status](https://badges.cranchecks.info/flavor/release/PSAboot.svg)](https://cran.r-project.org/web/checks/check_results_PSAboot.html)
<!-- badges: end -->

Package website: https://jbryer.github.io/PSAboot/
Poster: https://github.com/jbryer/PSAboot/blob/master/Poster/PSAboot_Poster.pdf


As the popularity of propensity score methods for estimating causal effects in observational studies increase, the choices researchers have for which methods to use has also increased. Estimated treatment effects may be sensitive to choice of method. One approach to test the sensitivity of method choice is to test the null hypothesis more than once using more than one method (Rosenbaum, 2012). With the wide availability of high power computers resampling methods such as bootstrapping (Efron, 1979) have become popular for providing more estimates of the sampling distribution. This paper introduces the `PSAboot` R package that provides functions for bootstrapping propensity score methods. It deviates from traditional bootstrapping methods by allowing for different sampling specifications for treatment and control groups, mainly to ensure the ratio of treatment-to-control observations are consistent. This approach can also be used in situations where there is imbalance between the number of treatment and control observations by allowing for up and/or down sampling. Lastly, by estimating balance statistics and treatment effects for each bootstrap sample we can compare the distributions across multiple propensity score methods to examine the relative performance of these methods.

Expand Down
39 changes: 23 additions & 16 deletions README.md
Original file line number Diff line number Diff line change
@@ -1,31 +1,38 @@

# <img src="man/figures/PSAboot.png" align="right" width="120" align="right" /> An R Package for Bootstrapping Propensity Score Analysis
# <img src="man/figures/PSAboot.png" align="right" width="120" align="right" /> Bootstrapping Propensity Score Analysis

<!-- badges: start -->

[![R-CMD-check](https://github.com/jbryer/PSAboot/actions/workflows/R-CMD-check.yaml/badge.svg)](https://github.com/jbryer/PSAboot/actions/workflows/R-CMD-check.yaml)
[![](https://img.shields.io/badge/devel%20version-1.3.6-blue.svg)](https://github.com/jbryer/PSAboot)
[![](https://img.shields.io/badge/devel%20version-1.3.7-blue.svg)](https://github.com/jbryer/PSAboot)
[![](https://www.r-pkg.org/badges/version/PSAboot)](https://cran.r-project.org/package=PSAboot)
[![CRAN
Status](https://badges.cranchecks.info/flavor/release/PSAboot.svg)](https://cran.r-project.org/web/checks/check_results_PSAboot.html)
<!-- badges: end -->

Package website: <https://jbryer.github.io/PSAboot/>
Poster:
<https://github.com/jbryer/PSAboot/blob/master/Poster/PSAboot_Poster.pdf>

As the popularity of propensity score methods for estimating causal
effects in observational studies increase, the choices researchers have
for which methods to use has also increased. Rosenbaum (2012) suggested
that there are benefits for testing the null hypothesis more than once
in observational studies. With the wide availability of high power
computers resampling methods such as bootstrapping (Efron, 1979) have
become popular for providing more stable estimates of the sampling
distribution. This paper introduces the `PSAboot` package for R that
provides functions for bootstrapping propensity score methods. It
deviates from traditional bootstrapping methods by allowing for
different sampling specifications for treatment and control groups,
mainly to ensure the ratio of treatment-to-control observations are
maintained. Additionally, this framework will provide estimates using
multiple methods for each bootstrap sample. Two examples are discussed:
the classic National Work Demonstration and PSID (Lalonde, 1986) study
and a study on tutoring effects on student grades.
for which methods to use has also increased. Estimated treatment effects
may be sensitive to choice of method. One approach to test the
sensitivity of method choice is to test the null hypothesis more than
once using more than one method (Rosenbaum, 2012). With the wide
availability of high power computers resampling methods such as
bootstrapping (Efron, 1979) have become popular for providing more
estimates of the sampling distribution. This paper introduces the
`PSAboot` R package that provides functions for bootstrapping propensity
score methods. It deviates from traditional bootstrapping methods by
allowing for different sampling specifications for treatment and control
groups, mainly to ensure the ratio of treatment-to-control observations
are consistent. This approach can also be used in situations where there
is imbalance between the number of treatment and control observations by
allowing for up and/or down sampling. Lastly, by estimating balance
statistics and treatment effects for each bootstrap sample we can
compare the distributions across multiple propensity score methods to
examine the relative performance of these methods.

## Installation

Expand Down
41 changes: 41 additions & 0 deletions man/boot.weighting.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

24 changes: 24 additions & 0 deletions man/calculate_ps_weights.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

2 changes: 1 addition & 1 deletion vignettes/PSAboot.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -243,7 +243,7 @@ hist(tutoring.boot)
boxplot(tutoring.boot)
```

```{r tutoring.matrixplot}
```{r tutoring.matrixplot, fig.width = 12, fig.height = 12}
matrixplot(tutoring.boot)
```

Expand Down

0 comments on commit 3f4d036

Please sign in to comment.