R/pem_surv.R

In csetraynor/rstanhaz: Bayesian Survival Models (rstanhaz)

#' Piecewise exponential
#'
#'
#' Semiparametric GAM Poisson model
#'
#' @param x data
#' @param surv_mod vector with covariates
#' @return mod
#' @seealso [coxph]
#' @keywords coxph
#'
#' @author Carlos S Traynor
#' @references
#'
#'  Terry M. Therneau and Patricia M. Grambsch (2000).
#'   Modeling Survival Data: Extending the Cox Model.
#'   Springer, New York. ISBN 0-387-98784-3.
#'@export pem_surv
pem_surv <- function(x, status = "status", time = "time", surv_form = NULL, prior = NULL, form = NULL,  ...) {
  
  if(!is.null(surv_form)) {
    input.formula <- as.formula( paste0(c( "status~1+offset(log_dtime)+s(tend)+", paste(surv_form, collapse = "+")), collapse = ""))
  } else{
    if(!is.null(form) ) {
      input.formula <- form
    } else {
      input.formula <- as.formula( "status~1+offset(log_dtime)+s(tend)")
    }
  }
  
  #prepare long format dataset
  long_x <- gen_long_dat(x, status = status, time = time)
  
  mgcv::gam( input.formula, 
            data = long_x,
            family=poisson() )
}


#' Additive hazard 
#'
#'
#' Additive risk model
#'
#' @param x data
#' @param surv_mod vector with covariates
#' @return mod
#' @seealso [coxph]
#' @keywords coxph
#'
#' @author Carlos S Traynor
#'@export ahaz_surv
ahaz_surv <- function(x, status = "status", time = "time", surv_form = NULL,  ...) {

  surv <- survival::Surv(x[[time]], x[[status]])
  
  X <- as.matrix(gam_matrix(x = x, vars = surv_form))
  
  
  fit <- ahaz::ahaz(surv, X)
  out <- list(fit, surv_form)
  names(out) <- c("mfit", "formula")
  out
}




#' Make GAM matrix
#'
#'
#' For PEM and Additive risk model
#'
#' @param x data
#' @param vars variables in vector form.
#' @return mod
#' @keywords coxph
#'
#' @author Carlos S Traynor
#'@export gam_matrix
gam_matrix <- function(x, vars){
  
  varis.obs <- gsub("factor|\\(|\\)|as.", "", vars)
  
  x <- x[varis.obs]
  
  X <- data.frame(id = seq_along(1: nrow(x)))
  
  for(i in seq_along(vars) ){
    
    if(is.numeric( x[[i]] )){
      matrixna <- x[i]
      X <- cbind(X, matrixna)
    } else{
      x[[i]] <- as.factor(x[[i]])
      coluna <- x[i]
      matrixna <- model.matrix(~ ., coluna)[ ,-1]
      matrixna <- as.data.frame(matrixna)
      if('matrixna' %in% colnames(matrixna) ){
        colnames(matrixna) <- paste0(colnames(coluna), "1")
      }
      X <- cbind(X, matrixna)
    }
  }
  
  X <- X[-match("id", colnames(X) )]
  X
}