library(readr)
library(mstate)
library(dplyr)
library(survival)
options(kableExtra.latex.load_packages = FALSE)
library(kableExtra)
library(lubridate)
library(tableone)
library(papeR)
library(rmarkdown)
library(Phase2.1MultStateModPackage)



# Prepare Data
my.data <- prep_data(params$dir_LPS ,
                     ymd(params$admin_censor))

# Set days from admission for length of stay estimates

end_limit <- 50

# Set number of bootstrap samples for confidence interval estimates

num_boot <- 1000


# Set x limit for plots

end_limit_plot <- 50

# Set Parameters for plot function

col_vec <- c("khaki1", "indianred1", "cornflowerblue","gray")

# Load data
#load(my.data)




# Set transition matrix for mstate
tra <- transMat(x = list(c(2,3,4), c(3,4), c(),c()),
                names = c("Non-Severe", "Severe", "Discharge","Death"))

# Add transition vectors
my.data$trans <- NA
for (i in 1:nrow(tra))
{
  for (j in 1:ncol(tra))
  {
    my.data$trans[which(my.data$from == i & my.data$to ==j)] <- tra[i, j]  
  }
}

# Create data frame with all possible transitions for when a patient is at risk
my.data_ext <- ext_mstate(my.data, tra)

# Create a data frome for table one and competing risks analyses
my.data_CR <- prep_data_CR(my.data)
statement1 <- c(paste("Site =", params$site_loc))
cat(statement1)
statement2 <- c(paste("Administrative censoring on", params$admin_censor))
cat(statement2)

Table One

my.data_table1 <- my.data_CR

my.data_table1$event <- factor(ifelse(my.data_table1$event == 0,  "censored", ifelse(my.data_table1$event == 1, "died", "discharged")),  ordered = FALSE)
#my.data_table1$event <- relevel(my.data_table1$event, ref="0")

my.data_table1$wave <- factor(my.data_table1$wave, labels = c("1: (01.01.20 - 31.07.20)", "2: (01.08.20 - 30.06.21)", "3: (01.07.21 - )"),  ordered = FALSE)

my.data_table1$age_cat <- relevel(my.data_table1$age_cat, ref="12-25")
my.data_table1$age_cat <- relevel(my.data_table1$age_cat, ref="00-11")

vars <- c("event", "sex", "age_cat", "wave")

tableOne <- CreateTableOne(vars = vars, data = my.data_table1)

kableone(tableOne, format = 'latex',  booktabs = T,linesep = c( "\\hline", "", "", "", "\\hline", "\\hline","", "", "", "", "", "", "\\hline", "", "", ""   ))

Full Cohort

Plots of Predicted Proportions

# Analysis -----------------------------------------------------------------------

# calculate transition probabilities
pt_full <- prob.tran.mat(my.data_ext, !is.na(my.data_ext$id), tra)

# ELOS and probability of long stayer
LOS_mat_full <- ELOS_mat(pt_full[[1]], end_limit, !is.na(my.data_ext$id), num_boot)

