Data Augmentation for Censored Survival Data

Suppose we have N individuals in the dataset where the time-to-event data for individual i is given by T_{i}. Of those N individuals, N_{obs} experience the event and N_{cen} are censored at \{C_{i}\}_{i=1}^{N_{cen}}. The vector representing the true event times, say \textbf{T}_{true} , can be split into \textbf{T} _{true}= [\textbf{T}_{obs,true}, \textbf{T}_{cen,true} ] where \textbf{T}_{obs,true} denotes the vector of (true ) death times for observed individuals and \textbf{T}_{cen,true} denotes the vectorof (true ) death times for censored individuals. This implies T_{i,obs,true} = T_{i} and T_{i,cen,true} \geq C_{i} .

In this simple example, suppose that the data comes from Weibull distribution (with shape \eta and scale \lambda). Our goal is augment the missing true death times for censored individuals. In other words, at each iteration q, q \in \{ 1, ..., Q\}, we wish to draw:

\textbf{T}_{i,cen,true}^{(q+1)} \sim p(\textbf{T}_{i,cen,true}| \textbf{T}_{i,cen,true} \geq C_{i}, \textbf{T}_{i,obs,true}, \eta^{(q)}, \lambda^{(q)})
\eta^{(q+1)} , \lambda^{(q+1)}\sim p(\eta , \lambda | \textbf{T}_{i,cen,true}^{(q+1)}, \textbf{T}_{i,obs,true})

Does the following code give us what we want? Please provide your thoughts and comments. Thanks!

data {
  int<lower=0> Nobs;                          // Number of individuals who experienced death before censoring
  int<lower=0> Ncen;                          // Number of censored individuals
  real<lower=0> ObsDeathTime[Nobs];           // Observed death times
  real<lower=0> ObsCenTime[Ncen];             // Observed Censoring times

parameters {
  real<lower = 0> eta;                           // shape parameter for ObsTime
  real<lower = 0> lambda;                        //  scale parameter for ObsTime
  vector<lower=0>[Ncen] T_miss_raw;             // vector for missing death times
transformed parameters{
   vector<lower=0>[Ncen] T_miss;
   for(j in 1:Ncen){
    T_miss[j] = T_miss_raw[j] + ObsCenTime[j];  // make sure that T_miss[j] is greater than or equal to C_i

model {
  // Priors
  eta ~ gamma(1,1);
  lambda ~ gamma(1,1);
  for(i in 1:Nobs){
    target+= weibull_lpdf(ObsDeathTime[i]|eta,lambda);    // observed death time
  for(j in 1:Ncen){
    target+= weibull_lpdf(T_miss[j]|eta,lambda);     // contribution of censored death time to the likelihood once we know what their values are