rNIX/R/cox.R

#' COX hazards
#'
#' @param x data frame
#' @param vars a vector of variables used in Cox's survival analysis
#'
#' @return model as returned by coxph
#'
#' @export 

cox.hazard<-function(x,vars){
	vlist=base::paste(vars,collapse='+')
	f<-stats::as.formula(base::paste('Surv(followup,censored)',vlist,sep='~'))
   s=survival::coxph(f,data=x)
   s
}

#' COX for univariate variables
#'
#' @param x data frame containing followup and censored columns
#' @param var the variable used to explain survival in format mode:name where mode drives unitChange selection, SD for m, 1 for r or logic and X for uX
#'
#' @return data frame with a single row for tested variable, with fields variableName,unitChange,HR,l95,u95,pWald,pLikelihood
#'
#' @export 

cox.univariate<-function(x,var){
   data=base::c(var)
   qv=base::strsplit(data,split=':')[[1]]
   mode=qv[1]
   v=qv[2]
   s=cox.hazard(x,v)
   s1<-base::summary(s)  
  
#calculate hazard ratio per std(m), given unit (u) or unit change (

#take std
   if (mode=='m'){
      x.dx<-stats::sd(x[,v],na.rm=TRUE)
   }

#use whatever follows u
   if (base::substr(mode,1,1)=='u'){
      x.dx=base::as.numeric(base::substring(mode,2))
   }
   
#take unit change
   if (mode=='logic'){
      x.dx=1
   }
   if (substr(mode,1,1)=='r'){
      x.dx=1
   }

#adjust output data
   sd=base::sqrt(s$var)*x.dx
   c=s1$coef[1]*x.dx
   l95=base::round(base::exp(c-1.96*sd),2)
   u95=base::round(base::exp(c+1.96*sd),2)
   HR=base::round(base::exp(c),2)
   pWald=base::round(s1$coefficients[1,5],2)
   pLikelihood=base::round(s1$logtest['pvalue'],3)
   unitChange=base::round(x.dx,2)
#return as data frame with a single row
   base::data.frame(varName=v,unitChange=unitChange,
      HR=HR,low95=l95,up95=u95,pWald=pWald,pLikelihood=pLikelihood)
  
}