[R-sig-ME] [FORGED] Re: Question about gls --- whoops!!!

Berwin A Turlach berwin.turlach at gmail.com
Sat Mar 24 09:06:17 CET 2018

G'day Rolf,

On Sat, 24 Mar 2018 12:39:20 +1300
Rolf Turner <r.turner at auckland.ac.nz> wrote:

> On 24/03/18 12:21, Douglas Bates wrote:
> > Your comment about S4 classes and methods is misdirected. The gls 
> > function is in the nlme package which uses S3 only.  
> Okay Doug; mea culpa.
> But then I must say that you have succeeded in obfuscating your code
> to a degree that I would not have thought possible without the
> assistance of S4 classes and methods! :-)

Even if this statement is spoken in jest, I think I have to disagree
with it.  The code just requires one to diligently follow the way S3
dispatches to figure out what is going on.

Actually, once you realise the rich class of models that the nlme
package supports via gls() and lme() (Pinheiro and Bates, 2000,
Springer-Verlag), will still allowing for extensions (see Galecki, 1994,
Communications in Statistics-Theory and Methods, 23(11): 3105-3119, or
Galecki and Burzykowski, 2012, Springer-Verlag), you will appreciate
that it cannot be a trivial piece of code.

Moreover, the variance parameters in these models have constraints on
them, variance matrices have to be positive semi-definite, the
correlation parameter in an AR(1) model has to be between -1 and 1, and
so forth.  Constrained optimisation is not trivial at the best of time,
non-linear optimisation under constraints is a real pain.  Thus,
reparametrisations that allow for unconstrained estimation are often
employed (see, among others, Bates and Watts, 1988, Wiley & Sons).
Such reparameterisations are widely used by the nlme package.

To figure out how to get from the unconstrained parameterisation which
was used to fit the model, and for which the values are stored in the
object returned by the fitting function, to the parametrisation of
interest, one just has to diligently follow the path that S3 uses when
dispatching to generic functions.  And, typically, it does not take to
long to figure out how to change from one parameterisation to another.
Though, to find out what the parameterisation is used as the
unconstrained parameterisation, one might have to dig quite deep, just
to be faced with a .C() call at the end and the realisation that it is
also necessary to reverse-engineer some C code.

Once you follow the path of how the print() or coef() functions obtain
the values that they return from those stored in the fitted model
object, you will start to appreciate how much thought went into the
design of nlme, how object oriented the design is, how it is
ensured that the various generics play nicely and correctly with
each other, which parts of the code are responsible for what action and
so forth....



More information about the R-sig-mixed-models mailing list