Comparing mixed effect models with the same number of degrees of freedom

Question

I have an experiment that I'll try to abstract here. Imagine I toss three white stones in front of you and ask you to make a judgment about their position. I record a variety of properties of the stones and your response. I do this over a number of subjects. I generate two models. One is that the nearest stone to you predicts your response, and the other is that the geometric center of the stones predicts your response. So, using lmer in R I could write.

mNear   <- lmer(resp ~ nearest + (1|subject), REML = FALSE)
mCenter <- lmer(resp ~ center  + (1|subject), REML = FALSE)

UPDATE AND CHANGE - more direct version that incorporates several helpful comments

I could try

anova(mNear, mCenter)

Which is incorrect, of course, because they're not nested and I can't really compare them that way. I was expecting anova.mer to throw an error but it didn't. But the possible nesting that I could try here isn't natural and still leaves me with somewhat less analytical statements. When models are nested naturally (e.g. quadratic on linear) the test is only one way. But in this case what would it mean to have asymmetric findings?

For example, I could make a model three:

mBoth <- lmer(resp ~ center + nearest + (1|subject), REML = FALSE)

Then I can anova.

anova(mCenter, mBoth)
anova(mNearest, mBoth)

This is fair to do and now I find that the center adds to the nearest effect (the second command) but BIC actually goes up when nearest is added to center (correction for the lower parsimony). This confirms what was suspected.

But is finding this sufficient? And is this fair when center and nearest are so highly correlated?

Is there a better way to analytically compare the models when it's not about adding and subtracting explanatory variables (degrees of freedom)?

Answer 1

Following ronaf's suggestion leads to a more recent paper by Vuong for a Likelihood Ratio Test on nonnested models. It's based on the KLIC (Kullback-Leibler Information Criterion) which is similar to the AIC in that it minimizes the KL distance. But it sets up a probabilistic specification for the hypothesis so the use of the LRT leads to a more principled comparison. A more accessible version of the Cox and Vuong tests is presented by Clarke et al; in particular see Figure 3 which presents the algorithm for computing the Vuong LRT test.

• Likelihood Ratio Tests for Model Selection and Non-nested Hypotheses (Vuong, 1999)
• Testing Nonnested Models of International Relations: Reevaluating Realism (Clarke et al, 2000)

It seems there are R implementations of the Vuong test in other models, but not lmer. Still, the outline mentioned above should be sufficient to implement one. I don't think you can obtain the likelihood evaluated at each data point from lmer as required for the computation. In a note on sig-ME, Douglas Bates has some pointers that might be helpful (in particular, the vignette he mentions).

---

Older

Another option is to consider the fitted values from the models in a test of prediction accuracy. The Williams-Kloot statistic may be appropriate here. The basic approach is to regress the actual values against a linear combination of the fitted values from the two models and test the slope:

• A Test for Discriminating Between Models (Atikinson, 1969)
• Growth and the Welfare State in the EU: A Causality Analysis (Herce et al, 2001)

The first paper describes the test (and others), while the second has an application of it in an econometric panel model.

---

When using lmer and comparing AICs, the function's default is to use the REML method (Restricted Maximum Likelihood). This is fine for obtaining less biased estimates, but when comparing models, you should re-fit with REML=FALSE which uses the Maximum Likelihood method for fitting. The Pinheiro/Bates book mentions some condition under which it's OK to compare AIC/Likelihood with either REML or ML, and these may very well apply in your case. However, the general recommendation is to simply re-fit. For example, see Douglas Bates' post here:

• How can I extract the AIC score from a mixed model object produced using lmer?

Answer 2

Still, you can compute confidence intervals for your fixed effects, and report AIC or BIC (see e.g. Cnann et al., Stat Med 1997 16: 2349).

Now, you may be interested in taking a look at Assessing model mimicry using the parametric bootstrap, from Wagenmakers et al. which seems to more closely resemble your initial question about assessing the quality of two competing models.

Otherwise, the two papers about measures of explained variance in LMM that come to my mind are:

• Lloyd J. Edwards, Keith E. Muller, Russell D. Wolfinger, Bahjat F. Qaqish and Oliver Schabenberger (2008). An R2 statistic for fixed effects in the linear mixed model, Statistics in Medicine, 27(29), 6137–6157.
• Ronghui Xu (2003). Measuring explained variation in linear mixed effects models, Statistics in Medicine, 22(22), 3527–3541.

But maybe there are better options.

Answer 3

there is a paper by d.r.cox that discusses testing separate [unnested] models. it considers a few examples, which do not rise to the complexity of mixed models. [as my facility with R code is limited, i'm not quite sure what your models are.]

altho cox's paper may not solve your problem directly, it may be helpful in two possible ways.

1. you can search google scholar for citations to his paper, to see if subsequent such results come closer to what you want.

2. if you are of an analytical bent, you could try applying cox's method to your problem. [perhaps not for the faint-hearted.]

btw - cox does mention in passing the idea srikant broached of combining the two models into a larger one. he doesn't pursue how one would then decide which model is better, but he remarks that even if neither model is very good, the combined model might give an adequate fit to the data. [it's not clear in your situation that a combined model would make sense.]

Answer 4

I do not know R well enough to parse your code but here is one idea:

Estimate a model where you have both center and near as covariates (call this mBoth). Then mCenter and mNear are nested in mBoth and you could use mBoth as a benchmark to compare the relative performance of mCenter and mNear.