Skip to contents

Introduction

In this vignette we derive the explicit formula for the overall survival under piecewise constant hazards. This is implemented in the function PWCsurvOS().

So we are in the situation where the hazards for the transitions are piecewise constant, i.e. all three hazard functions λ01(t)\lambda_{01}(t) (stable to progression), λ02(t)\lambda_{02}(t) (stable to death) and λ12(t)\lambda_{12}(t) (progression to death) are step functions. Say the start points of the constant hazard pieces for λ01(t)\lambda_{01}(t) are 0t01(1)<<t01(k01)0 \equiv t_{01}^{(1)} < \dotsb < t_{01}^{(k_{01})}, k011k_{01} \geq 1, with corresponding constant positive hazards h01(1),,h01(k01)h_{01}^{(1)}, \dotsc, h_{01}^{(k_{01})}. Obviously we use here the smallest set of pieces, i.e. neighboring hazards are required to be different, h01(j)h01(j+1)h_{01}^{(j)} \neq h_{01}^{(j+1)}. This holds analogously for the hazard functions of the other two state transitions.

We denote the cumulative hazards similarly as Λ01(t)\Lambda_{01}(t), Λ02(t)\Lambda_{02}(t) and Λ12(t)\Lambda_{12}(t). Note that these are piecewise linear, with the slope changes occurring at the times of hazard changes.

Overall survival calculation

Now we want to calculate the overall survival (OS) survival function induced by the piecewise constant hazard model. We start from

SOS(t)=SPFS(t)+0tSPFS(u)λ01(u)exp(Λ12(u)Λ12(t))du S_{\text{OS}}(t) = S_{\text{PFS}}(t) + \int_0^t S_{\text{PFS}}(u) \lambda_{01}(u)\exp(\Lambda_{12}(u) - \Lambda_{12}(t))\, du where SOS(t)S_{\text{OS}}(t) is the survival function for OS, and SPFS(t)S_{\text{PFS}}(t) is the survival function for PFS with the closed form

SPFS(t)=exp(Λ01(t)Λ02(t)). S_{\text{PFS}}(t) = \exp(- \Lambda_{01}(t) - \Lambda_{02}(t)). Hence we can rewrite the integral from above as

exp(Λ12(t))0texp(Λ12(u)Λ01(u)Λ02(u))λ01(u)du \exp(- \Lambda_{12}(t)) \int_0^t \exp(\Lambda_{12}(u) - \Lambda_{01}(u) - \Lambda_{02}(u))\lambda_{01}(u)\, du So overall we now have

$$ S_{\text{OS}}(t) = S_{\text{PFS}}(t) + \exp(- \Lambda_{12}(t)) \cal{I}(t) $$ and we can rewrite the integral $$ \cal{I}(t) := \int_0^t \exp(\Lambda_{12}(u) - \Lambda_{01}(u) - \Lambda_{02}(u))\lambda_{01}(u)\, du $$ in terms of the unique starting time points 0t(1)<<t(k)0 \equiv t_{(1)} < \dotsb < t_{(k)}, chosen such that the set {t(1),,t(k)}\{t_{(1)}, \dotsc, t_{(k)}\} is the smallest super set of all state specific transition starting points {t01(1),,t01(k01)}\{t_{01}^{(1)}, \dotsc, t_{01}^{(k_{01})}\}, {t02(1),,t02(k02)}\{t_{02}^{(1)}, \dotsc, t_{02}^{(k_{02})}\} and {t12(1),,t12(k12)}\{t_{12}^{(1)}, \dotsc, t_{12}^{(k_{12})}\}, as

$$\begin{align} \cal{I}(t) = &\int_{t_{(1)}}^{t_{(2)}} \exp(a_{(1)} + b_{(1)}(u - t_{(1)}))h_{01(1)}\,du \\ &+ \dotsb + \\ &\int_{t_{(l)}}^{t} \exp(a_{(l)} + b_{(l)}(u - t_{(l)}))h_{01(l)}\,du, \end{align}$$ where:

  • t(l)t_{(l)} is the start of the last integral part and defined as the maximum starting point that is smaller than tt
  • h01(j)=λ01(t(j))h_{01(j)} = \lambda_{01}(t_{(j)}), j=1,,kj=1, \dotsc, k are the hazard values for the stable to progression transition within the unique time intervals
  • a(j)=Λ12(t(j))Λ01(t(j))Λ02(t(j))a_{(j)} = \Lambda_{12}(t_{(j)}) - \Lambda_{01}(t_{(j)}) - \Lambda_{02}(t_{(j)}), j=1,,kj=1, \dotsc, k are the intercepts
  • b(j)=h12(j)h01(j)h02(j)b_{(j)} = h_{12(j)} - h_{01(j)} - h_{02(j)}, j=1,,kj=1, \dotsc, k are the slopes, again using the hazard values within the unique time intervals for the specific transitions

Note that this is essentially just because of Λ(t)=Λ(s)+h(ts) \Lambda(t) = \Lambda(s) + h (t-s) when there is a constant hazard λ(t)h\lambda(t) \equiv h and two time points s<ts<t.

We can then easily derive a closed form for each integral part, j=1,,lj = 1, \dotsc, l, where t(l+1)tt_{(l+1)}\equiv t is the end point of the last integral:

$$\begin{align} \cal{I}_{j} &= \int_{t_{(j)}}^{t_{(j+1)}} \exp(a_{(j)} + b_{(j)}(u - t_{(j)}))h_{01(j)}\,du\\ &= \exp(a_{(j)} - b_{(j)}t_{(j)})h_{01(j)} \int_{t_{(j)}}^{t_{(j+1)}} \exp(b_{(j)}u)\,du\\ &= \exp(a_{(j)} - b_{(j)}t_{(j)})h_{01(j)} \left[ b_{(j)}^{-1}\exp(b_{(j)}u) \right]_{t_{(j)}}^{t_{(j+1)}}\\ &= \frac{h_{01(j)}}{b_{(j)}} \exp(a_{(j)} - b_{(j)}t_{(j)}) (\exp(b_{(j)}t_{(j+1)}) - \exp(b_{(j)}t_{(j)}))\\ &= \frac{h_{01(j)}}{b_{(j)}} \exp(a_{(j)}) (\exp(b_{(j)}(t_{(j+1)} - t_{(j)})) - 1) \end{align}$$

Note that if it should happen that b(j)=0b_{(j)} = 0, the integral simplifies further to $$ \cal{I}_{j} = \int_{t_{(j)}}^{t_{(j+1)}} \exp(a_{(j)})h_{01(j)}\,du = \exp(a_{(j)})h_{01(j)}(t_{(j+1)} - t_{(j)}). $$

Implementation

The above formula is implemented in the function PWCsurvOS(). Note that there are a few modifications compared to above exposition:

  1. In order to be efficient, we look at a vector of time points tt directly from the beginning. Then, find the unique time points and order them, and perform piecewise integration between the sorted unique time points.
  2. In order to avoid numerical overflow with large time points, we move the term exp(Λ12(t))\exp(- \Lambda_{12}(t)) inside the integral, similarly as in the starting point of the derivation.