gls() vs. lme() in the nlme package

Question

In the nlme package there are two functions for fitting linear models (lme and gls).

1. What are the differences between them in terms of the types of models that can be fit, and the fitting process?
2. What is the design rational for having two functions to fit linear mixed models where most other systems (e.g. SAS SPSS) only have one?

Update: Added bounty. Interested to know differences in the fitting process, and the rational.

Answer 1

Interesting question.

In principle the only difference is that gls can't fit models with random effects, whereas lme can. So the commands

fm1 <- gls(follicles ~ sin(2*pi*Time) + cos(2*pi*Time),Ovary,correlation=corAR1(form=~1|Mare))

and

lm1=lme(follicles~sin(2*pi*Time)+cos(2*pi*Time),Ovary,correlation=corAR1(form=~1|Mare))

ought to give the same result but they don't. The fitted parameters differ slightly.

Answer 2

From Pinheiro & Bates 2000, Section 5.4, p250:

The gls function is used to fit the extended linear model, using either maximum likelihood, or restricted maximum likelihood. It can be veiwed as an lme function without the argument random.

For further details, it would be instructive to compare the lme analysis of the orthodont dataset (starting on p147 of the same book) with the gls analysis (starting on p250). To begin, compare

---

orth.lme <- lme(distance ~ Sex * I(age-11), data=Orthodont)
summary(orth.lme)

Linear mixed-effects model fit by REML
 Data: Orthodont 
       AIC     BIC    logLik
  458.9891 498.655 -214.4945

Random effects:
 Formula: ~Sex * I(age - 11) | Subject
 Structure: General positive-definite
                      StdDev    Corr                
(Intercept)           1.7178454 (Intr) SexFml I(-11)
SexFemale             1.6956351 -0.307              
I(age - 11)           0.2937695 -0.009 -0.146       
SexFemale:I(age - 11) 0.3160597  0.168  0.290 -0.964
Residual              1.2551778                     

Fixed effects: distance ~ Sex * I(age - 11) 
                          Value Std.Error DF  t-value p-value
(Intercept)           24.968750 0.4572240 79 54.60945  0.0000
SexFemale             -2.321023 0.7823126 25 -2.96687  0.0065
I(age - 11)            0.784375 0.1015733 79  7.72226  0.0000
SexFemale:I(age - 11) -0.304830 0.1346293 79 -2.26421  0.0263
 Correlation: 
                      (Intr) SexFml I(-11)
SexFemale             -0.584              
I(age - 11)           -0.006  0.004       
SexFemale:I(age - 11)  0.005  0.144 -0.754

Standardized Within-Group Residuals:
        Min          Q1         Med          Q3         Max 
-2.96534486 -0.38609670  0.03647795  0.43142668  3.99155835 

Number of Observations: 108
Number of Groups: 27

---

orth.gls <- gls(distance ~ Sex * I(age-11), data=Orthodont)
summary(orth.gls)

Generalized least squares fit by REML
  Model: distance ~ Sex * I(age - 11) 
  Data: Orthodont 
       AIC      BIC    logLik
  493.5591 506.7811 -241.7796

Coefficients:
                          Value Std.Error  t-value p-value
(Intercept)           24.968750 0.2821186 88.50444  0.0000
SexFemale             -2.321023 0.4419949 -5.25124  0.0000
I(age - 11)            0.784375 0.1261673  6.21694  0.0000
SexFemale:I(age - 11) -0.304830 0.1976661 -1.54214  0.1261

 Correlation: 
                      (Intr) SexFml I(-11)
SexFemale             -0.638              
I(age - 11)            0.000  0.000       
SexFemale:I(age - 11)  0.000  0.000 -0.638

Standardized residuals:
        Min          Q1         Med          Q3         Max 
-2.48814895 -0.58569115 -0.07451734  0.58924709  2.32476465 

Residual standard error: 2.256949 
Degrees of freedom: 108 total; 104 residual

---

Notice that the estimates of the fixed effects are the same (to 6 decimal places), but the standard errors are different, as is the correlation matrix.