Since the Cox PH model is a non-linear model, we would expect the effect to be non-collapsible. i.e., the marginal and conditional effects differ.
I did some calculation for a setting where the baseline hazard is constant, and the covariate is normally distributed and treatment is binary. In this case, would the effect actually be collapsible? See the calculation below.
Assume a treatment Z is binary (0, 1), a covariate x has standard normal distribution N(0, 1), and the baseline hazard is exponential with constant hazard $\lambda_0$, the conditional effect is $\beta$.
Suppose the hazard function is $\lambda=\lambda_0 exp(X+\beta Z)$.
The marginal hazard for $Z=1$ would be: $$\int_{-\infty}^\infty 1/{\sqrt(2\pi)}\exp(x+\beta)\exp(-1/2 x^2) dx = \exp(\beta+0.5) $$
The marginal hazard for $Z=0$ would be: $$\int_{-\infty}^\infty 1/{\sqrt(2\pi)}\exp(x)\exp(-1/2 x^2) = \exp(0.5)$$
The ratio of the two would be $\exp(\beta)$, which is the same as the conditional hazard ratio.
Not sure if my understanding is correct.
Best Answer
Maths is hard, so I would always start by checking this sort of thing with simulation
It doesn't seem to be collapsible. So what's wrong with the maths?
Well, here are smooth estimates of the hazard rate for the two $Z$ groups
They aren't constant. The marginal hazard decreases with time (for any mixture of exponentials), so your calculations can't be right.
What I think you've calculated is $E[\exp(X+\beta Z)|Z=1]$ and $E[\exp(X+\beta Z)|Z=0]$. Those aren't the marginal hazard (except at time zero) because they don't account for how the population at risk changes over time.
In fact, the hazard ratio at time zero, estimated from the smoother, is very close to
exp(1)
, but it decreases over time.