Cox Proportional Hazards Model

Definition

Cox proportional hazards model assume that covariates affect the Hazard Function.

Let $T_1,T_2,\dots ,T_n$ be i.i.d. survival time, and $C_1,C_2,\dots ,C_n$ be i.i.d. censoring time. We can observe $(Y_i,{\delta }_i),i=1,2,\dots ,n$, where $Y_i=\min (T_i,C_i)$ and ${\delta }_i=I(T_i\le C_i)$ is censoring indicator, and have covariates ${\mathbf{x}}_i=(x_{i1},\dots ,x_{ip})$. Then, the Cox proportional hazards model is defined as $\lambda (t:\mathbf{x})={\lambda }_0(t)\exp ({\mathbf{x}}^{⊺}\mathbf{\beta })$ where ${\lambda }_0$ is called the baseline hazard function, i.e. hazard at $\mathbf{x}=0$

Conditional Likelihood

Let $Y_1<Y_2<\cdots <Y_n$ (no ties case), and ${\mathcal{R}}_i$ be the risk set. For each uncensored time $Y_i$, $P\{\text{a death in}[y_i,y_i+\Delta )|{\mathcal{R}}_i\}\approx \sum _{j\in {\mathcal{R}}_i}{\lambda }_0(y_j)\exp ({\mathbf{x}}_j^{⊺}\mathbf{\beta })\Delta$ Therefore, $P\{\text{a death of }i\text{ at time}{y}_i|\text{one death in }{\mathcal{R}}_i\text{ at time }{y}_i\}=\frac{{\lambda }_0(y_i)\exp ({\mathbf{x}}_{\mathbf{i}}^{⊺}\mathbf{\beta })\Delta }{\sum _{j\in {\mathcal{R}}_i}{\lambda }_0(y_i)\exp ({\mathbf{x}}_j^{⊺}\mathbf{\beta })\Delta }\approx \frac{\exp ({\mathbf{x}}_{\mathbf{i}}^{⊺}\mathbf{\beta })}{\sum _{j\in {\mathcal{R}}_i}\exp ({\mathbf{x}}_j^{⊺}\mathbf{\beta })}$ Taking the product of these conditional probabilities gives a conditional likelihood $L_C(\mathbf{\beta })={\prod }_{i:I_u}\frac{\exp ({\mathbf{x}}_{\mathbf{i}}^{⊺}\mathbf{\beta })}{\sum _{j\in {\mathcal{R}}_i}\exp ({\mathbf{x}}_j^{⊺}\mathbf{\beta })}={\prod }_{i=1}^n{\left[\frac{\exp ({\mathbf{x}}_{\mathbf{i}}^{⊺}\mathbf{\beta })}{\sum _{j\in {\mathcal{R}}_i}\exp ({\mathbf{x}}_j^{⊺}\mathbf{\beta })}\right]}^{{\delta }_i}$ where $I_u$ is the indicator set for uncensored samples.

The $L_C(\mathbf{\beta })$ is not a likelihood. However, Cox suggested treating the conditional likelihood as an ordinary likelihood to find the Maximum Likelihood Estimation.

Since there’s no analytic solution for the MLE, iterative methods such as Newton–Raphson method is used to estimate the coefficient $\mathbf{\beta }$.

The hazard ratio $\exp ({\hat{\beta }}_j)$ represents the relative change in Hazard Rate for a one-unit increase in the covariate $x_j$.

Goodness-of-Fit Test

For testing the null hypothesis $H_0:\mathbf{\beta }=0$, Cox suggested the Rao Test. ${\left(\frac{\partial \ln L_C(0)}{\partial \mathbf{\beta }}\right)}^{⊺}{\left(-\frac{{\partial }^2\ln L_C(0)}{\partial {\mathbf{\beta }}^2}\right)}^{-1}\left(\frac{\partial \ln L_C(0)}{\partial \mathbf{\beta }}\right)\stackrel{a}{\sim }{\chi }^2(p)$

Asymptotic Normality of MLE

$\hat{\mathbf{\beta }}\stackrel{D}{\to }{N}_p(\beta ,I^{-1}(\mathbf{\beta }))$ where $I(\mathbf{\beta })$ is the observed Fisher Information

Estimation of Survival Function

