Friday, August 31, 2012

Conditional survival times

Using this paper explaining how to generate conditional survival times and the survival library in R I created some example.

Here's an example for the Weibull distribution an a conditional T>2 event time against the non-conditional case.

The survival function is given by
$$ 1-F(t) = S(t) = \exp(-H(t))$$
So, by inverting the formulae we can generate event times using the cumulative hazard function
$$ T = H^{-1}(-\log U)$$
Now, substituting in the conditional hazard function

$$ S(t|t>T_k) = \exp(-a(H(t) - H(T_k)))$$
And rearranging as before gives
$$ T = H^{-1}(-{\log U}/{a} + H(T_k)) - T_k.$$

The code looks like this

library(survival)

n <- 10000
nbreak <- 100
lamb <- 2
r <- 1    # hazard adjustment
tmin <- 1    # conditional minimum event time e.g. intervention time
nu <- 2


# Weibull
invH_Wei <- function(t){
(1/lamb * t)^(1/nu)
}

T <- invH_Wei(-log(rv_unif)/r)


H_Wei <- function(t){
lamb * (t^nu)
}

T2 <- invH_Wei(-log(rv_unif)/r + H_Wei(tmin)) - tmin

Tcond <- T2 + tmin

plot(survfit(Surv(Tcond)~1))
lines(survfit(Surv(T)~1))

No comments:

Post a Comment