Checking PH assumption using scaled Schoenfeld residual tests

My cox model violates the PH assumption (according to cox.zph test), and I'm trying to interpret the scaled Schoenfeld residual plot, but I've never seen a plot with points lying in such a flat manner. I'm wondering what the beta(t) means and if it's possible that almost all the points, across time, have almost the same values...

I'd appreciate any help!

• What is plotted here? Is Beta(t) the df-beta for a Cox time-to-event model with one covariate (sex)? What code gave you this graph? – AdamO Dec 20 '19 at 16:40
• model<-coxph(Surv(time, status) ~ sex, data) then ggcoxzph(cox.zph(model)) was the code! – llbia Dec 20 '19 at 16:43

The Schoenfeld residuals take the difference between the scaled covariate value(s) for the i-th observed failure and what is expected by the model.

So for your model, you have a single binary covariate sex which takes values 0 or 1. Supposing there is 0 association and supposing sex is balanced in the sample, the "expected" covariate would be 0.5, that is, any given event is equally likely to have belonged to either a male or a female at any time. This amount can be offset by the number of males or females in the overall sample or whether they are right censored and entered into the sample in a staggered fashion. The lack of variability and association in the covariate is what's essentially being reflected in the plot.

For the smoothing spline that trends somewhat upward, we are clued into a possible accelerated failure time effect ($$\beta(t) \ne 0$$). A tried-and-true sensitivity analysis would be to inspect a Kaplan-Meier plot for binary sex and inspect whether the survival curves cross in the sample.

We can simulate this using the simple result that an exponential random variable has a constant hazard whereas a Weibull model has a linear hazard that is not flat when the shape is not 1.

Take:

set.seed(123)
group <- rep(0:1, each=1000)
t1 <- rexp(1000)
t2 <- rweibull(1000, 2, 1/gamma(1+1/2))
t <- c(t1, t2)
par(mfrow=c(1,2))
plot(survfit(Surv(t) ~ group), col=1:2, xlab='Survival time', ylab='Surviving proportion')
fit <- coxph(Surv(t) ~ group)
summary(fit)
proph <- cox.zph(fit)
proph
plot(proph)


Gives:

> summary(fit)
Call:
coxph(formula = Surv(t) ~ group)

n= 2000, number of events= 2000

coef exp(coef) se(coef)     z Pr(>|z|)
group 0.10060   1.10584  0.04675 2.152   0.0314 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

exp(coef) exp(-coef) lower .95 upper .95
group     1.106     0.9043     1.009     1.212

Concordance= 0.465  (se = 0.007 )
Likelihood ratio test= 4.64  on 1 df,   p=0.03
Wald test            = 4.63  on 1 df,   p=0.03
Score (logrank) test = 4.63  on 1 df,   p=0.03


i.e. a marginally significant hazard ratio because the expected failure time is close to exactly equal. However, the exponential model shows more earlier failures and this tends to power Cox models.

The failure time averaged hazard ratio comes out to 1.1. It should be noted this is the correct interpretation, there's nothing wrong with the estimate it's telling us exactly what's in the data. If you use the robust=TRUE option, the 95% CI is correct and there's no reason to discredit the Cox model except that possible answers to the general question: which group survives better, may have better, more nuanced answers. Thus is the result of any sensitivity analysis, rather than "invalidating models".

As we see, the proportional hazards test is highly significant.

> proph
rho chisq        p
group 0.325   230 6.19e-52


And the plot is similar to yours. My example shows more extreme trends because I have exaggerated the accelerated failure time effect in the Weibull model. This is mostly to illustrate the reason why we see a temporal trend in the smoothing, and you also see the "parallel groups" for covariate values both high and low.

• Thank you so much for the answer!!!!I'm using retrospective data in which participants were followed up about five years after the initial survey, and at that point, some people had experienced the event of interest, and they recorded (retrospectively) the month that event occurred. But I realized now that this means for those who hasn't experience the event ("survivors"), "time" variable would be how many months have passed since the last survey and would not have as much variability (the minimum being like 52 months...). I guess this isn't data I could use Cox!! Thanks ALOT – llbia Dec 20 '19 at 19:37