R/WCH3df.R

Defines functions WCH3df

Documented in WCH3df

#' Weighted cumulative hazard estimator for three gap times distribution 
#' function.
#' @description Provides estimates for three gap times distribution function 
#' based on Weighted cumulative hazard estimator (WCH).
#' @usage WCH3df(object, x, y, z)
#' @param object An object of class multidf. 
#' @param x The first time for obtaining estimates for the three gap times
#' distribution function.
#' @param y The second time for obtaining estimates for the three gap times 
#' distribution function.
#' @param z The third time for obtaining estimates for the three gap times 
#' distribution function.
#' @return Vector with the Weighted cumulative hazard estimates for three gap 
#' times distribution function.
#' 
#' @references Wang, M.C. and Wells, M.T. (1998). Nonparametric Estimation of 
#' successive duration times under dependent censoring, Biometrika 85, 561-572.
#' @seealso \code{\link{KMW3df}}, \code{\link{LIN3df}} and \code{\link{LDM3df}}.
#' 
#' @examples
#' data("bladder5state")
#' 
#' b4state<-multidf(gap1=bladder5state$y1, event1=bladder4state$d1, 
#'                  gap2=bladder5state$y2, event2=bladder4state$d2,
#'                  gap3=bladder5state$y3, status=bladder4state$d3)
#'                  
#' head(b4state)[[1]]
#' 
#' WCH3df(b4state, x=13, y=20, z=40)
#' 
#' b4<-multidf(gap1=bladder4$t1, event1=bladder4$d1, 
#'             gap2=bladder4$t2-bladder4$t1, event2=bladder4$d2,
#'             gap3=bladder4$t3-bladder4$t2, status=bladder4state$d3)
#'               
#' WCH3df(b4, x=13, y=20, z=40)
#' 
#' @author Gustavo Soutinho and Luis Meira-Machado

WCH3df <-
function(object, x, y, z)
{
 obj <- object[[1]]
 time1 <- obj$time1
 time2 <- obj$time2 - obj$time1
 time3 <- obj$time - obj$time2 
 time <- obj$time
 event1 <- obj$event1
 event2 <- obj$event2
 status <- obj$status

  Stime <- time[status == 0]
  Stime <- c(0, Stime, max(time))
  Stime <- sort(unique(Stime))
  cen.surv <- sapply(Stime, function (u) KM(time, 1 - status, u))	
  
  w <- sort(unique(time3[time3 <= z & time3 > 0]))
  nw <- length(w)
  dv <- rep(0, nw)
  dv2 <- rep(0, nw)
  
  pn <- lapply(w, function (u) which(time1 <= x & 
                                       time2 <= y & status == 1 & time3 == u))
  pd <- lapply(w, function (u) which(time1 <= x & 
                                       time2 <= y & event2 == 1 & time3 >= u))
  
  if (z == 0) res <- 0
  else {
  for (i in 1:nw) {
    m1 <- length(pn[[i]])
    m2 <- length(pd[[i]])
    if (m1 == 0) Gn <- 0
    else {
      Gn <- vector(length = m1)
      for (j in 1:m1) { p0 <- max(which(Stime <= time1[pn[[i]]][j] + 
                                          time2[pn[[i]]][j] + w[i]))
      Gn[j] <- 1/cen.surv[p0]}}
    
    if (m2 == 0) Gd <- 1
    else {
      Gd <- vector(length = m2)
      for (j in 1:m2) { p0 <- max(which(Stime <= time1[pd[[i]]][j] + 
                                          time2[pd[[i]]][j] + w[i]))
      Gd[j] <- 1/cen.surv[p0]}}
    
    dv[i] <- ifelse(sum(Gd) == 0, 1 , 1 - sum(Gn) / sum(Gd))
  }
    
    #db3 <- multidf(time1=obj$time1, event1=obj$event1,
    #               time=obj$time2,status=obj$event2)
    
    db3 <- multidf(gap1=obj$time1, event1=obj$event1,
                   gap2=obj$time2-obj$time1,status=obj$event2)
    
  res <- (1 - prod(dv)) * WCHdf(object=db3, x, y)
  }
  
  
  if (res < 0) res <- 0
  if (res > 1) res <- 1
  return(res)
}

Try the survivalREC package in your browser

Any scripts or data that you put into this service are public.

survivalREC documentation built on Aug. 9, 2023, 5:09 p.m.