I'm working with acute myeloid leukemia (AML) data. I want to study if there are differences in Overall Survival between the patients who received an allogeneic transplant and those who did not. I have performed the Proportional Hazards Assumption of a Cox Regression test and plotted the residuals, and we can observe that the assumption is not satisfied. Here is the code:
library(survival)
library(survminer)
OS_AlloHCT <- coxph(Surv(SGmesos, Estadoactual) ~ AlloTPH, data = AML)
test.ph <- cox.zph(OS_AlloHCT)
test.ph
# chisq df p
# AlloTPH 93.3 1 <2e-16
# GLOBAL 93.3 1 <2e-16
ggcoxzph(test.ph)
The plot of the residuals leads us to conclude that the Proportional Hazards Assumption is not met. For this reason, and due to clinical considerations, I have decided to treat this variable as a time-dependent covariate. To do this, I have created new variables in the database: txtime, start, and stop.
FechaAloTPH indicates the transplant date, Fechadiagnóstico is the diagnosis date, and Fechaúltimavisita is the last follow-up visit. txtime represents the difference, in days, between the diagnosis date (Fechadiagnóstico) and the transplant date (FechaAloTPH). It needs to be a numeric variable:
AML$txtime <- difftime(AML$FechaAloTPH, AML$Fechadiagnóstico,"days")
AML$txtime <- as.numeric(AML$txtime)
AloTPH is a binary variable that takes the value 1 if the patient has received a bone marrow transplant and 0 if not. If the patient has received a bone marrow transplant and the days of overall survival are greater than the days until the transplant (which is logical), the start date will be set to txtime (the number of days until the transplant date). Otherwise, the start date will be set to 0.
If the patient has received a bone marrow transplant and the days of overall survival are greater than the days until the transplant (again, which is logical), the stop date will be set to SGdies (the number of days until the last visit). Otherwise, the stop date will be set to the value of txtime.
AML$start <- ifelse(AML$AloTPH ==1 & (AML$SGdies > AML$txtime),AML$txtime,0)
AML$stop <- ifelse(AML$AloTPH ==1 & (AML$SGdies > AML$txtime),AML$SGdies, AML$txtime)
AML$stop <- ifelse(is.na(AML$stop),AML$SGdies, AML$stop)
Here, we can see the variables required for the analysis:
All the variables have been previously described, and Estadoactual takes the value 1 if the patient is dead and 0 if they are alive.
I have created the Simon-Makuch plot:
library(survival)
library(survminer)
library(ggplot2)
model <- survfit(Surv(start,stop,Estadoactual)~AloTPH, data=AML)
ggsurvplot(model)
The Cox regression model, with the AloTPH variable treated as a time-dependent covariate, is as follows:
AML$Surv <- Surv(AML$start,AML$stop,AML$Estadoactual)
AML$AloTPH<- relevel(as.factor(AML$AloTPH), ref = "1")
cox_OS_AloTPH<- coxph(Surv~AloTPH,data=AML)
summary(cox_OS_AloTPH)
# coef exp(coef) se(coef) z Pr(>|z|)
# AloTPH0 -0.1831 0.8327 0.1082 -1.691 0.0908 .
#---
#Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# exp(coef) exp(-coef) lower .95 upper .95
# AloTPH0 0.8327 1.201 0.6736 1.029
# Concordance= 0.518 (se = 0.009 )
# Likelihood ratio test= 2.82 on 1 df, p=0.09
# Wald test = 2.86 on 1 df, p=0.09
# Score (logrank) test = 2.87 on 1 df, p=0.09
I have discovered an alternative method for creating the Simon-Makuch plot and conducting the Mantel-Byar test:
#install.packages("aod")
library(aod)
#install.packages("RcmdrPlugin.EZR")
library(RcmdrPlugin.EZR)
library(Rcmdr)
TempTD <- AML
Mantel.Byar(Group = "AloTPH", Event = TempTD$Estadoactual,
StartTime = TempTD$start, StopTime = TempTD$stop,
method = c("SAS", "Tominaga"), plot=0, landmark=0)
What do you think of this analysis? Is it correct?