Subtleties in the interpretation of hazard contrasts

Abstract

The hazard ratio is one of the most commonly reported measures of treatment effect in randomised trials, yet the source of much misinterpretation. This point was made clear by Hernán (Epidemiology (Cambridge, Mass) 21(1):13–15, 2010) in a commentary, which emphasised that the hazard ratio contrasts populations of treated and untreated individuals who survived a given period of time, populations that will typically fail to be comparable—even in a randomised trial—as a result of different pressures or intensities acting on different populations. The commentary has been very influential, but also a source of surprise and confusion. In this note, we aim to provide more insight into the subtle interpretation of hazard ratios and differences, by investigating in particular what can be learned about a treatment effect from the hazard ratio becoming 1 (or the hazard difference 0) after a certain period of time. We further define a hazard ratio that has a causal interpretation and study its relationship to the Cox hazard ratio, and we also define a causal hazard difference. These quantities are of theoretical interest only, however, since they rely on assumptions that cannot be empirically evaluated. Throughout, we will focus on the analysis of randomised experiments.

Appendix A

1.1 A.1 Binary frailty variable

Let Z be binary, for instance $P(Z=0.2)=0.2$ and $P(Z=1.2)=0.8$, low and high risk groups, then $E(Z)=1$, and similar results as in the Gamma-distribution case (Sect. 3.2) are obtained. In this binary case, the Laplace transform is

$$\begin{aligned} \phi _Z(u)=0.2e^{-0.2u}+0.8e^{-1.2u}. \end{aligned}$$

Therefore

$$\begin{aligned} \text{ HR}_Z(t)\frac{g_Z(e^{-\Lambda (t,a=1)})}{g_Z(e^{-\Lambda (t,a=0)})} =\frac{\lambda (t;a=1)}{\lambda (t;a=0)}, \end{aligned}$$

where

$$\begin{aligned} g_Z(u)=\{D\log {(\phi _Z)}\}\{\phi _Z^{-1}(u)\}, \end{aligned}$$

(17)

which is an increasing function. Hence, if we take $\beta _1<0$ and $\beta _2=0$, then again

$$\begin{aligned} \text{ HR}_Z(t)<1 \end{aligned}$$

for all t.

1.2 A.2 Frailty model arising by marginalization

The DGP (4) given by $\lambda (t;A,Z)$ with $\nu =\infty$ induces the marginal Cox model for the observed data where we do not get access to Z. It is tempting to interpret $\text{ HR}_{Z}(t)$ as a causal hazard ratio, but this only holds under further untestable assumptions as shown in Sect. 4.1. Below, we formulate a more general DGP involving an additional frailty variable so that both the marginal Cox model and the model $\lambda (t;A,Z)$ are correctly specified.

Rename Z to $Z_1$. We then show that similarly we can also pick a DGP $\lambda (t;A,Z_1,Z_2)$ so that it marginalizes to $\lambda (t;A,Z_1)$ that further marginalizes to $\lambda (t;A)$, the latter being the Cox model. For ease of calculations, let

$$\begin{aligned} \lambda (t;A,Z_1,Z_2)=Z_2\lambda ^*(t;A,Z_1) \end{aligned}$$

with $Z_2$ being Gamma distributed with mean and variance equal to 1, and independent of $Z_1$ and A. Similar calculations as those in Sect. 3.1 gives the following expression

