# How to interpret interaction results of subgroup analysis in cox model in R

Context

I am using survival::coxph() to fit the cox model. I want to get HR for exposure factors under different sex.

I am using the dataset lung from the survival package to simplify the problem I am having.

Say that I want to get HR of ph.karno (Continuous variable) in different sex (Categorical variable). I can do this in two ways.

Method 1: Use the interaction term of sex and ph.karno to get the HR of ph.karno across sex.

Method 2: Filter the dataset containing only males (sex == 1) or females (sex == 1) separately. Then fit the cox model using these two datasets separately.

Both of the above methods can get HR of ph.karno in different sex. According to the results, the results obtained by the two methods are not consistent (see Reproducible code).

Question

How can I explain the results of the two methods?

Which method obtains the correct HR of ph.karno in different sex.

Reproducible code

library(survival)

# method 1: get ph.karno effect in different sex by adding an interact term

fit <- coxph(Surv(time, status) ~ factor(sex):ph.karno, data = lung)

interact_result = fit\$coefficients

# method 2: get ph.karno effect in different sex by subsetting differnt data

lung = lung
lung_1 = lung[lung$$sex == 1,] lung_2 = lung[lung$$sex == 2,]

fit1 = coxph(Surv(time, status) ~ ph.karno, data = lung_1)
fit2 = coxph(Surv(time, status) ~ ph.karno, data = lung_2)

subgroup_result = c(sex_1 = fit1$$coefficients, sex_2 = fit2$$coefficients)

# comparing ph.karno effect in different sexe using the two methods described above

rbind(interact_result, subgroup_result)

#                 factor(sex)1:ph.karno factor(sex)2:ph.karno
# interact_result          -0.013675278           -0.02030161
# subgroup_result          -0.009673906           -0.03744254
$$$$


Either method could be correct, depending on your assumption about the baseline survival function.

With method 1 you are implicitly assuming a single shared baseline hazard for both sexes.

With method 2 you are allowing different baseline hazards for each sex. You get the same results from the following stratified model.

fit3 <- coxph(Surv(time, status) ~ ph.karno:strata(sex), data = lung)


Thus the best method depends on whether you think that the proportional hazard (PH) assumption holds for sex.

In general, the best way to proceed is to start with a combined model that includes the full interaction term sex*ph.karno to get individual coefficients for each of the predictors along with the interaction. In more complicated models you can get into trouble trying to interpret coefficients when some individual coefficients are omitted.

Then evaluate whether proportional hazards hold. If it doesn't hold for sex then you might consider proceeding with a model stratified by sex`.