Handling correlated fixed effects

Question

I have been performing Cox mixed effects regressions using coxme in R. These models have survival as the response, up to 6 interacting fixed effects and 2 random effects in the full model. I am concerned that multicollinearity may be occuring. For example I use the yields of several crops; if these are influenced by climatic variation then it is possible for them to covary with each other. Here's an example of the model, imagine that crop1 and crop3 are correlated, how important is this?

coxme(Survival ~ S * Crop1 * Crop2 * Crop3 * TreatmentA * TreatmentB * TreatmentC + (1|site) + (1|year), data = df1)

I have tried to explore this using VIF (car package) but it could be problematic with a Cox regression because there is no intercept.

• Do I need to worry about correlated fixed effects? How much correlation is too much correlation?

• What is the correct way to deal with correlated fixed effects?

• Should I discard the approach using multiple crops, and switch to an approach where I would test it with three models of (where crop is either crop1, crop2, or crop3):

coxme(Survival ~ S * Crop * TreatmentA * TreatmentB * TreatmentC + (1|site) + (1|year), data = df1)

Crops 1-3 have 45, 45, and 61 "levels" each, treatments are all binary, $\rho_{1,2} = −0.5; \rho_{1,2} = 0.4; \rho_{1,3} = -0.4$. There are 2800 individuals, with a minimum of 10 individuals per level of each crop.

---

A clear (keep it simple for me ;D) & well referenced answer will be rewarded with a bounty

Answer 1

First, note that multicollinearity is primarily a function of the predictors in a regression. To gauge multicollinearity per se you could calculate VIFs based on regression of each predictor against all of the others, as explained on the Wikipedia page. For models fit by maximum (partial) likelihood like coxme(), that measure of multicollinearity does not necessarily have the strict relationship with coefficient variance inflation seen in standard linear regression. The basic idea that high predictor multicollinearity can mean high variance of coefficient estimates is nevertheless the same.*

Second, your multicollinearity calculation would necessarily include not only the 7 individual fixed predictor variables but also all of the 2-way, 3-way,..., 7-way interactions among them. As your data don't seem to be orthogonal, that will tend to lead to high multicollinearity.

Third, multicollinearity isn't necessarily a problem, particularly when it arises from interaction terms. See this question and the link in its answer for more details. So I don't see any advantage to treating each of the Crops separately, as that will lose information that their combinations of values might provide.

Yes, with multicollinearity your regression coefficients will have higher errors than they would in their absence. Nevertheless, multicollinearity seems to be inherent in your data. So accept it and deal with it. One way to gauge this problem with coefficients would be to repeat the model building process on multiple bootstrap samples of your data and examine variability of the coefficients among those models.

Ridge regression is a well accepted way to deal with multicollinearity, as it nicely combines information among correlated predictors while maintaining all in the model. The coxme help page (page 4) shows an example of how to incorporate "Shrinkage effects (equivalent to ridge regression)." I have no experience with this, however, and it looks like you would have to specify and scale all of the interaction terms yourself. Other than that, I'm not aware of any R functions, at least, that do both mixed effects and ridge regression for Cox models.

Finally, think very carefully about the way you have set up the model, based on your knowledge of the underlying subject matter. Do you really need all of those interaction terms? It's easy to include them all in your formula, but that might not be what you really want. Even without multicollinearity the interaction terms take away degrees of freedom from other predictors and might make it harder for you to document the significance of important predictors.

---

*The vif() method in the R car package is based directly on the coefficient variance-covariance matrix produced during the maximum-likelihood fitting of the model, as shown here.