My data have non-proportional hazards with clear separation. Should I handle it via stratified Cox regression or using separate Log-rank test within subsets?
I will use R only to illustrate.
I want to split the follow up into two periods, T<100 days and T>= 100 days. For the Cox analysis I can use the survSplit(), but survdiff doesn't work with it.
lung$ph.karno_cat <- ifelse(lung$ph.karno < 100, "A", "B")
lung.split <- survSplit(Surv(time, status) ~ ., data= lung, cut=c(180, 350), episode= "tgroup", id="id")
> survdiff(Surv(tstart, time, status) ~ ph.karno_cat:strata(tgroup), data=lung.split)
Error in survdiff(Surv(tstart, time, status) ~ ph.karno_cat:strata(tgroup), :
Right censored data only
Having just one categorical covariate, can I just use the p-values from the coxph() per each stratum? They test whether the coefficient is non-zero and should correspond to the log-rank, isn't it?
> coef(summary(coxph(Surv(tstart, time, status)
~ ph.karno_cat:strata(tgroup), data=lung.split) ))
coef exp(coef) se(coef) z Pr(>|z|)
ph.karno_catA:strata(tgroup)tgroup=1 0.82367716 2.278864 0.5170404 1.59306146 0.1111464
ph.karno_catB:strata(tgroup)tgroup=1 NA NA 0.0000000 NA NA
ph.karno_catA:strata(tgroup)tgroup=2 0.65219820 1.919756 0.4709701 1.38479747 0.1661144
ph.karno_catB:strata(tgroup)tgroup=2 NA NA 0.0000000 NA NA
ph.karno_catA:strata(tgroup)tgroup=3 -0.03294059 0.967596 0.3683481 -0.08942787 0.9287419
ph.karno_catB:strata(tgroup)tgroup=3 NA NA 0.0000000 NA NA
Or should I use the survdiff on filtered data?
> survdiff(Surv(time, status) ~ ph.karno_cat, data=lung, subset = time < 180)
Call:
survdiff(formula = Surv(time, status) ~ ph.karno_cat, data = lung,
subset = time < 180)
n=67, 1 observation deleted due to missingness.
N Observed Expected (O-E)^2/E (O-E)^2/V
ph.karno_cat=A 63 57 59.17 0.0795 2.72
ph.karno_cat=B 4 4 1.83 2.5704 2.72
Chisq= 2.7 on 1 degrees of freedom, p= 0.1
> survdiff(Surv(time, status) ~ ph.karno_cat, data=lung, subset = time >= 180 & time < 350)
Call:
survdiff(formula = Surv(time, status) ~ ph.karno_cat, data = lung,
subset = time >= 180 & time < 350)
N Observed Expected (O-E)^2/E (O-E)^2/V
ph.karno_cat=A 72 46 40.3 0.82 4.07
ph.karno_cat=B 13 5 10.7 3.07 4.07
Chisq= 4.1 on 1 degrees of freedom, p= 0.04
and
> survdiff(Surv(time, status) ~ ph.karno_cat, data=lung, subset = time >= 350)
Call:
survdiff(formula = Surv(time, status) ~ ph.karno_cat, data = lung,
subset = time >= 350)
N Observed Expected (O-E)^2/E (O-E)^2/V
ph.karno_cat=A 63 43 43.08 0.000143 0.000839
ph.karno_cat=B 12 9 8.92 0.000689 0.000839
Chisq= 0 on 1 degrees of freedom, p= 1
The p-values are different, but quite close to the p-values for the Cox model. Which option would you prefer and why?
EDIT: I used the statistical package only to illustrate the problem. I am not asking for anyone to write code for me and I don't ask for datasets or debugging.