Survival Analysis with causaldef

Data Preparation

In a competing risks setting, we need to distinguish between the event of interest (Relapse) and the competing event (Death). We create separate indicators for proper specification.

# Create indicator for Primary Event (Relapse)
hct_outcomes$relapse <- ifelse(hct_outcomes$event_status == "Relapse", 1, 0)

# Create indicator for Competing Event (Death)
# Patients who die without relapse are removed from the risk set for relapse
# but are NOT merely censored in the same way as loss-to-follow-up.
hct_outcomes$death_cr <- ifelse(hct_outcomes$event_status == "Death", 1, 0)

# Check the distribution of events
# Check the distribution of events
table(hct_outcomes$event_status)
#> 
#> Censored    Death  Relapse 
#>      123      131      546

# Visualize Competing Risks
coords <- data.frame(
  name = c("Covariates", "CondIntensity", "Relapse", "Death"),
  x = c(0, -1.5, 1.5, 1.5),
  y = c(1, 0, 0.5, -0.5)
)
edges <- data.frame(
  from = c("Covariates", "Covariates", "Covariates", "CondIntensity", "CondIntensity"),
  to =   c("CondIntensity", "Relapse", "Death", "Relapse", "Death")
)
plot_dag(coords, edges, title = "Competing Risks Structure")

Specification

We using causal_spec_survival() to fully define the structure, including the competing event.

spec_hct <- causal_spec_survival(
  data = hct_outcomes,
  treatment = "conditioning_intensity",
  time = "time_to_event",
  event = "relapse",           # Primary event of interest
  competing_event = "death_cr", # Competing event
  covariates = c("age", "disease_status", "kps", "donor_type"),
  estimand = "ATE"
)
#> ✔ Created survival causal specification: n=800, 546 events

print(spec_hct)
#> 
#> -- Survival Causal Specification ----------------------------------------
#> 
#> * Treatment: conditioning_intensity ( binary )
#> * Time: time_to_event (max = 86.08 )
#> * Event: relapse ( 546 events)
#> * Covariates: age, disease_status, kps, donor_type 
#> * Sample size: 800 
#> * Estimand: ATE 
#> * Competing event: death_cr

Deficiency Estimation

We compare unadjusted analysis (Naive) with Inverse Probability of Treatment Weighting (IPTW) to see if we can reduce the causal deficiency.

# Estimate deficiency
# Note: In survival settings, this measures how well we can recover the 
# interventional survival curves.
results_hct <- estimate_deficiency(
  spec_hct,
  methods = c("unadjusted", "iptw"),
  n_boot = 0 # Set > 0 for bootstrap confidence intervals
)
#> ℹ Inferred treatment value: Reduced
#> ℹ Estimating deficiency: unadjusted
#> ℹ Estimating deficiency: iptw

print(results_hct)
#> 
#> -- Deficiency Proxy Estimates (PS-TV) ------
#> 
#>      Method  Delta SE            Quality
#>  unadjusted 0.3030  - Insufficient (Red)
#>        iptw 0.0644  -   Caution (Yellow)
#> Note: delta is a propensity-score TV proxy (overlap/balance diagnostic).
#> 
#> Best method: iptw (delta = 0.0644 )

Clinical Interpretation of Deficiency:

The deficiency \(\delta\) quantifies how much information is lost by using observational data instead of a randomized trial:

For this HCT Study:

• Unadjusted analysis (high \(\delta\)): Observed distribution differs significantly from target population (High TV distance) -> Biased comparison
• IPTW adjusted (lower \(\delta\)): Reweights patients to balance baseline characteristics (Low TV distance) -> Closer to “what-if-randomized” comparison

Clinical Significance:

• \(\delta\) < 0.05: Excellent control – observational analysis approximates RCT-quality evidence
• 0.05 \(\leq\) \(\delta\) < 0.10: Adequate – reasonable for clinical decision support, acknowledge limitations
• \(\delta\) \(\geq\) 0.10: Caution – substantial bias risk, consider additional covariates or sensitivity analysis

Why This Matters:
In HCT, patients receiving myeloablative conditioning are typically younger and healthier (confounding by indication). IPTW attempts to create a “synthetic randomization” by upweighting underrepresented patient types. Low \(\delta\) validates this approach.

Quantifying Uncertainty: Policy Regret Bounds

