Cox Proportional Hazards Regression

The Cox proportional hazards (PH) model (Cox 1972) specifies the hazard as h(t | X) = h0(t) · exp(β'X), where h0(t) is an arbitrary, unestimated baseline hazard and β are log hazard ratios. Its central trick is the partial likelihood: by conditioning on the risk set at each observed event time, h0(t) cancels and the estimator depends only on the rank order of event times, not their spacing. That is what makes Cox "semiparametric" — no functional form is imposed on how risk evolves with time, only on how covariates multiply it. It is the workhorse of pharmacoepidemiology and HEOR because it accommodates right-censoring, left-truncation (delayed entry), and time-varying covariates with mature, regulator-accepted software. This entry covers the standard model; pushing exposure or covariates inside the partial likelihood at each event time is its own concept (`standard-cox-time-dependent`).

Core conceptual distinction

The quantity Cox estimates is a conditional (covariate-specific) hazard ratio, and three properties of it are routinely misunderstood.

1. Conditional vs marginal, and non-collapsibility. Even with no confounding, the HR adjusted for a prognostic covariate is generally not equal to the crude HR — the HR is non-collapsible (Stensrud & Hernán 2020). Adding a strong predictor moves the HR away from the null even when that predictor is not a confounder, so "the HR changed when I adjusted" is not evidence of confounding.
2. The HR is a ratio of instantaneous rates among survivors, so under PH it still mixes a true causal effect with selection: as the more-protected arm depletes its susceptibles, the surviving risk sets become non-comparable, and a constant true effect can produce a time-varying observed HR.
3. A single HR is a weighted average over follow-up when PH is violated; the weights depend on the censoring distribution, making the number partly an artifact of study duration. For decision-makers, an absolute or marginal summary — cumulative incidence, restricted mean survival time (RMST), or a survival difference at a fixed horizon — is often the more interpretable estimand and should be pre-specified alongside (or instead of) the HR.

Pros, cons, and trade-offs

• vs fully parametric survival models (Weibull, exponential, generalized gamma): Cox is robust because it never commits to a shape for h0(t), and it is the easiest survival model to communicate. Cost: it is less efficient when a parametric form is correct, it does not directly yield absolute survival curves (you need a Breslow baseline estimate), and it cannot extrapolate beyond observed follow-up. Prefer parametric for HTA survival extrapolation (`survival-extrapolation-hta-rwe`) and lifetime cost-effectiveness, where you must project past the data.
• vs Poisson / negative-binomial rate models: Cox uses exact event-time ordering and handles smoothly varying baseline risk; Poisson assumes piecewise-constant rates and is simpler for aggregate person-time and recurrent counts. Prefer Poisson/NB (`poisson-negative-binomial-count-models`) for incidence-rate reporting and recurrent events when individual event timing is coarse.
• vs cause-specific / Fine-Gray competing-risks models: standard Cox censoring a competing event silently targets the cause-specific hazard and, if read as "risk," overstates the cumulative incidence of the event of interest whenever competing mortality differs by arm. Prefer the explicit competing-risks framing (`competing-risks-cause-specific-fine-gray-rwe`) in elderly, oncology, and end-stage populations.
• vs g-methods (MSM/IPTW-weighted Cox, g-estimation, clone-censor-weight): plain Cox cannot handle time-varying confounders affected by prior treatment, and conditioning on post-baseline mediators induces collider bias. Prefer g-methods (`marginal-structural-models-g-methods`, `g-estimation-structural-nested-models`, `clone-censor-weight-per-protocol`) for sustained-strategy or per-protocol estimands under time-varying confounding; Cox remains the analytic engine that runs after the weighting.

When NOT to use — and when it is actively misleading or dangerous

• Strong, structured PH violation (crossing or converging hazards). A reported single HR near 1 can hide an early harm that reverses to a late benefit. Do not paper over this; model time-by-covariate interaction, stratify the baseline hazard, switch to RMST, or report time-segmented HRs. Stensrud & Hernán (2020) argue the relevant question is usually not "does PH hold?" but "what estimand do I actually want?"
• Competing risks ignored. Censoring death to study a non-fatal endpoint in an elderly claims cohort, then interpreting 1 − KM as risk, overstates incidence when mortality is differential by arm — a classic, dangerous error. Use cause-specific and subdistribution analyses.
• Time-varying confounding by indication. When a post-baseline lab (eGFR, HbA1c) both responds to treatment and drives the next treatment decision and the outcome, standard Cox cannot remove the bias; adjusting for it opens a mediator/collider path. This is the g-methods boundary.
• Immortal time and misaligned time zero. If follow-up starts before exposure can occur (e.g., "ever had the procedure" with follow-up clocked from diagnosis), exposed person-time is inflated and the HR is biased toward benefit. Fix the design (new-user, time zero at the index event), not the model.
• Informative censoring. Loss to follow-up correlated with prognosis (sicker patients disenroll or switch) violates the independent-censoring assumption; require IPCW or a competing-risks/sensitivity framing (`attrition-and-loss-to-follow-up-rwe`).

Data-source operational depth

• Claims (FFS vs MA vs commercial): time zero is the index fill (`fill_date`) or, for procedures, the `service_date` on the CPT/HCPCS claim; outcomes are validated code algorithms (e.g., death via discharge status or a linked mortality index, MI via 1-IP or 2-OP rules). The dominant failure mode is enrollment-driven observability: Medicare Advantage encounter data are incomplete and MA-only person-time can lack the FFS claims that define both exposure and outcomes — restrict to enrollees with continuous Parts A/B/D (or full commercial medical+pharmacy) and exclude MA-only spans, or "no event" is unobserved censoring masquerading as survival. In elderly claims, competing risk of death differs by exposure, so a cause-specific-only Cox misleads. Procedure studies invite immortal time if the index date is not pinned to the actual CPT date; stockpiling and 90-day mail-order distort `days_supply` and any on-treatment window.
• EHR: the event is the order/administration, not the dispensing; link to fills to confirm initiation. Rich labs/vitals enable severity adjustment but arrive at irregular times — naive last-observation-carried-forward into a time-fixed Cox can bias, and a lab measured because the patient deteriorated is a collider. Visit-driven capture makes loss to follow-up informative; define observation windows explicitly.
• Registry: cleanest for staging, biomarkers, and adjudicated recurrence/progression (SEER-Medicare, disease-specific), but typically weak for complete pharmacy exposure and full mortality — link to claims and a death index. Often used to validate claims-based survival algorithms.
• Linked claims–EHR–vital records: the ideal substrate (severity + completeness + reliable mortality), but linkage selection and order/fill/service-date discrepancies must be reconciled before time-zero assignment, or immortal time creeps back in at the seam.

Worked claims example

Question: time to first hospitalized heart-failure event, second-generation sulfonylurea vs DPP-4 inhibitor, among adults with type 2 diabetes in a commercial + Medicare FFS database.

1. Eligibility: age ≥18, ≥2 diabetes diagnoses, and 365 days of continuous A/B/D (or commercial medical+pharmacy) enrollment with no MA-only spans before the first study fill.
2. Washout / new-user: no fill of any sulfonylurea or DPP-4 inhibitor in the 365-day lookback.
3. Time zero: the date of the first qualifying fill; assign `arm` from the dispensed NDC.
4. Outcome and time: `time_to_event` = days from `index_date` to the first validated HF hospitalization (1-IP primary-position algorithm); `event = 1` at that date.
5. Censoring: at disenrollment, end of data, and the competing event of non-HF death — analyze death as a competing risk, not as administrative censoring, and report both the cause-specific HR (from Cox) and the subdistribution HR / cumulative incidence.
6. Adjustment: fit `coxph(Surv(time_to_event, event) ~ arm + covariates)` on a propensity-matched set with baseline covariates measured only in [index_date − 365, index_date]; test PH with weighted Schoenfeld residuals (Grambsch & Therneau 1994), and if the HF hazards converge over follow-up, report RMST difference at 3 years alongside the HR rather than a single number.

Interpreting the output

A Cox model returns: adjusted HR = 0.75 (95% CI 0.60–0.94) for DPP-4 inhibitor vs sulfonylurea, time to first hospitalized heart-failure event, propensity-matched cohort.

Formal interpretation. The estimated hazard ratio of 0.75 means that, at any instant during follow-up, patients in the DPP-4 group who are still event-free have, on average, 75% of the instantaneous rate of heart-failure hospitalization observed in the sulfonylurea group, conditional on the matched covariates. Under the proportional-hazards assumption this ratio is assumed constant across follow-up time; if Schoenfeld-residual tests indicate time-varying hazards, the summary HR should be supplemented or replaced by a restricted mean survival time difference or a time-split analysis.

Practical interpretation. The data are consistent with a lower rate of heart-failure events on DPP-4 therapy, but the HR is not "25% fewer hospitalizations": it describes an instantaneous rate ratio among those still event-free, and its translation to absolute risk depends on the baseline hazard and length of follow-up. Pair this HR with an absolute measure — such as the RMST difference or cumulative incidence at 36 months — to communicate the magnitude of benefit in patient-meaningful terms.

Decision diagram

flowchart TD
  Q[Time-to-event question with a defensible time zero] --> PH{Proportional hazards<br/>plausible? Schoenfeld test}
  PH -->|Yes, PH holds| CR{Competing event<br/>differential by arm?}
  PH -->|No, hazards cross/converge| NPH[Stratify baseline OR time-by-covariate interaction<br/>report RMST / segmented HR]
  CR -->|No| TVC{Time-varying confounding<br/>affected by prior treatment?}
  CR -->|Yes| FG[Cause-specific Cox PLUS Fine-Gray / cumulative incidence]
  TVC -->|No| COX[Standard or PS-weighted Cox<br/>report HR + absolute summary]
  TVC -->|Yes| GM[g-methods: MSM-weighted Cox / g-estimation / clone-censor-weight]

flowchart LR
  U[Unmeasured healthy behaviors / frailty] -->|confounds| X[Treatment arm]
  U -->|drives event hazard| T[Event time]
  X --> T
  U -.->|depletes faster in higher-risk arm,<br/>making risk sets non-comparable| RS[Surviving risk set over time]
  RS -.->|induces time-varying observed HR<br/>even under constant true effect| T
  style U fill:#ffcccc

Worked example

Scenario

A researcher wants to know whether patients newly starting a DPP-4 inhibitor reach their first hospitalized heart-failure event more slowly than patients newly starting a second-generation sulfonylurea. Eight patients are enrolled on the day they pick up their first prescription (time zero). Four are in the DPP-4 arm and four are in the sulfonylurea arm. The study tracks each patient until either they are hospitalized for heart failure (the event) or their insurance coverage ends and we can no longer observe them (censored). We want to estimate the hazard ratio comparing the two arms.

Dataset

One row per patient. time_to_event is the number of days from first fill until hospitalization or last observable day. event = 1 means hospitalization occurred; event = 0 means the patient left follow-up before the event (censored).

FIG. 1 — DESIGN TIMELINE

Steps

Result

In the illustrative large-scale claims study described in this concept (type 2 diabetes patients, commercial + Medicare FFS database, DPP-4 inhibitor vs second-generation sulfonylurea, time to first hospitalized heart-failure event), the adjusted hazard ratio is 0.71 (95% CI 0.61–0.83), meaning patients on a DPP-4 inhibitor experienced hospitalized heart failure at 29% lower rate at any given point in follow-up compared to patients on a sulfonylurea, after adjusting for baseline differences between the groups.

Trade-offs

Runnable example

import pandas as pd
from lifelines import CoxPHFitter

# cohort: analysis-ready table described in the header (NO toy data created here).
cohort = pd.read_parquet("cohort.parquet")

covariates = ["arm", "age", "sex", "prior_hf", "ckd", "baseline_utilization"]
model_df = cohort[["time_to_event", "event"] + covariates].dropna()

cph = CoxPHFitter()
cph.fit(model_df, duration_col="time_to_event", event_col="event",
        robust=True)                      # robust SEs; required if weights/clustering are added later
print(cph.summary[["coef", "exp(coef)", "exp(coef) lower 95%",
                   "exp(coef) upper 95%", "p"]])   # exp(coef) = adjusted HR

# Proportional-hazards check via weighted Schoenfeld residuals (Grambsch & Therneau 1994).
# A small p-value for `arm` => PH violated for the exposure: stratify, add a time interaction,
# or report RMST / cumulative incidence instead of a single averaged HR.
cph.check_assumptions(model_df, p_value_threshold=0.05, show_plots=False)

library(survival)

# cohort: analysis-ready data.frame described in the header.
cohort <- readRDS("cohort.rds")
cohort$arm <- relevel(factor(cohort$arm), ref = "0")

fit <- coxph(
  Surv(time_to_event, event) ~ arm + age + sex + prior_hf + ckd + baseline_utilization,
  data   = cohort,
  ties   = "efron",        # Efron handling of tied event times (default; preferred over Breslow)
  robust = TRUE            # robust variance; required once IPTW/IPCW weights are added
)
summary(fit)               # exp(coef) = adjusted hazard ratio with 95% CI

# Weighted Schoenfeld PH test + diagnostic plot (Grambsch & Therneau 1994).
zph <- cox.zph(fit)
print(zph)                 # global + per-term p-values; p<0.05 for arm => PH violated
# If PH fails for arm: strata(stage) for nuisance violators, or add tt(arm) for a time-varying
# effect, or summarize with RMST: survival::rmean via summary(survfit(...), rmean=1095).

/* Cumulative incidence by arm (and competing-risk CIF if a competing event is coded). */
proc lifetest data=work.cohort plots=cif(test);
  time time_to_event*event(0);
  strata arm;
run;

/* Primary Cox model: adjusted (cause-specific) hazard ratio with PH diagnostics. */
proc phreg data=work.cohort;
  class arm (ref='0') sex (ref='0');
  model time_to_event*event(0) = arm age sex prior_hf ckd baseline_utilization / ties=efron rl;
  hazardratio 'Treatment effect' arm / diff=ref;
  assess ph / resample seed=123;          /* weighted Schoenfeld residual PH test (Grambsch-Therneau) */
  output out=work.resid ressch=sch_arm;   /* per-event Schoenfeld residuals for the exposure */
run;

/* Competing-risks (subdistribution) HR when a competing event exists (event=2 = competing death). */
proc phreg data=work.cohort;
  class arm (ref='0');
  model time_to_event*event(0) = arm age sex / eventcode=1 ties=efron;  /* Fine-Gray */
  hazardratio arm / diff=ref;
run;

Citations

1. [1]Cox DR. Regression models and life-tables. Journal of the Royal Statistical Society Series B (Methodological). 1972;34(2):187-220.
2. [2]Grambsch PM, Therneau TM. Proportional hazards tests and diagnostics based on weighted residuals. Biometrika. 1994;81(3):515-526.
3. [3]Stensrud MJ, Hernán MA. Why test for proportional hazards? JAMA. 2020;323(14):1401-1402.

1. [4]Hernán MA, Robins JM. Using big data to emulate a target trial when a randomized trial is not available. American Journal of Epidemiology. 2016;183(8):758-764.