[R-sig-ME] Generalized least squares program

[Douglas Bates]

In general update applied to a fitted model object modifies the
original call then evaluates the modified call.  Hence there will be
no gain in speed.  The refit function in the lme4 package replaces the
response and refits the model with the new response.  If there was
such a function for a gls object that would be what you would want to
use.  However, there is no refit method for gls objects.

On Wed, Jul 16, 2008 at 4:23 PM, Gang Chen <gangchen at mail.nih.gov> wrote:
> No, update() doesn't make much difference in my case either. I've been
> reluctant to dive into the source code, but I may have to do so if no better
> choices turn up.
>
> Thanks again,
> Gang
>
>
> On Jul 16, 2008, at 4:32 PM, Steven McKinney wrote:
>
>> Hi Gang,
>>
>> Sorry, I should have said "update()" instead
>> of "predict()".
>>
>> However, a timing test on a small example shows no
>> considerable gain to updating an existing object
>> compared to refitting a new object.  This may be different
>> for your scenario.  If not, then it's an exercise of
>> diving into the code base to figure out how to adapt
>> the update.gls() function or related fitting
>> functions so you can run your extra fits without
>> recomputing everything.
>>
>>
>>> require("nlme")
>>>
>>> Ovary <- Ovary
>>> set.seed(123)
>>> Ovary$jitterFollicles <- Ovary$follicles + sample(((-2):2), nrow(Ovary),
>>> replace = TRUE)
>>>
>>> fm1 <- gls(follicles ~ sin(2*pi*Time) + cos(2*pi*Time), Ovary,
>>
>> +            correlation = corAR1(form = ~ 1 | Mare))
>>>
>>> fm2 <- update(fm1, model. = jitterFollicles ~ .)
>>>
>>> fm3 <-  gls(jitterFollicles ~ sin(2*pi*Time) + cos(2*pi*Time), Ovary,
>>
>> +            correlation = corAR1(form = ~ 1 | Mare))
>>>
>>> all.equal(fitted(fm2), fitted(fm3))
>>
>> [1] TRUE
>>>
>>> system.time(
>>
>> +             for(i in 1:100) {
>> +               fm2 <- update(fm1, model. = jitterFollicles ~ .)
>> +             }
>> +             )
>>   user  system elapsed
>>  12.022   5.238  17.264
>>>
>>> system.time(
>>
>> +             for(i in 1:100) {
>> +               fm3 <-  gls(jitterFollicles ~ sin(2*pi*Time) +
>> cos(2*pi*Time), Ovary,
>> +                           correlation = corAR1(form = ~ 1 | Mare))
>> +             }
>> +             )
>>   user  system elapsed
>>  11.898   5.227  17.126
>>>
>>
>> Steve McKinney
>>
>>
>> -----Original Message-----
>> From: Gang Chen [mailto:gangchen at mail.nih.gov]
>> Sent: Wed 7/16/2008 11:34 AM
>> To: Steven McKinney
>> Cc: r-sig-mixed-models at r-project.org
>> Subject: Re: [R-sig-ME] Generalized least squares program
>>
>> Thanks for the suggestion!
>>
>> It seems to me predict() is the opposite of what I'm looking for.
>> With an extended linear model
>>
>> Y = X*b + e, e ~ N(0, R s^2), R s^2 is the  covariance matrix of the
>> residual term e
>>
>> both X and the structure of R remain the same across all the loops in
>> my case, and Y is different and so are the estimates of b and R s^2
>> from one iteration to another.
>>
>> Gang
>>
>>
>> On Jul 16, 2008, at 1:29 PM, Steven McKinney wrote:
>>
>>> Hi Gang,
>>>
>>> Have you tried the predict() function?
>>> You should be able to hand off new data
>>> to your saved fitted model object.
>>> See
>>> ?predict.gls
>>> for an example.
>>>
>>> HTH
>>>
>>> Steve McKinney
>>>
>>> -----Original Message-----
>>> From: r-sig-mixed-models-bounces at r-project.org on behalf of Gang Chen
>>> Sent: Wed 7/16/2008 8:21 AM
>>> To: r-sig-mixed-models at r-project.org
>>> Subject: [R-sig-ME] Generalized least squares program
>>>
>>> I'm using function gls in nlme to run some extended linear model with
>>> generalized least squares many many times in a loop with the same
>>> design matrix X and same residual covariance structure but different
>>> values for the response variable. A lot of time is wasted in
>>> reformulating the same model structures and overhead matrix
>>> calculations. Such a waste makes the runtime range from one day to
>>> one month. Are there some slots available in gls that would allow me
>>> to avoid the unnecessary waste? If not, is there a counterpart of gls
>>> for extended linear model with generalized least squares available in
>>> lme4 that would give such slots?
>>>
>>> Thanks,
>>> Gang
>>>
>>> _______________________________________________
>>> R-sig-mixed-models at r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>
>>
>>
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>