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).


How can I explain the results of the two methods?

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

Reproducible code


# 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

1 Answer 1


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.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.