I am using mixed effects Cox models in the R package coxme, with the model

SurvObj ~ Sex*NE + (1|Year)

where Sex is a categorical fixed effect with two levels each (M and F), NE is a continuous predictor of particular interest, and Year the sampling year used as a random effect.

Here's dummy data/script with longer lifespan in Sex = M, and a positive effect of NE on lifespan:

dumDF = data.frame(
    "Sex" = rep(c("M","F"), each = 1000), 
    "Year" = sample(c("A","B"), replace = T), 
    "NE" = rnorm(2000,0,5), 
    "Status" = sample(c(1,2,2,2), 2000, replace = T), 
    "Life" = rnorm(2000, 50, 5))

dumDF$Life = ifelse(dumDF$Sex=="M", dumDF$Life -10 + dumDF$NE, dumDF$Life+dumDF$NE)

dumDF$SurvObj = with(dumDF, Surv(Life, Status == 2))
coxdum <- coxme(SurvObj ~ Sex*NE + (1|Year), data = dumDF)


I would like to get some guidance on how to interpret the result/output:

> coxdum
Cox mixed-effects model fit by maximum likelihood
  Data: dumDF
  events, n = 1480, 2000
  Iterations= 1 7 
                    NULL Integrated    Fitted
Log-likelihood -9793.495  -9068.067 -9068.067

                    Chisq df p     AIC     BIC
Integrated loglik 1450.86  4 0 1442.86 1421.66
 Penalized loglik 1450.86  3 0 1444.86 1428.96

Model:  SurvObj ~ Sex * NE + (1 | Year) 
Fixed coefficients
               coef exp(coef)    se(coef)      z       p
SexM     1.92364026 6.8458338 0.063091867  30.49 0.0e+00
NE      -0.17453933 0.8398438 0.009085403 -19.21 0.0e+00
SexM:NE -0.05491849 0.9465623 0.011945520  -4.60 4.3e-06

Random effects
 Group Variable  Std Dev Variance
 Year  Intercept 2e-02   4e-04   

What I'd like to do is predict the coefficient for the range of NE, specific to either sex. I have done similarly so previously using a different sort of model, but the problem I have here is that the coxme does not give me an intercept for the group where Sex = F and it is necessary for all predictions. How can I predict the coefficent for either group/how can I get the intercept (b1)?

$ coef_{F} = b_1 + b_3 \times x$

$ coef_{M} = (b_1 + b_2) + (b_3 + b_4) \times x$

where $b_1$ is the F specific intercept, $b_2$ is the M specific intercept ("SexM") relative to $b_1$, $b_3$ is the F specific slope ("NE") and $b_4$ is the M specific slope ("SexM:NE") relative to the F specific slope.

Also how does one interpret the number (coef), given that SurvObj is the response, would a negative coefficient for NE mean that there is negative relationship between NE and survival (higher NE = short lifespan/higher rate of mortality) or negative relationship between NE and mortality (higher NE = long lifespan/lower rate of mortality)? Following from my dummy data, where I know NE is positively associated with lifespan, then I believe the coefficient would represent mortality rate, because NE is negative in the example, thus mortality rate decreases with increasing NE.


1 Answer 1


A Cox proportional hazards model does not provide an intercept like the glmm models in the page to which you linked. It is a semi-parametric model, a substantially "different sort of model" from those.

The reference for a Cox model is effectively the empirical hazard function over time, at the covariate values for the baseline class. Conceptually, you could think of that entire empirical hazard function as the "intercept." The coefficient for a predictor variable documents how that predictor affects the overall slope of the log of the hazard function, assuming proportional effects on hazards. Significance tests are properly performed on the coefficients ($\hat\beta$) and their standard errors, but in clinical use at least practical attention is more often on $exp(\hat\beta)$, which represent the hazard ratios relative to the baseline class.

You might find the documentation in the survival package or other introductory explanations to be more helpful for these basic concepts than what is provided in the coxme package.

  • $\begingroup$ thanks, what I am planning on doing is testing for sex differences in the slopes - I can still get these slopes (in this example it's simply comparing b3 with (b3+b4)) and to do this I am using the coefficients, and randomising NE in the data 1000 times, then taking the ratio of randomised slopes which are steeper than the actual slope as a psuedo p value... I'll take a look at the documentation you noted, thanks very much for your answer! $\endgroup$
    – rg255
    Aug 1, 2016 at 17:08
  • $\begingroup$ For sex differences in slopes with female as the baseline level of that factor, all you need to do is test the coefficient reported for the male level; that already is expressed relative to the value for female. In your terminology, b3 is subsumed in the baseline hazard, all you need to test is b4. $\endgroup$
    – EdM
    Aug 1, 2016 at 20:11
  • $\begingroup$ thanks @EdM - you're right, essentially I am testing whether b4 is significant... I have been using the reconstructed slope because in my actual data I have an additional two-level factor fixed effect, thus I have four slopes of interest. Would you agree that if I construct the slopes as I have done in my previous question I can then compare the slopes to one another and test by randomisation? $\endgroup$
    – rg255
    Aug 2, 2016 at 8:08
  • $\begingroup$ If you have two categorical factors each of which has two levels, then there will only be 2 coefficients from them to test in the Cox model, or 3 if you include their interaction. The baseline hazard will be for the reference value of all factors (and at a value of 0 for any continuous predictors in your model). It's not completely clear what you mean by "test by randomisation" but if you mean bootstrap resampling that would be fine. The covariance matrix of the estimated coefficient values also could allow for significance tests. $\endgroup$
    – EdM
    Aug 2, 2016 at 11:31

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.