[FORGED] Re: Question about gls --- whoops!!!

G'day Rolf,

On Sat, 24 Mar 2018 11:30:17 +1300

Nope!  ***NOT*** bingo!!!  Louise just pointed out to me that what I
got by applying as.vector(), i.e. 0.3533897, is not the same as
[...] 0.1748786 [...] Bit of a duhhh on my part.

:-)

But it still mystifies me why as.vector() doesn't give the same
numeric value as is obtained from just printing the object in
question.  

Because, as you said, you are printing the object. Essentially,
as.vector() just returns the value stored in the object will all
attributes stripped.  Printing on the other hand might call some
function, and here it does:

R> class(fit.gls$modelStruct$corStruct)
[1] "corAR1"    "corStruct"
R> getAnywhere("print.corAR1")
no object named ?print.corAR1? was found
R> getAnywhere("print.corStruct")
A single object matching ?print.corStruct? was found
It was found in the following places
  registered S3 method for print from namespace nlme
  namespace:nlme
with value
[...]

Also looking at the code of nlme:::coef.corStruct I would
have thought that this function would simply throw an error when
unconstrained is set to FALSE.  (Just type "nlme:::coef.corStruct"
and you'll see what I mean.)

If I set foo <- nlme:::coef.corStruct and then do

     foo(fit.gls$modelStruct$corStruct)

I do indeed get an error.  WTF?

What error did you get?  I do not get an error:

R> foo <- nlme:::coef.corStruct
R> foo(fit.gls$modelStruct$corStruct)
[1] 0.3533897

Or did you mean

R> foo(fit.gls$modelStruct$corStruct, unconstrained=FALSE)
Error in foo(fit.gls$modelStruct$corStruct, unconstrained = FALSE) : 
  do not know how to obtain parameters of ?corAR1? object

There would appear to be no way for the results to be different (look

at the code) --- but different they are.

Well, you are missing that 

R> class(fit.gls$modelStruct$corStruct)
[1] "corAR1"    "corStruct"

And S3 will try to dispatch initially on the first class, then the
second class, and so forth...

And coef.corAR1 exists in the nlme package:

R> nlme:::coef.corAR1
function (object, unconstrained = TRUE, ...) 
{
    if (unconstrained) {
        if (attr(object, "fixed")) {
            return(numeric(0))
        }
        else {
            return(as.vector(object))
        }
    }
    aux <- exp(as.vector(object))
    aux <- c((aux - 1)/(aux + 1))
    names(aux) <- "Phi"
    aux
}
<bytecode: 0x2935e60>
<environment: namespace:nlme>