# LOS_mat_full[[1]] %>%
#   kbl() %>%
#   kable_styling()
# plot of transition probabilities for patients starting normal ward
plot(pt_full[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25,  
      main= "Patients Starting in Non-Severe")

# # plot of transition probabilities for patients starting normal ward
# plot(pt_full[[1]], from = 2,ord = c(4,2,1,3),type= "filled", cols = col_vec,
#      lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25,  legend = c("","Severe","Discharge","Death"),
#  main= "Patients Starting in Severe")

Full Cohort

Length of stay estimates

LOS_mat_full[[1]] %>%
  kbl(booktabs = T, linesep = c("", "", "\\hline", "", "")) %>%
  kable_styling(latex_options = c("hold_position"))

\newpage

Sex

Plots of Predicted Proportions

# calculate transition probabilities
pt_male <- prob.tran.mat(my.data_ext,my.data_ext$sex == 'male', tra)

# ELOS and probability of long stayer
LOS_mat_male <- ELOS_mat(pt_male[[1]], end_limit, my.data_ext$sex == 'male', num_boot)


# calculate transition probabilities
pt_female <- prob.tran.mat(my.data_ext,my.data_ext$sex == 'female', tra)

# ELOS and probability of long stayer
LOS_mat_female <- ELOS_mat(pt_female[[1]], end_limit, my.data_ext$sex == 'female', num_boot)


LOS_mat_sex <- cbind(LOS_mat_male[[1]], LOS_mat_female[[1]])

colnames(LOS_mat_sex) <- c("Male", "Female")

# LOS_mat_sex %>%
#   kbl() %>%
#   kable_styling()
par(mfrow=c(1,2))

# plot of transition probabilities for patients starting normal ward
plot(pt_male[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25,  
 main= "Male Patients")

# plot of transition probabilities for patients starting normal ward
plot(pt_female[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25,  
   main= "Female Patients")

Sex

Length of stay estimates

LOS_mat_sex %>%
  kbl(booktabs = T, linesep = c("", "", "\\hline", "", "")) %>%
  kable_styling(latex_options = c("hold_position"))

Age Categories

Plots of Predicted Proportions

# calculate transition probabilities
p_age1 <- prob.tran.mat(my.data_ext,my.data_ext$age_cat == "00-11", tra)

# ELOS and probability of long stayer
LOS_mat_age1 <- ELOS_mat(p_age1[[1]], end_limit, my.data_ext$age_cat == "00-11",
                         num_boot)
# calculate transition probabilities
p_age2 <- prob.tran.mat(my.data_ext,my.data_ext$age_cat == "12-25", tra)

# ELOS and probability of long stayer
LOS_mat_age2 <- ELOS_mat(p_age2[[1]], end_limit, my.data_ext$age_cat == "12-25", 
                         num_boot)
# calculate transition probabilities
p_age3 <- prob.tran.mat(my.data_ext,my.data_ext$age_cat == "26-49", tra)


# ELOS and probability of long stayer
LOS_mat_age3 <- ELOS_mat(p_age3[[1]], end_limit, my.data_ext$age_cat == "26-49",
                         num_boot)
# calculate transition probabilities
p_age4 <- prob.tran.mat(my.data_ext,my.data_ext$age_cat == "50-69", tra)

# ELOS and probability of long stayer
LOS_mat_age4 <- ELOS_mat(p_age4[[1]], end_limit, my.data_ext$age_cat == "50-69",
                         num_boot)
# calculate transition probabilities
p_age5 <- prob.tran.mat(my.data_ext,my.data_ext$age_cat == "70-79", tra)


# ELOS and probability of long stayer
LOS_mat_age5 <- ELOS_mat(p_age5[[1]], end_limit, my.data_ext$age_cat == "70-79",
                         num_boot)
# calculate transition probabilities
p_age6 <- prob.tran.mat(my.data_ext,my.data_ext$age_cat == "80+", tra)

# ELOS and probability of long stayer
LOS_mat_age6 <- ELOS_mat(p_age6[[1]], end_limit,my.data_ext$age_cat == "80+",
                         num_boot)
LOS_mat_age <- cbind(LOS_mat_age1[[1]], LOS_mat_age2[[1]], LOS_mat_age3[[1]], LOS_mat_age4[[1]], LOS_mat_age5[[1]], LOS_mat_age6[[1]])

colnames(LOS_mat_age) <- c("00-11", "12-25", "26-49", "50-69", "70-79", "80+")

# LOS_mat_age %>%
#   kbl() %>%
#   kable_styling()
par(mfrow=c(1,2))
plot(p_age1[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25, main= "Age 0-11")
plot(p_age2[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25, main= "Age 12-25")

par(mfrow=c(1,2))
plot(p_age3[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25, main= "Age 26-49")
plot(p_age4[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25, main= "Age 50-69")

par(mfrow=c(1,2))
plot(p_age5[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25, main= "Age 70-79")
plot(p_age6[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25, main= "Age 80+")

Age Categories

Length of stay estimates

LOS_mat_age %>%
  kbl(booktabs = T,linesep = c("", "", "\\hline", "", "")) %>%
  kable_styling(latex_options = c("hold_position","scale_down"))

\newpage

Calendar Time

Plots of Predicted Proportions

# calculate transition probabilities
p_wave1 <- prob.tran.mat(my.data_ext,my.data_ext$wave == 1, tra)

# ELOS and probability of long stayer
LOS_mat_wave1 <- ELOS_mat(p_wave1[[1]], end_limit, my.data_ext$wave == 1,
                          num_boot)
# calculate transition probabilities
p_wave2 <- prob.tran.mat(my.data_ext,my.data_ext$wave == 2, tra)

# ELOS and probability of long stayer
LOS_mat_wave2 <- ELOS_mat(p_wave2[[1]], end_limit, my.data_ext$wave == 2,
                          num_boot)
# calculate transition probabilities
p_wave3 <- prob.tran.mat(my.data_ext,my.data_ext$wave == 3, tra)

# ELOS and probability of long stayer
LOS_mat_wave3 <- ELOS_mat(p_wave3[[1]], end_limit, my.data_ext$wave == 3,
                          num_boot)
LOS_mat_wave <- cbind(LOS_mat_wave1[[1]], LOS_mat_wave2[[1]], LOS_mat_wave3[[1]])

colnames(LOS_mat_wave) <- c("Wave 1", "Wave 2", "Wave 3")

# LOS_mat_wave %>%
#   kbl() %>%
#   kable_styling()
par(mfrow=c(1,3))

plot(p_wave1[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25, main= "Wave 1 (01.01.20 - 31.07.20)")


plot(p_wave2[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25, main= "Wave 2 (01.08.20 - 30.06.21)")


plot(p_wave3[[1]], from = 1,ord = c(4,2,1,3),type= "filled", cols = col_vec,
     lwd= 2, xlab = "Days Since Admission", ylab = "Predicted Probabilities", xlim = c(0,end_limit_plot), cex.lab = 1.25, main= "Wave 3 (01.07.21 - )")

Calendar time

Length of stay estimates

LOS_mat_wave %>%
  kbl(booktabs = T, linesep = c("", "", "\\hline", "", "")) %>%
  kable_styling(latex_options = c("hold_position"))

\newpage

Competing Risks Model

Cause Specific Cox Regression

Death

Cox_Death <- coxph(Surv( exit, event == 1) ~ sex + age_cat + wave ,  data= my.data_CR, method = "breslow") 
prty_Cox_Death <- prettify(summary(Cox_Death))

prty_Cox_Death <- prty_Cox_Death %>% mutate_if(is.numeric, ~round(., 3))

prty_Cox_Death <- prty_Cox_Death[,c(1, 3:5, 8:9)]

kbl(prty_Cox_Death, booktabs = T, linesep = c("\\hline", "", "", "","", "\\hline", "",""))

Competing Risks Model

Cause Specific Cox Regression

Discharge

Cox_Dis <- coxph(Surv( exit, event == 2) ~ sex + age_cat + wave ,  data= my.data_CR, method = "breslow") 

prty_Cox_Dis <- prettify(summary(Cox_Dis))

prty_Cox_Dis <- prty_Cox_Dis %>% mutate_if(is.numeric, ~round(., 3))

prty_Cox_Dis <- prty_Cox_Dis[,c(1, 3:5, 8:9)]

kbl(prty_Cox_Dis, booktabs = T,linesep = c("\\hline", "", "", "","", "\\hline", "",""))

\newpage

Competing Risks Model

Fine and Gray Model

Death

# event has to be factor for finegray function

my.data_CR$event_FG <- as.factor(my.data_CR$event)


my.data_FG_dea <- finegray(Surv(exit, event_FG) ~ . ,  data= my.data_CR, 
                           etype = 1)




Cox_Death_FG <- coxph(Surv(fgstart, fgstop, fgstatus) ~ sex + age_cat + wave, ties = 'breslow', data = my.data_FG_dea, weight= fgwt)

prty_Cox_Death_FG <- prettify(summary(Cox_Death_FG))

prty_Cox_Death_FG <- prty_Cox_Death_FG %>% mutate_if(is.numeric, ~round(., 3))

prty_Cox_Death_FG <- prty_Cox_Death_FG[,c(1, 3:5, 9:10)]

kbl(prty_Cox_Death_FG, booktabs = T, linesep = c("\\hline", "", "", "","", "\\hline", "",""))

Competing Risks Model

Fine and Gray Model

Discharge

my.data_FG_dis <- finegray(Surv(exit, event_FG ) ~ id + age_cat + sex + wave ,  data= my.data_CR, etype = 2)


Cox_Dis_FG <- coxph(Surv(fgstart, fgstop, fgstatus) ~ sex + age_cat + wave, ties = 'breslow', data = my.data_FG_dis, weight= fgwt )

prty_Cox_Dis_FG <- prettify(summary(Cox_Dis_FG))

prty_Cox_Dis_FG <- prty_Cox_Dis_FG %>% mutate_if(is.numeric, ~round(., 3))

prty_Cox_Dis_FG <- prty_Cox_Dis_FG[,c(1, 3:5, 9:10)]

kbl(prty_Cox_Dis_FG, booktabs = T,linesep = c("\\hline", "", "", "","", "\\hline", "",""))


dyhazard/Phase2.1MultStateModPackage documentation built on Dec. 20, 2021, 2:17 a.m.