The RMST difference provides the observed treatment effect. However, this is biased if treatment assignment was confounded.

The deficiency \(\delta\) we calculated earlier tells us how much we can trust this estimate.

• If \(\delta \approx 0\) (after IPTW), we can interpret the weighted RMST difference as causal (up to sampling error).
• If \(\delta > 0\), policy_regret_bound() converts this information gap into decision-relevant regret bounds on the same scale as the utility (here: months of RMST up to the horizon).

# Calculate policy regret bounds (transfer penalty + minimax floor)
# Utility range is [0, horizon] (months of survival).
# The bound is tied to the previously estimated IPTW deficiency proxy.
bound <- policy_regret_bound(results_hct, utility_range = c(0, horizon), method = "iptw")
#> ℹ Transfer penalty: 1.5457 (delta = 0.0644)
print(bound)
#> 
#> -- Policy Regret Bounds -------------------------------------------------
#> 
#> * Deficiency delta: 0.0644 
#> * Delta mode: point 
#> * Delta method: iptw 
#> * Delta selection: pre-specified method 
#> * Utility range: [0, 24]
#> * Transfer penalty: 1.5457 (additive regret upper bound)
#> * Minimax floor: 0.7728 (worst-case lower bound)
#> 
#> Note: this is a plug-in bound using a deficiency proxy rather than an identified exact deficiency.
#> 
#> Interpretation: Transfer penalty is 6.4 % of utility range given delta

Clinical Interpretation (Transfer Penalty and Minimax Safety Floor):

The bounds answer: “How much regret could confounding add to my decision, and what regret is unavoidable without stronger assumptions?”

Formulas:

• Transfer penalty = \(M \times \delta\) (additive worst-case inflation term in the regret upper bound)
• Minimax safety floor = \((M/2) \times \delta\) (irreducible worst-case regret lower bound)

where \(M\) is the utility range (here, the RMST horizon in months) and \(\delta\) is the deficiency proxy.

What It Means:

• The transfer penalty is the conservative quantity for deployment: it upper-bounds how much worse your deployed policy could be (in regret) due to the observational–interventional information gap.
• The minimax safety floor is a fundamental limit: if this floor is already clinically unacceptable, no algorithm can guarantee safe decisions from the observational data alone.

Clinical Decision Framework (use Transfer Penalty):

1. If Transfer Penalty < 1 month:
   [PASS] Proceed with confidence – Observational evidence is decision-relevant at this horizon.

2. If 1 \(\leq\) Transfer Penalty < 3 months:
   [CAUTION] Proceed with acknowledgment – Often acceptable, but document the bound: “Based on observational evidence with transfer penalty of X months…”

3. If Transfer Penalty \(\geq\) 3 months:
   [STOP] High uncertainty – Consider:

   • Improving adjustment (add covariates, better models)
   • Conducting a prospective RCT
   • Using shared decision-making (acknowledge uncertainty to patients)
   • Limiting recommendations to subgroups with better data

Example:
If \(\delta\) = 0.02 (excellent IPTW performance) and \(M\) = 24 months:
Transfer Penalty = \(24 \times 0.02\) = 0.48 months
Minimax Safety Floor = \((24/2) \times 0.02\) = 0.24 months

Interpretation: If the observed RMST difference is 4 months, then even the conservative transfer-penalty ratio is roughly \(4/0.48 \approx 8:1\), suggesting the decision is robust at this horizon (subject to the modeling assumptions).

Reporting to Clinicians:
“Our observational analysis suggests a 4-month survival advantage for myeloablative conditioning (95% CI: [2.1, 5.9]). After adjustment, the deficiency proxy is \(\delta \approx 0.02\) at a 24-month horizon, implying a transfer penalty of about 0.5 months and a minimax safety floor of about 0.25 months. These bounds are well below the observed benefit, supporting the recommendation, though an RCT would eliminate this residual uncertainty.”

Connection to Deficiency

The key insight:

• Low \(\delta\) (successful adjustment) -> Low transfer penalty / minimax floor (trustworthy decisions)
• High \(\delta\) (poor adjustment) -> High transfer penalty / minimax floor (risky decisions)

The deficiency \(\delta\) and regret bounds work together:

1. \(\delta\) tells you how well you approximated an RCT (information quality)
2. Transfer penalty / minimax floor translate that gap into decision risk on your utility scale

Both should be reported to provide complete transparency about observational evidence quality.