Under the Cox proportional hazards model, $S(t;\mathbf{x})=S_0(t{)}^{\exp ({\mathbf{x}}^{⊺}\mathbf{\beta })}$ To estimate $S(t;\mathbf{x})$, we can use $\hat{\mathbf{\beta }}$ for $\mathbf{\beta }$ but we still need to estimate $S_0(t)$, ${\lambda }_0(t)$, or ${\Lambda }_0(t)$.

Breslow suggested the estimators of ${\lambda }_0(t)$ and $S_{0,B}(t)$ as ${\hat{\lambda }}_{0,B}=\frac{1}{(y_{u_i}-y_{u_{i-1}})\sum _{j\in {\mathcal{R}}_i}\exp ({\mathbf{x}}_j^{⊺}\hat{\mathbf{\beta }})}$ If $Y_{u_{i-1}}<t<Y_{u_i}$

${\hat{S}}_{0,B}(t)={\prod }_{i:y_i\le t}\left(1-\frac{{\delta }_i}{\sum _{j\in {\mathcal{R}}_i}\exp ({\mathbf{x}}_j^{⊺}\hat{\mathbf{\beta }})}\right)$ If $\beta =0$, then ${\hat{S}}_{0,B}$ is the Kaplan-Meier Estimator

It has a few drawbacks

• ${\hat{S}}_{0,B}(t)\ne \exp (-{\hat{\Lambda }}_0(t))$
• ${\hat{S}}_{0,B}(t)$ can take negative values.

Tsiatis suggested a non-negative version of ${\hat{S}}_{0,B}(t)$ ${\hat{S}}_{0,T}(t)=\exp (-{\hat{\Lambda }}_{0,T}(t))$ where ${\hat{\Lambda }}_{0,T}(t)=\sum _{y_i\le t}\frac{{\delta }_i}{\sum _{j\in {\mathcal{R}}_i}\exp ({\mathbf{x}}_j^{⊺}\hat{\mathbf{\beta }})}$

Link suggested using the linear smooth of ${\hat{\Lambda }}_{0,T}(t)$.

Discrete on Grouped Data

When data is discrete or grouped, there are ties at each failure. Denote the ordered discrete failure time by $Y_1<Y_2<\cdots <Y_r$ and let ${\mathcal{R}}_i$ be the risk set at $Y_i^{-}$, ${\mathcal{D}}_i$ be the death set at $Y_i$, and $d_i=\#(D_i)$.

Cox suggested combining the all possible permutations. However, it is computationally infeasible. $L_C={\prod }_{i=1}^r\left[\frac{\prod _{j\in {\mathcal{D}}_i}{\psi }_j}{\sum _{{\mathcal{D}}_i^{\ast }}\prod _{j\in {\mathcal{D}}_i^{\ast }}{\psi }_j}\right]$ where ${\psi }_j=\exp ({\mathbf{x}}_j^{⊺}\mathbf{\beta })$, and ${\mathcal{D}}_i^{\ast }$ is the size $d_i$ subset of ${\mathcal{R}}_i$

Peto suggested an alternative likelihood that instead of all possible permutations, use the same contribution. $L_C={\prod }_{i=1}^r\left[\frac{\prod _{j\in {\mathcal{D}}_i}{\psi }_j}{(\sum _{j\in {\mathcal{R}}_i}{\psi }_j{)}^{d_i}}\right]$

Time Dependent Covariates

In the case, the covariate depends on time. We observe ${\mathbf{x}}_i(t)$ and the conditional likelihood defined as $L_C(\mathbf{\beta })={\prod }_{i:I_u}\frac{\exp ({\mathbf{x}}_i(y_i{)}^{⊺}\mathbf{\beta })}{\sum _{j\in {\mathcal{R}}_i}\exp ({\mathbf{x}}_j(y_i{)}^{⊺}\mathbf{\beta })}$ where $I_u$ is the indicator set for uncensored samples, and ${\mathcal{R}}_i$ is the risk set.

Facts

Any two individuals have hazard functions that are constant multiples of the one another.

The Survival Function of the Cox proportional hazard model is a family of Lehmann alternatives. $\forall S\in \mathcal{S},\exists S_0\in \mathcal{S}s.t.S=S_0^{\gamma }$ where $\gamma \in {\mathbb{R}}^{+}$.

If $p=1$, $x_i=I(i\in I_1)$, where $I_1$ is the indicator set for sample 1, and there are no ties, then Cox test is exactly equal to the Mantel-Haenszel Test.