adding files

README.md

@@ -1,2 +1,50 @@
-# survival_models
-A repository of code and resources for survival models
+Examples of plotting survival curves and time to event analysis
+
+
+* Fitting a Kaplan-Meier and a Cox proportional hazards model
+
+	* survival_analysis_example.Rmd
+
+* Fitting a piecewise exponential proportional hazards model
+
+	* piecewise_exponential_proportional_hazards_model.R
+
+* Fitting a survival model using a GLM Poisson model
+
+	* simple_glm_survival.R
+
+* Fitting survival trees
+
+	* survival_trees.R
+
+* Cheatsheets and other resources
+
+	* cheatsheets folder
+
+* Adapted from 
+
+
+     https://stat.ethz.ch/R-manual/R-devel/library/survival/html/survfit.formula.html
+
+
+* Adapted from
+
+
+    https://rviews.rstudio.com/2017/09/25/survival-analysis-with-r/
+
+
+* Adapted from 
+
+
+    https://stats.stackexchange.com/questions/348922/piecewise-exponential-model-how-to-fit
+
+
+* Additional resources
+
+
+    www.stat.columbia.edu/~madigan/W2025/notes/survival.pdf
+
+    http://data.princeton.edu/wws509/notes/
+
+    https://rpkgs.datanovia.com/survminer/index.html
+

piecewise_exponential_proportional_hazards_model.R

@@ -0,0 +1,53 @@
+##################################################################################################
+# Piecewise exponential proportional hazards survival model
+# 
+# Adapted from:
+# https://stats.stackexchange.com/questions/348922/piecewise-exponential-model-how-to-fit
+#
+# Theory:
+#  http://data.princeton.edu/wws509/notes/
+#  
+# Installation:
+#     install.packages("pammtools")
+#
+##################################################################################################
+
+##################################################
+# Load library
+##################################################
+library(pammtools)
+library(survival)
+#library(ggp)
+
+#########################
+# load example data
+#########################
+data(tumor)
+head(tumor)
+
+##################################################
+## transform to piece-wise exponential data (PED)
+##################################################
+ped_tumor <- as_ped(Surv(days, status) ~ ., data = tumor, cut = seq(0, 3000, by = 100))
+# the interval column indicates the interval, ped_status indicates the status
+# in the respective interval
+filter(ped_tumor, id == 1) %>% head() %>% select(1:8)
+
+###########################################################################
+# Fit the Piece-wise-exponential Additive Model (PAM) instead of PEM
+###########################################################################
+pam <- mgcv::gam(ped_status ~ s(tend) + age + sex + complications, 
+                 data = ped_tumor, family = poisson(), offset = offset)
+summary(pam)
+plot(pam)
+
+#########################
+# compare to Cox 
+#########################
+cox <- coxph(Surv(days, status)~ age + sex + complications, data=tumor)
+cbind(coef(pam)[2:4], coef(cox))
+
+cox_fit <- survival::survfit(cox)
+plot(cox_fit)
+
+

simple_glm_survival.R

@@ -0,0 +1,105 @@
+#############################################################################
+# Simple code to perform survival analysus using GLM Poisson model
+#
+# Adapted from:
+#   https://data.princeton.edu/wws509/r/c7s1
+#   https://rdrr.io/github/datashield/dsBaseClient/man/ds.glmSLMA.html
+#
+# Other resources:
+#   why intercept is taken from Poisson GLM model:
+#   http://sphweb.bumc.bu.edu/otlt/MPH-Modules/BS/R/R7_LogisticRegression-Survival/R7_LogisticRegression-Survival_print.html
+#   https://stats.stackexchange.com/questions/8117/does-cox-regression-have-an-underlying-poisson-distribution/8118#8118
+#   https://stats.stackexchange.com/questions/115479/calculate-incidence-rates-using-poisson-model-relation-to-hazard-ratio-from-cox
+#   https://bmcmedresmethodol.biomedcentral.com/articles/10.1186/1471-2288-12-34
+#   https://data.princeton.edu/wws509/r/c7s1
+# 
+#############################################################################
+
+####################
+# Load library
+####################
+library(survival)
+library(metafor)
+library(ggplot2)
+library(survminer)
+
+
+#######################
+# Load data
+#######################
+data("veteran")
+head(veteran)
+
+#########################################
+# Build Cox-proportional hazards model
+#########################################
+surv_object <- survival::Surv(time = veteran$time, event = veteran$status)
+cxph <- survival::coxph(formula = surv_object ~ trt + diagtime + age, data = veteran)
+summary(cxph)
+survminer::ggforest(model = cxph, data = veteran)
+
+# calculate HR and se
+# TODO: compute a 2x2 table for use in metafor
+#     https://www.openepi.com/TwobyTwo/TwobyTwo.htm
+summary(cxph)$coef
+
+# get hazard ratio and se
+hr_1 = summary(cxph)$coef[1,1]
+se_1 = summary(cxph)$coef[1,2]
+
+
+
+
+#####################################
+# Repeat using glm() Poisson model
+#####################################
+
+veteran$logsurv <- log(veteran$time)
+
+glm_hazard <- glm(formula = status ~ trt + diagtime + age,  # +  offset(veteran$logsurv),
+                  family = 'poisson', 
+                  offset = veteran$logsurv, 
+                  data = veteran
+                  )
+
+# Poisson model hazard ratio
+#	the intercept is the hazard ratio
+# http://sphweb.bumc.bu.edu/otlt/MPH-Modules/BS/R/R7_LogisticRegression-Survival/R7_LogisticRegression-Survival_print.html
+# https://stats.stackexchange.com/questions/8117/does-cox-regression-have-an-underlying-poisson-distribution/8118#8118
+
+summary(glm_hazard)$coef
+summary(glm_hazard)$coef[1,1]
+summary(glm_hazard)$coef[1,2]
+
+exp(summary(glm_hazard)$coef[1,1])
+exp(summary(glm_hazard)$coef[1,2])
+
+# Compare to Cox model hazard ratios
+summary(cxph)$coef
+summary(cxph)$coef[1,1]
+summary(cxph)$coef[1,2]
+
+
+# https://data.princeton.edu/wws509/r/c7s1
+pchisq(deviance(glm_hazard), glm_hazard$df.residual, lower.tail=FALSE)
+
+############################################
+# calculate hazard ratios
+############################################
+# hazard ratios from Poisson model
+exp(coef(glm_hazard[1])) # - 1
+# hazard ratios form Cox proportional hazards model
+exp(coef(cxph))
+
+
+# Other resources and explanations
+# why intercept:
+# http://sphweb.bumc.bu.edu/otlt/MPH-Modules/BS/R/R7_LogisticRegression-Survival/R7_LogisticRegression-Survival_print.html
+# https://stats.stackexchange.com/questions/8117/does-cox-regression-have-an-underlying-poisson-distribution/8118#8118
+# https://stats.stackexchange.com/questions/115479/calculate-incidence-rates-using-poisson-model-relation-to-hazard-ratio-from-cox
+# https://bmcmedresmethodol.biomedcentral.com/articles/10.1186/1471-2288-12-34
+# https://data.princeton.edu/wws509/r/c7s1
+
+
+
+