Compare the fit of two mixed effects models estimated with robust standard errors

[JaimeAndres Castaneda]

Dear all,

I'm using Stata 14.

I'm trying to compare the fit of two mixed effects models estimated with robust standard errors.

I know that lrtest doesn't work when the models are estimated with robust standard errors. Below I provide a code that exemplifies this:

mixed performance i.profit i.feedback i.behavior || id: R.profit, vce(r)
estimates store A
mixed performance i.profit i.feedback || id: R.profit, vce(r)
estimates store B
lrtest A B

(Where performance captures some performance metric from my sample of participants, profit is a factor variable that captures three profit levels (e.g., low, medium, high), feedback is a factor variable that captures three feedback types (e.g., simple, normal, complex), behavior is a factor variable that captures four behavioral patterns from my sample of participants (e.g., pattern A, B, C, D) and id identifies my participants)

Stata would show the following error:

LR test likely invalid for models with robust vce
r(498);

I know I can make lrtest run by using the option force:

lrtest A B, force

But I understand that this approach could be questionable.

If I estimate those two models without robust standard errors and compare their fit with lrtest, the command works fine.

Is there a command analogous to lrtest or a set of commands that allow to compare the fit of two mixed effects models estimated with robust standard errors?

Thanks!

[Clyde Schechter]

Well, the likelihood ratio is not, in any case a comparison of goodness of fit of two models. It is a joint significance test of the coefficients of the variables that differ between the models. In your case, it is a significance test for the variable behavior. So if you want a significance test for the variable behavior, and you can't have a likelihood ratio test because you need robust standard errors, use the one that's already in the regression output table from the regression that includes i.behavior. If you want a name for it, it's called a Wald test. When non-robust errors are used, the likelihod ratio and Wald tests usually give very similar, though typically not identical, results.

As far as actual goodness of fit is concerned, the model with more variables always fits better than the model with fewer. The issue is not whether adding i.behavior improves fit--the question should be whether the improvement of fit outweighs the potential overfitting that comes with adding more variables. Many people use the Akaike or Bayes information criteria to make this assessment (AIC or BIC). You can get those statistics by running -estat ic- after each of the models and seeing the differences between the models. The model with the lower AIC or BIC is considered to be better overall. (Unfortunately, sometimes the AIC and BIC give opposing verdicts on this question, but usually they agree.) People often also adopt some minimum threshold of difference between models on these statistics, and differences smaller than that are informally considered to mean the models are roughly equivalent overall.

[JaimeAndres Castaneda]

Hello Clyde,

Thank you for the clarification about what the likelihood ratio test does.

And thank you also for the information about the AIC and BIC and the command -estat ic-. I heard about the AIC and BIC before and have used the AIC in other research. Before starting this post, I was thinking to use the AIC, but I was wondering whether it was right to compute the AIC using the log pseudolikelihood that Stata shows in the regression output table. An intuition based on my basic econometrics background was telling me that maybe the distinction between a likelihood and a pseudolikelihood could make a difference when computing an AIC. I did an initial search and couldn't find much to guide a deeper and more focused search on that (then I decided to come to this forum).

It turns out that the log likelihood reported by -estat ic- is the same as the log pseudolikelihood reported in the regression output table. I don't know if this an artifact of my data and/or the use of mixed effects estimation or it's something normal (that I shouldn't be worried about).

With my data, both the AIC and BIC lead to the same conclusion.

Thanks again!

[Clyde Schechter]

StataCorp is pretty good about programming commands so that they do not attempt to calculate a statistic (either throwing an error message or returning a missing value for the statistic) when it is not legitimate to do so. So I think you are OK here.

I am going out a bit on a limb here because I am not certain of what follows, but I believe that when Stata estimates -mixed- models with robust standard errors, it first does a regular maximum likelihood estimation, and then subsequently adjusts the standard errors using the sandwich (Huber-White) estimator. Why they choose to call the likelihood a pseudo-likelihood in the output is not clear to me, but it is, indeed, the actual likelihood for the model, and statistics like AIC and BIC derived from it are still OK. The AIC and BIC are not used for hypothesis testing and the validity of the standard errors derived from the likelihood function, or lack thereof, are not an issue. (This contrasts with the situation for a likelihood ratio test: by using the robust standard errors, you are stating that you do not believe that the usual standard errors derived from the information matrix, which is a second derivative of the likelihood function, are not valid, and so tests that correspond to that calculation are not valid. Consequently, Stata will not give you a likelihood-ratio test comparing two models in this situation.)

[JaimeAndres Castaneda]

Hello Clyde,

Thanks again for your reply and the clarifications!