$$\begin{aligned} \lambda (t;A,Z_1,Z_2)=Z_1Z_2\lambda _0(t)e^{\beta A} \exp {\bigl [\Lambda _0(t)e^{\beta A}+Z_1\{ \exp {\{\Lambda _0(t)e^{\beta A}}-1\}\bigr ]} \end{aligned}$$

(18)

Hence, if the DGP is governed by (18) then model (4), with $\nu =\infty$, and the marginal Cox model, only conditioning on A, are also correctly specified.

1.3 A.3 Hazard differences

Assume (5) and (6). Then

$$\begin{aligned} \lim _{h\rightarrow 0}P(t\le T^1< t+h|T^0\ge t,T^1\ge t)=\psi (t) +\frac{E\{\omega (t,Z)e^{-2\int _0^t\omega (s,Z)\, ds}\}}{E\{e^{-2\int _0^t\omega (s,Z)\, ds}\}} \end{aligned}$$

and

$$\begin{aligned} \lim _{h\rightarrow 0}P(t\le T^0< t+h|T^0\ge t,T^1\ge t) =\frac{E\{\omega (t,Z)e^{-2\int _0^t\omega (s,Z)\, ds}\}}{E\{e^{-2\int _0^t\omega (s,Z)\, ds}\}} \end{aligned}$$

and therefore

$$\begin{aligned} \psi (t)= & {} \lim _{h\rightarrow 0}P(t\le T^1< t+h|T^0\ge t,T^1\ge t)\\&-\lim _{h\rightarrow 0}P(t\le T^0< t+h|T^0\ge t,T^1\ge t). \end{aligned}$$

1.4 A.4 $\text{ HR }(t)$ under copulas not obeying (11)

We outline now how to generate samples from the copula $C^{*}$ to gain further insight, Nelsen (2006), pp. 40–41,

• 1. Generate a standard uniformly distributed random variable u;

• 2. Set

  $$\begin{aligned} v = \left\{ \begin{array}{ll} u&{} \quad \text {if}\ u\le \xi { or}u\ge 1-\xi \\ 1-u &{} \quad \text {if}\ \xi<u<1-\xi \end{array}\right. \end{aligned}$$

• 3. Set $T^0=\Lambda _0^{-1}\{-\log {(u)}\}$ and $T^1=\Lambda _0^{-1}\{-\log {(v)}/\phi \}$ (assuming $\Lambda _0$ is strictly increasing).

As an illustration, we took $\Lambda _0(t)=t$, $\phi =2/3$ and $\xi =0.251$ and generated $n=2,000,000$ samples. For these simulated data we obtained an estimated Cox HR of 0.666 and Spearman’s $\rho$ was estimated to equal 0.753, which fits nicely with the theoretical values. As mentioned, if we follow the path $(\exp {\{-\Lambda _0(t)\}},\exp {\{-\phi \Lambda _0(t)\}})$ (starting at (1,1) when $t=0$) in $[0,1]^2$ in such a way that we pass through the upper triangle of $J_2$ then $\text{ HR }(t)$ will be equal to 0 before we reach $J_2$ but then it jumps to $k(u,v)\phi$, where $k(u,v)>1$. To illustrate this, we estimated $\text{ HR }(t)$ non-parametrically in time-points t until we cross the line $v=1-u$ in $J_2$. This was done by taking the observed frequency of $T^1$-failures in the time interval $[t,t+\delta ]$, $\delta =0.01$, among those where both $T^0\ge t$ and $T^1\ge t$; and similarly for $T^0$; and then taking the ratio of the two. This is seen in Fig. 3, where the theoretical values $\phi \exp {\{-t(\phi -1) \}}$ are superimposed. The dashed line corresponds to $\phi =2/3$. We see that the theoretical values fit nicely to the estimated ones. For this copula it can be shown that $\text{ HR }(t)$ is always below 1, but as we have just seen it can also be larger than $\phi$! Hence, to have that $\text{ HR }(t)<\phi$ we need more than just positive correlation between $T^0$ and $T^1$ as this example illustrates.

1.5 A.5 Monotonicity and $\text{ HR }(t)$

Assume that $T^1\ge T^0$ a.s. and that the distribution of $(T^0,T^1)$ is absolutely continuous with density $f(t_0,t_1)$ with respect to the Lebesgue measure. Monotonicity is attained if $f(t_0,t_1)=0$ for $t_0>t_1$. We now look at

$$\begin{aligned} z(h)&=P(T^0\ge t, t\le T^1\le t+h)\\&=P(T^0\ge t, t\le T^1\le t+h, T^0\le T^1)\\&=\int _t^{t+h}\int _t^{t_1}f(t_0,t_1)dt_0dt_1. \end{aligned}$$

Hence

$$\begin{aligned} \lim _{h\rightarrow 0}\frac{z(h)}{h}=\lim _{h\rightarrow 0}\frac{G(t+h) -G(t)}{h}=G^{\prime }(t)=\int _t^{t}f(t_0,t)dt_0=0, \end{aligned}$$

where

$$\begin{aligned} G(u)=\int _0^{u}\left\{ \int _t^{t_1}f(t_0,t_1)dt_0\right\} dt_1. \end{aligned}$$

But this means $\lambda _{T^1}(t|T^0\ge t,T^1\ge t)=0$, and therefore $\text{ HR }(t)=0$. To avoid having $\lambda _{T^1}(t|T^0\ge t,T^1\ge t)$ equal to zero we need for instance $T^0$ given $T^1$ to have a discrete distribution, but this seems unappealing.