marCoxph: Marginal Cox Model with Propensity Score Weighting

In PSsurvival: Propensity Score Methods for Survival Analysis

Marginal Cox Model with Propensity Score Weighting

Description

Main user interface for estimating marginal hazard ratios using propensity score weighting. Supports binary and multiple treatment groups with various weighting schemes (ATE, ATT, overlap) and optional trimming. Variance can be estimated via bootstrap or robust sandwich estimator.

Usage

marCoxph(
  data,
  ps_formula,
  time_var,
  event_var,
  reference_level,
  estimand = "ATE",
  att_group = NULL,
  trim = NULL,
  delta = NULL,
  alpha = NULL,
  variance_method = "bootstrap",
  boot_level = "full",
  B = 100,
  parallel = FALSE,
  mc.cores = 2,
  seed = NULL,
  ps_control = list(),
  robust = TRUE
)

Arguments

data	Data frame containing treatment, survival outcome, and covariates.
ps_formula	Formula for propensity score model: treatment ~ covariates.
time_var	Character string specifying the time-to-event variable name.
event_var	Character string specifying the event indicator variable name. Should be coded as 1=event, 0=censored.
reference_level	Treatment level to use as reference in Cox model. MANDATORY. Must be one of the treatment levels.
estimand	Target estimand: "ATE" (average treatment effect), "ATT" (average treatment effect on the treated), or "overlap" (overlap weighting). Default "ATE".
att_group	Target group for ATT. Required if estimand = "ATT".
trim	Trimming method: "symmetric" or "asymmetric". Default NULL (no trimming).
delta	Threshold for symmetric trimming (e.g., 0.1). Required if trim = "symmetric".
alpha	Percentile for asymmetric trimming (e.g., 0.05). Required if trim = "asymmetric".
variance_method	Variance estimation method: "bootstrap" (default) or "robust". "bootstrap" resamples the entire analysis pipeline. "robust" uses the sandwich variance estimator from coxph() without bootstrap.
boot_level	Bootstrap sampling level: "full" (default) or "strata". "full" resamples from entire dataset (standard for observational studies). "strata" resamples within each treatment group preserving group sizes (useful when treatment assignment follows a stratified or fixed-ratio design). Only used if variance_method = "bootstrap".
B	Number of bootstrap iterations. Default 100. Used only if variance_method = "bootstrap".
parallel	Logical. Use parallel bootstrap computation? Default FALSE.
mc.cores	Number of cores for parallel bootstrap. Default 2.
seed	Random seed for bootstrap reproducibility. Default NULL.
ps_control	Control parameters for propensity score model. Default list().
robust	Logical. Use robust (sandwich) variance in Cox model fitting? Default TRUE. When TRUE, coxph() is called with robust = TRUE.

Details

**Analysis Workflow:** 1. Extract treatment variable from ps_formula. 2. Estimate propensity scores using multinomial logistic regression (or logistic for binary treatment). 3. Calculate propensity score weights based on estimand and optional trim. 4. Fit marginal Cox model Surv(time, event) ~ treatment with weights. 5. Estimate variance via bootstrap (resampling full pipeline) or robust sandwich estimator.

**Variance Estimation:** - bootstrap: Resamples data (full or stratified), re-estimates PS and weights, re-fits Cox model. Provides bootstrap SE for log hazard ratios. - robust: Uses robust sandwich variance from coxph() directly. No bootstrap performed (faster but may be less accurate with extreme weights).

**Trimming:** - Symmetric: Crump extension for multiple treatments (Yoshida et al., 2019). - Asymmetric: Sturmer extension for multiple treatments (Yoshida et al., 2019). - Not supported with overlap weights (already bounded [0,1]).

Value

Object of class "marCoxph" containing:

coxph_fitted	Fitted coxph model object.
logHR_est	Named vector of estimated log hazard ratios. Names are formatted as "treatment_var:level" (e.g., "Z:B" for treatment Z, level B vs reference).
logHR_se_robust	Named vector of robust standard errors from coxph.
logHR_se_bootstrap	Named vector of bootstrap standard errors. NULL if variance_method = "robust".
n_coxph_fitted	Named vector of sample sizes per treatment group used in Cox model fitting (after trimming).
events_coxph_fitted	Named vector of event counts per treatment group used in Cox model fitting (after trimming).
variance_method	Variance method used: "bootstrap-full", "bootstrap-strata", or "robust".
estimand	Target estimand used.
att_group	Target group for ATT (NULL if not applicable).
trim_method	Trimming method (NULL if no trimming).
delta	Symmetric trimming threshold (NULL if not applicable).
alpha	Asymmetric trimming threshold (NULL if not applicable).
treatment_var	Name of treatment variable.
treatment_levels	Sorted unique treatment values.
reference_level	Reference level used in Cox model.
n_levels	Number of treatment groups.
n	Number of complete cases used in analysis.
ps_result	Propensity score estimation results.
weight_result	Weight estimation results.
boot_result	Bootstrap results (NULL if variance_method = "robust"). Contains: boot_samples, boot_allocation, n_success_by_group, B.

References

Li, F., & Li, F. (2019). Propensity score weighting for causal inference with multiple treatments. The Annals of Applied Statistics, 13(4), 2389-2415.

Yoshida, K., et al. (2019). Multinomial extension of propensity score trimming methods: A simulation study. American Journal of Epidemiology, 188(3), 609-616.

Examples

# Example 1: Binary treatment with overlap weighting
data(simdata_bin)
result1 <- marCoxph(
  data = simdata_bin,
  ps_formula = Z ~ X1 + X2 + X3 + B1 + B2,
  time_var = "time",
  event_var = "event",
  reference_level = "A",
  estimand = "overlap"
)
summary(result1)

# Example 2: Multiple treatments with ATT and robust variance
data(simdata_multi)
result2 <- marCoxph(
  data = simdata_multi,
  ps_formula = Z ~ X1 + X2 + X3 + B1 + B2,
  time_var = "time",
  event_var = "event",
  reference_level = "C",
  estimand = "ATT",
  att_group = "C",
  variance_method = "robust"
)
summary(result2)

PSsurvival documentation built on Dec. 9, 2025, 9:07 a.m.