Cross Validated is a question and answer site for people interested in statistics, machine learning, data analysis, data mining, and data visualization. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

Reading the Linear Mixed Model (LMM) literature I am aware that fitting a model using REML provides better estimates of variance parameters than fitting via ML. However, we should not compare nested models fitted with REML that have different fixed effects.

Recently, I have been fitting some models using GLS via the gls() function in the nlme package for R. The default fitting method for that function is REML. Do the same principals of REML vs ML for LMM also apply to GLS?

Specifically, I am fitting a model with and without a linear trend with a correlation structure in the residuals:

m1 <- gls(Response ~ Time, data = foo, correlation = corAR1(form ~ Time))
m0 <- gls(Response ~ 1, data = foo, correlation = corAR1(form ~ Time))

In the above, I should fit the models using ML as they have different fixed effects. Is this correct?

Secondly, consider two GLS models that differ in the correlation structure:

m1 <- gls(Response ~ Time, data = foo, correlation = corARMA(form ~ Time, p = 1))
m2 <- gls(Response ~ Time, data = foo, correlation = corARMA(form ~ Time, p = 2))

What fitting method should ideally be used here? REML or ML? Here my intuition would say fit via REML as we are estimating (co)variance parameters. Is my intuition correct or have I got this all mixed up?

share|improve this question
up vote 5 down vote accepted

Your intuition is correct, the same principles apply. I looked in Pinheiro/Bates section 5.4, where gls is introduced, but it doesn't say so explicitly, so you'll just have to trust me, I guess. :)

In Chapter 2 they go through the theory of REML and ML and you'll notice that none of the theory depends on there being any random effects, and that actually, you could write any random effect model using just correlation structure instead and fit with gls, though for complex random effects it would be quite complex. The simplest example is that a random intercept model is equivalent to a compound symmetry model.

share|improve this answer
Many thanks for the answer and the pointer to the theory in PB. I had looked at 5.4 but like you didn't see anything there. Will take a look at the earlier chapter. – Gavin Simpson Jun 19 '11 at 9:21

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

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