Parametric models (Chapters 17–22) require us to specify the shape of the baseline hazard \(h_0(t)\). In practice, we may not be confident in this choice, and a misspecified baseline can bias the estimated effects of the explanatory variables \(\hat{\bbeta}\).
In 1972, Sir David Cox proposed a semi-parametric approach that estimates \(\bbeta\) without specifying \(h_0(t)\) at all. The baseline hazard is treated as an unspecified nuisance function, and inference on \(\bbeta\) is based on a partial likelihood that eliminates \(h_0(t)\).
The explanatory-variable effects are the same as in the parametric PH model (Chapter 22). What is new is that \(h_0(t)\) is left completely unspecified.
Cox’s key insight was to construct a likelihood for \(\bbeta\) that does not involve \(h_0(t)\).
Let \(t_{(1)} < t_{(2)} < \cdots < t_{(r)}\) be the ordered event times (ties handled separately). At each event time \(t_{(j)}\), define the risk set \[ \mathcal{R}(t_{(j)}) = \{\text{individuals still at risk immediately before } t_{(j)}\}. \]
Given that an event occurs at \(t_{(j)}\), the conditional probability that it is individual \(i_{(j)}\) who experiences the event (rather than any other at-risk individual) is \[ \frac{h_0(t_{(j)})\exp(\bbeta^T\mathbf{x}_{i_{(j)}})} {\sum_{\ell \in \mathcal{R}(t_{(j)})} h_0(t_{(j)})\exp(\bbeta^T\mathbf{x}_\ell)} = \frac{\exp(\bbeta^T\mathbf{x}_{i_{(j)}})} {\sum_{\ell \in \mathcal{R}(t_{(j)})} \exp(\bbeta^T\mathbf{x}_\ell)}. \]
The baseline hazard \(h_0(t_{(j)})\) cancels in the numerator and denominator.
Proof: See the supplementary notes proof_partial_likelihood.pdf for a full justification of this conditioning argument.
This can be written compactly as \[ L(\bbeta) = \prod_{i=1}^n \left[ \frac{\exp(\bbeta^T\mathbf{x}_i)} {\sum_{\ell \in \mathcal{R}(t_i)} \exp(\bbeta^T\mathbf{x}_\ell)} \right]^{\delta_i}. \]
Only the \(\delta_i = 1\) terms contribute to the product.
\[ \ell(\bbeta) = \sum_{i=1}^n \delta_i \left\{ \bbeta^T\mathbf{x}_i - \log\sum_{\ell \in \mathcal{R}(t_i)}\exp(\bbeta^T\mathbf{x}_\ell) \right\}. \tag{23.1}\]
The maximum partial likelihood estimator (MPLE) \(\hat{\bbeta}\) maximises Equation 23.1. Under regularity conditions, \(\hat{\bbeta}\) has the same asymptotic properties as an MLE: \[ \hat{\bbeta} \;\dot\sim\; \mathrm{N}_p\!\left(\bbeta,\; \mathcal{I}(\bbeta)^{-1}\right), \] where \(\mathcal{I}(\bbeta) = -\partial^2\ell/\partial\bbeta\partial\bbeta^T\) is the observed information matrix.
The partial likelihood above assumes no ties. With ties, the most common approximation is Breslow’s approximation: \[ L(\bbeta) \approx \prod_{j=1}^r \frac{\exp(\bbeta^T\mathbf{s}_{(j)})} {\left(\displaystyle\sum_{\ell \in \mathcal{R}(t_{(j)})} \exp(\bbeta^T\mathbf{x}_\ell)\right)^{d_j}}, \] where \(d_j\) is the number of deaths at \(t_{(j)}\) and \(\mathbf{s}_{(j)}\) is the sum of the explanatory-variable vectors of those who die.
library(survival)
# Basic Cox model
cox_fit <- coxph(Surv(t, delta) ~ x1 + x2, data = mydata,
ties = "breslow")
summary(cox_fit)
# Output includes:
# - Estimated coefficients beta_hat
# - exp(coef): hazard ratios
# - 95% CI for hazard ratios
# - Wald, likelihood ratio, and score tests for overall significance
library(survival)
# Load data
h <- read.table("https://www.richardpmann.com/MATH3701/Datasets/hypernephroma.dat", header = TRUE)
h$age <- as.factor(h$age)
levels(h$age) <- c("<60", "60-70", ">70")
# Fit Cox model to nephrectomy patients only
h2 <- h[h$nephrectomy == 1, ]
cox_fit <- coxph(Surv(time, status) ~ age, data = h2, ties = "breslow")
summary(cox_fit)23.1. A small clinical trial compares Drug A (\(x=0\)) and Drug B (\(x=1\)). Survival times (days) are:
The Cox proportional hazards model is \(h(t, x) = e^{\beta x} h_0(t)\).
At each event time \(t_{(j)}\), the partial likelihood contribution is \(e^{\beta x_{(j)}} / \sum_{k \in R(t_{(j)})} e^{\beta x_k}\), where \(x_{(j)}\) is the value of the explanatory variable for the subject who had the event. The risk set \(R(t_{(j)})\) contains all subjects with \(t_k \geq t_{(j)}\). Multiply contributions over all event times.
Differentiate \(\ell(\beta) = \log L(\beta)\) with respect to \(\beta\) and set to zero.
23.2. Continuing from Exercise 23.1 with \(\hat\beta = -\log 3\).
Compute the approximate variance of \(\hat\beta\) using the observed information \(-\ell''(\hat\beta)\).
Test \(H_0: \beta = 0\) (drugs equally effective) against a two-sided alternative at the 5% significance level using a Wald test.
Compute \(\ell''(\beta) = \mathrm{d}^2\ell/\mathrm{d}\beta^2\) from \(\ell'(\beta) = 1 - 4e^\beta/(1+e^\beta)\). Evaluate at \(\hat\beta\) and take the negative reciprocal for the variance. The Wald statistic is \(\hat\beta/\widehat{\mathrm{se}}(\hat\beta)\).