From 4e2a3ca6bc45f9b9e93ca5aa95efd53a2329127d Mon Sep 17 00:00:00 2001 From: Malcolm Barrett Date: Mon, 23 Feb 2026 20:21:22 -0500 Subject: [PATCH] upkeep --- DESCRIPTION | 13 ++- LICENSE | 2 +- LICENSE.md | 2 +- R/ps_calibration.R | 12 ++- R/ps_trunc.R | 9 ++ README.Rmd | 171 ++++--------------------------- README.md | 216 ++++++---------------------------------- man/is_ps_calibrated.Rd | 2 +- man/ps_calibrate.Rd | 10 ++ man/ps_trunc.Rd | 9 ++ 10 files changed, 103 insertions(+), 343 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 8b12350..b6ed0f8 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -4,13 +4,12 @@ Version: 0.0.0.9000 Authors@R: person("Malcolm", "Barrett", , "malcolmbarrett@gmail.com", role = c("aut", "cre", "cph"), comment = c(ORCID = "0000-0003-0299-5825")) -Description: A comprehensive toolkit for propensity score analysis in - causal inference. Supports calculating propensity score weights for - multiple causal estimands across binary, continuous, and categorical - exposures. Provides methods for handling extreme propensity scores - through trimming, truncation, and calibration. Includes inverse - probability weighted estimators that correctly account for propensity - score estimation uncertainty. +Description: Calculates propensity score weights for multiple causal + estimands across binary, continuous, and categorical exposures. + Provides methods for handling extreme propensity scores through + trimming, truncation, and calibration. Includes inverse probability + weighted estimators that correctly account for propensity score + estimation uncertainty. License: MIT + file LICENSE URL: https://r-causal.github.io/propensity/, https://github.com/r-causal/propensity diff --git a/LICENSE b/LICENSE index 188d777..9feedf9 100644 --- a/LICENSE +++ b/LICENSE @@ -1,2 +1,2 @@ -YEAR: 2022 +YEAR: 2026 COPYRIGHT HOLDER: propensity authors diff --git a/LICENSE.md b/LICENSE.md index bfcb661..5bf4666 100644 --- a/LICENSE.md +++ b/LICENSE.md @@ -1,6 +1,6 @@ # MIT License -Copyright (c) 2022 propensity authors +Copyright (c) 2026 propensity authors Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal diff --git a/R/ps_calibration.R b/R/ps_calibration.R index 47ba76c..cf53ef0 100644 --- a/R/ps_calibration.R +++ b/R/ps_calibration.R @@ -131,6 +131,16 @@ pava_weighted <- function(x, y, w = rep(1, length(x))) { #' smoothing was applied). Use [is_ps_calibrated()] to test whether an #' object has been calibrated. #' +#' @references +#' Platt, J. (1999). Probabilistic outputs for support vector machines and +#' comparisons to regularized likelihood methods. *Advances in Large Margin +#' Classifiers*, 61--74. +#' +#' Zadrozny, B., & Elkan, C. (2002). Transforming classifier scores into +#' accurate multiclass probability estimates. *Proceedings of the Eighth ACM +#' SIGKDD International Conference on Knowledge Discovery and Data Mining*, +#' 694--699. \doi{10.1145/775047.775151} +#' #' @seealso [is_ps_calibrated()] to test for calibrated scores; #' [ps_trim()] and [ps_trunc()] for alternative approaches to extreme #' propensity scores; [wt_ate()] and other weight functions that accept @@ -395,7 +405,7 @@ ps_calib_meta <- function(x) { #' #' is_ps_calibrated(ps) #' -#' calibrated <- ps_calibrate(ps, exposure) +#' calibrated <- ps_calibrate(ps, exposure, smooth = FALSE) #' is_ps_calibrated(calibrated) #' #' @export diff --git a/R/ps_trunc.R b/R/ps_trunc.R index 521a6e6..db8560a 100644 --- a/R/ps_trunc.R +++ b/R/ps_trunc.R @@ -59,6 +59,15 @@ #' metadata including `method`, `lower_bound`, `upper_bound`, and #' `truncated_idx` (positions of modified values). #' +#' @references +#' Crump, R. K., Hotz, V. J., Imbens, G. W., & Mitnik, O. A. (2009). Dealing +#' with limited overlap in estimation of average treatment effects. +#' *Biometrika*, 96(1), 187--199. +#' +#' Walker, A. M., Patrick, A. R., Lauer, M. S., et al. (2013). A tool for +#' assessing the feasibility of comparative effectiveness research. +#' *Comparative Effectiveness Research*, 3, 11--20. +#' #' @seealso [ps_trim()] for removing (rather than bounding) extreme values, #' [ps_refit()] for refitting the propensity model after trimming, #' [is_ps_truncated()], [is_unit_truncated()], [ps_trunc_meta()] diff --git a/README.Rmd b/README.Rmd index 93a3ef8..2f5cb2d 100644 --- a/README.Rmd +++ b/README.Rmd @@ -21,7 +21,16 @@ set.seed(123) [![Codecov test coverage](https://codecov.io/gh/r-causal/propensity/graph/badge.svg)](https://app.codecov.io/gh/r-causal/propensity) -propensity provides a toolkit for propensity score analysis in causal inference. It supports multiple causal estimands across binary, categorical, and continuous exposures, handles extreme propensity scores through trimming, truncation, and calibration, and estimates causal effects with valid standard errors via inverse probability weighting. +## Overview + +propensity makes it easy to calculate propensity score weights and use them to estimate causal effects. It supports: + +- Six estimands for binary exposures (ATE, ATT, ATU, ATO, ATM, and entropy weights) +- Binary, categorical, and continuous exposures +- Trimming, truncation, and calibration for extreme propensity scores +- Inverse probability weighted estimation with standard errors that account for propensity score estimation + +You can learn more in `vignette("propensity")`. ## Installation @@ -38,11 +47,9 @@ You can install the development version of propensity from [GitHub](https://gith pak::pak("r-causal/propensity") ``` -## Quick Start - -Estimate a causal effect in three steps: +## Usage -```{r quick-start, message = FALSE, warning = FALSE} +```{r usage, message = FALSE, warning = FALSE} library(propensity) # Simulate data with a confounder, binary exposure, and binary outcome @@ -63,159 +70,25 @@ outcome_mod <- glm(y ~ z, data = dat, family = binomial(), weights = wts) ipw(ps_mod, outcome_mod) ``` -`ipw()` performs inverse probability weighted estimation. It accounts for uncertainty in the estimated propensity scores when computing standard errors, producing valid confidence intervals and p-values. - -## Multiple Estimands - -propensity supports six causal estimands for binary exposures: - -```{r estimands, eval = FALSE} -ps <- fitted(ps_mod) +`ipw()` uses linearization to account for uncertainty in the estimated propensity scores when computing standard errors. -wt_ate(ps, z) # Average Treatment Effect -wt_att(ps, z) # Average Treatment Effect on the Treated -wt_atu(ps, z) # Average Treatment Effect on the Untreated -wt_ato(ps, z) # Average Treatment Effect, Overlap population -wt_atm(ps, z) # Average Treatment Effect, Matched population -wt_entropy(ps, z) # Entropy-weighted ATE -``` +## Estimands -Choose your estimand based on your research question: +Each weight function targets a different population: | Estimand | Target population | Function | |----------|-------------------|----------| | **ATE** | Entire population | `wt_ate()` | | **ATT** | Treated units | `wt_att()` | | **ATU** | Untreated units | `wt_atu()` (alias: `wt_atc()`) | -| **ATO** | Overlap population -- units with the most equipoise between groups | `wt_ato()` | -| **ATM** | Matched population -- mimics 1:1 matching | `wt_atm()` | -| **Entropy** | Entropy-balanced compromise between ATE and overlap | `wt_entropy()` | - -ATO and ATM produce naturally bounded weights, making them good alternatives when ATE weights are highly variable. - -## Flexible Input - -Weight functions accept propensity scores in several forms: - -```{r inputs, eval = FALSE} -# GLM object -- extracts fitted values and exposure automatically (recommended) -wt_ate(ps_mod) - -# Numeric vector of propensity scores -wt_ate(ps, z) - -# Data frame of class probabilities (e.g., from multinomial models) -ps_df <- data.frame(control = 1 - ps, treated = ps) -wt_ate(ps_df, z) -``` - -## Handling Extreme Weights - -Propensity scores near 0 or 1 produce extreme weights that inflate variance. propensity offers four complementary strategies: - -1. **Switch estimand** -- `wt_ato()` and `wt_atm()` produce bounded weights by design. -2. **Trim** -- remove observations with extreme scores (sets them to `NA`). -3. **Truncate** -- winsorize extreme scores to a fixed range. -4. **Calibrate** -- adjust scores so they better reflect true treatment probabilities. - -```{r extreme-weights, message = FALSE, warning = FALSE} -ps <- fitted(ps_mod) - -# Diagnose: inspect the weight distribution -summary(wt_ate(ps, z)) - -# 1. Switch estimand -- overlap weights are bounded by design -summary(wt_ato(ps, z)) - -# 2. Trim -- remove observations with extreme propensity scores -ps_trimmed <- ps_trim(ps, method = "adaptive") -summary(wt_ate(ps_trimmed, z)) - -# 3. Truncate -- bound extreme propensity scores to a range -ps_truncated <- ps_trunc(ps, lower = 0.05, upper = 0.95) -summary(wt_ate(ps_truncated, z)) - -# 4. Calibrate -- adjust scores to better reflect true probabilities -ps_calibrated <- ps_calibrate(ps, z) -summary(wt_ate(ps_calibrated, z)) -``` - -After trimming, refit the propensity score model on the retained subset so the scores reflect the trimmed population: - -```{r refit, eval = FALSE} -ps_refitted <- ps_refit(ps_trimmed, ps_mod) -wts_refitted <- wt_ate(ps_refitted, z) -``` - -## Advanced Features - -### Weight Stabilization - -Stabilized weights can reduce variance. Stabilization is supported for `wt_ate()` and `wt_cens()`: - -```{r stabilization, message = FALSE} -summary(wt_ate(ps, z, stabilize = TRUE)) -``` - -### Continuous Exposures - -For continuous treatments, weights use the ratio of the marginal to conditional density. Stabilization is strongly recommended: - -```{r continuous, eval = FALSE} -# Fit a model for the continuous exposure -continuous_exposure <- rnorm(n, mean = 0.5 * x1) -dat$a <- continuous_exposure -ps_continuous <- glm(a ~ x1, data = dat, family = gaussian()) - -# Density-ratio weights (stabilization strongly recommended) -wts_continuous <- wt_ate(ps_continuous, stabilize = TRUE) -``` - -### Categorical Exposures - -For multi-level treatments, supply a matrix or data frame of class probabilities with one column per treatment level: - -```{r categorical, eval = FALSE} -# Multinomial propensity scores (one column per treatment level) -ps_matrix <- matrix(c(0.3, 0.5, 0.2), ncol = 3, nrow = n, byrow = TRUE) -categorical_exposure <- factor(sample(1:3, n, replace = TRUE)) +| **ATO** | Overlap population | `wt_ato()` | +| **ATM** | Matched population | `wt_atm()` | +| **Entropy** | Entropy-balanced population | `wt_entropy()` | -wt_ate(ps_matrix, categorical_exposure, exposure_type = "categorical") - -# For ATT with categorical exposures, specify the focal level -wt_att(ps_matrix, categorical_exposure, .focal_level = "2") -``` - -### Censoring Weights - -`wt_cens()` calculates inverse probability of censoring weights for survival or longitudinal analyses. These address informative censoring -- not treatment assignment: - -```{r censoring, eval = FALSE} -# Model the probability of being uncensored -cens_mod <- glm(uncensored ~ x1 + x2, data = dat, family = binomial()) - -# Censoring weights (uses the same formula as wt_ate()) -wts_cens <- wt_cens(cens_mod) - -# Combine with treatment weights for a doubly-weighted analysis -wts_combined <- wt_ate(ps_mod) * wts_cens -``` - -### Calibration Methods - -`ps_calibrate()` supports two calibration methods: - -```{r calibration, eval = FALSE} -# Logistic calibration (default) -- fits a logistic regression of exposure on -# predicted propensity scores -ps_calibrate(ps, z, method = "logistic") - -# Isotonic regression calibration -- fits a monotone step function; useful for -# non-smooth relationships with large samples -ps_calibrate(ps, z, method = "isoreg") -``` +ATO and ATM weights are bounded by construction, making them a good alternative when ATE weights are highly variable. -## Learn More +## Learn more +- [Causal Inference in R](https://www.r-causal.org/) -- A book on causal inference methods in R +- `vignette("propensity")` -- Getting started with propensity score weighting - [propensity package documentation](https://r-causal.github.io/propensity/) -- Full reference and articles -- [Causal Inference in R](https://www.r-causal.org/) -- A comprehensive guide to causal inference methods in R diff --git a/README.md b/README.md index 73f8f47..e647b66 100644 --- a/README.md +++ b/README.md @@ -10,11 +10,19 @@ coverage](https://codecov.io/gh/r-causal/propensity/graph/badge.svg)](https://app.codecov.io/gh/r-causal/propensity) -propensity provides a toolkit for propensity score analysis in causal -inference. It supports multiple causal estimands across binary, -categorical, and continuous exposures, handles extreme propensity scores -through trimming, truncation, and calibration, and estimates causal -effects with valid standard errors via inverse probability weighting. +## Overview + +propensity makes it easy to calculate propensity score weights and use +them to estimate causal effects. It supports: + +- Six estimands for binary exposures (ATE, ATT, ATU, ATO, ATM, and + entropy weights) +- Binary, categorical, and continuous exposures +- Trimming, truncation, and calibration for extreme propensity scores +- Inverse probability weighted estimation with standard errors that + account for propensity score estimation + +You can learn more in `vignette("propensity")`. ## Installation @@ -33,9 +41,7 @@ You can install the development version of propensity from pak::pak("r-causal/propensity") ``` -## Quick Start - -Estimate a causal effect in three steps: +## Usage ``` r library(propensity) @@ -74,187 +80,31 @@ ipw(ps_mod, outcome_mod) #> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ``` -`ipw()` performs inverse probability weighted estimation. It accounts -for uncertainty in the estimated propensity scores when computing -standard errors, producing valid confidence intervals and p-values. - -## Multiple Estimands - -propensity supports six causal estimands for binary exposures: - -``` r -ps <- fitted(ps_mod) - -wt_ate(ps, z) # Average Treatment Effect -wt_att(ps, z) # Average Treatment Effect on the Treated -wt_atu(ps, z) # Average Treatment Effect on the Untreated -wt_ato(ps, z) # Average Treatment Effect, Overlap population -wt_atm(ps, z) # Average Treatment Effect, Matched population -wt_entropy(ps, z) # Entropy-weighted ATE -``` - -Choose your estimand based on your research question: - -| Estimand | Target population | Function | -|----|----|----| -| **ATE** | Entire population | `wt_ate()` | -| **ATT** | Treated units | `wt_att()` | -| **ATU** | Untreated units | `wt_atu()` (alias: `wt_atc()`) | -| **ATO** | Overlap population – units with the most equipoise between groups | `wt_ato()` | -| **ATM** | Matched population – mimics 1:1 matching | `wt_atm()` | -| **Entropy** | Entropy-balanced compromise between ATE and overlap | `wt_entropy()` | - -ATO and ATM produce naturally bounded weights, making them good -alternatives when ATE weights are highly variable. - -## Flexible Input - -Weight functions accept propensity scores in several forms: - -``` r -# GLM object -- extracts fitted values and exposure automatically (recommended) -wt_ate(ps_mod) - -# Numeric vector of propensity scores -wt_ate(ps, z) - -# Data frame of class probabilities (e.g., from multinomial models) -ps_df <- data.frame(control = 1 - ps, treated = ps) -wt_ate(ps_df, z) -``` - -## Handling Extreme Weights - -Propensity scores near 0 or 1 produce extreme weights that inflate -variance. propensity offers four complementary strategies: - -1. **Switch estimand** – `wt_ato()` and `wt_atm()` produce bounded - weights by design. -2. **Trim** – remove observations with extreme scores (sets them to - `NA`). -3. **Truncate** – winsorize extreme scores to a fixed range. -4. **Calibrate** – adjust scores so they better reflect true treatment - probabilities. - -``` r -ps <- fitted(ps_mod) - -# Diagnose: inspect the weight distribution -summary(wt_ate(ps, z)) -#> Min. 1st Qu. Median Mean 3rd Qu. Max. -#> 1.380 1.806 1.956 2.000 2.154 2.902 - -# 1. Switch estimand -- overlap weights are bounded by design -summary(wt_ato(ps, z)) -#> Min. 1st Qu. Median Mean 3rd Qu. Max. -#> 0.2752 0.4463 0.4889 0.4912 0.5358 0.6554 - -# 2. Trim -- remove observations with extreme propensity scores -ps_trimmed <- ps_trim(ps, method = "adaptive") -summary(wt_ate(ps_trimmed, z)) -#> Min. 1st Qu. Median Mean 3rd Qu. Max. -#> 1.380 1.806 1.956 2.000 2.154 2.902 +`ipw()` uses linearization to account for uncertainty in the estimated +propensity scores when computing standard errors. -# 3. Truncate -- bound extreme propensity scores to a range -ps_truncated <- ps_trunc(ps, lower = 0.05, upper = 0.95) -summary(wt_ate(ps_truncated, z)) -#> Min. 1st Qu. Median Mean 3rd Qu. Max. -#> 1.380 1.806 1.956 2.000 2.154 2.902 +## Estimands -# 4. Calibrate -- adjust scores to better reflect true probabilities -ps_calibrated <- ps_calibrate(ps, z) -summary(wt_ate(ps_calibrated, z)) -#> Min. 1st Qu. Median Mean 3rd Qu. Max. -#> 1.401 1.804 1.956 2.000 2.156 2.881 -``` - -After trimming, refit the propensity score model on the retained subset -so the scores reflect the trimmed population: +Each weight function targets a different population: -``` r -ps_refitted <- ps_refit(ps_trimmed, ps_mod) -wts_refitted <- wt_ate(ps_refitted, z) -``` +| Estimand | Target population | Function | +|-------------|-----------------------------|--------------------------------| +| **ATE** | Entire population | `wt_ate()` | +| **ATT** | Treated units | `wt_att()` | +| **ATU** | Untreated units | `wt_atu()` (alias: `wt_atc()`) | +| **ATO** | Overlap population | `wt_ato()` | +| **ATM** | Matched population | `wt_atm()` | +| **Entropy** | Entropy-balanced population | `wt_entropy()` | -## Advanced Features - -### Weight Stabilization - -Stabilized weights can reduce variance. Stabilization is supported for -`wt_ate()` and `wt_cens()`: - -``` r -summary(wt_ate(ps, z, stabilize = TRUE)) -#> Min. 1st Qu. Median Mean 3rd Qu. Max. -#> 0.7036 0.9017 0.9817 0.9999 1.0732 1.4536 -``` - -### Continuous Exposures - -For continuous treatments, weights use the ratio of the marginal to -conditional density. Stabilization is strongly recommended: - -``` r -# Fit a model for the continuous exposure -continuous_exposure <- rnorm(n, mean = 0.5 * x1) -dat$a <- continuous_exposure -ps_continuous <- glm(a ~ x1, data = dat, family = gaussian()) - -# Density-ratio weights (stabilization strongly recommended) -wts_continuous <- wt_ate(ps_continuous, stabilize = TRUE) -``` - -### Categorical Exposures - -For multi-level treatments, supply a matrix or data frame of class -probabilities with one column per treatment level: - -``` r -# Multinomial propensity scores (one column per treatment level) -ps_matrix <- matrix(c(0.3, 0.5, 0.2), ncol = 3, nrow = n, byrow = TRUE) -categorical_exposure <- factor(sample(1:3, n, replace = TRUE)) - -wt_ate(ps_matrix, categorical_exposure, exposure_type = "categorical") - -# For ATT with categorical exposures, specify the focal level -wt_att(ps_matrix, categorical_exposure, .focal_level = "2") -``` - -### Censoring Weights - -`wt_cens()` calculates inverse probability of censoring weights for -survival or longitudinal analyses. These address informative censoring – -not treatment assignment: - -``` r -# Model the probability of being uncensored -cens_mod <- glm(uncensored ~ x1 + x2, data = dat, family = binomial()) - -# Censoring weights (uses the same formula as wt_ate()) -wts_cens <- wt_cens(cens_mod) - -# Combine with treatment weights for a doubly-weighted analysis -wts_combined <- wt_ate(ps_mod) * wts_cens -``` - -### Calibration Methods - -`ps_calibrate()` supports two calibration methods: - -``` r -# Logistic calibration (default) -- fits a logistic regression of exposure on -# predicted propensity scores -ps_calibrate(ps, z, method = "logistic") - -# Isotonic regression calibration -- fits a monotone step function; useful for -# non-smooth relationships with large samples -ps_calibrate(ps, z, method = "isoreg") -``` +ATO and ATM weights are bounded by construction, making them a good +alternative when ATE weights are highly variable. -## Learn More +## Learn more +- [Causal Inference in R](https://www.r-causal.org/) – A book on causal + inference methods in R +- `vignette("propensity")` – Getting started with propensity score + weighting - [propensity package documentation](https://r-causal.github.io/propensity/) – Full reference and articles -- [Causal Inference in R](https://www.r-causal.org/) – A comprehensive - guide to causal inference methods in R diff --git a/man/is_ps_calibrated.Rd b/man/is_ps_calibrated.Rd index ca659a1..b8a56ca 100644 --- a/man/is_ps_calibrated.Rd +++ b/man/is_ps_calibrated.Rd @@ -22,7 +22,7 @@ exposure <- rbinom(100, 1, ps) is_ps_calibrated(ps) -calibrated <- ps_calibrate(ps, exposure) +calibrated <- ps_calibrate(ps, exposure, smooth = FALSE) is_ps_calibrated(calibrated) } diff --git a/man/ps_calibrate.Rd b/man/ps_calibrate.Rd index 00edeaa..93e749f 100644 --- a/man/ps_calibrate.Rd +++ b/man/ps_calibrate.Rd @@ -102,6 +102,16 @@ if (rlang::is_installed("mgcv")) { cal_smooth <- ps_calibrate(ps, exposure) } +} +\references{ +Platt, J. (1999). Probabilistic outputs for support vector machines and +comparisons to regularized likelihood methods. \emph{Advances in Large Margin +Classifiers}, 61--74. + +Zadrozny, B., & Elkan, C. (2002). Transforming classifier scores into +accurate multiclass probability estimates. \emph{Proceedings of the Eighth ACM +SIGKDD International Conference on Knowledge Discovery and Data Mining}, +694--699. \doi{10.1145/775047.775151} } \seealso{ \code{\link[=is_ps_calibrated]{is_ps_calibrated()}} to test for calibrated scores; diff --git a/man/ps_trunc.Rd b/man/ps_trunc.Rd index 79fdb0c..d9a39d4 100644 --- a/man/ps_trunc.Rd +++ b/man/ps_trunc.Rd @@ -118,6 +118,15 @@ wt_ate(ps_t, .exposure = z) # Inspect which observations were truncated is_unit_truncated(ps_t) +} +\references{ +Crump, R. K., Hotz, V. J., Imbens, G. W., & Mitnik, O. A. (2009). Dealing +with limited overlap in estimation of average treatment effects. +\emph{Biometrika}, 96(1), 187--199. + +Walker, A. M., Patrick, A. R., Lauer, M. S., et al. (2013). A tool for +assessing the feasibility of comparative effectiveness research. +\emph{Comparative Effectiveness Research}, 3, 11--20. } \seealso{ \code{\link[=ps_trim]{ps_trim()}} for removing (rather than bounding) extreme values,