Re: [R] Function to Sort and test AIC for mixed model lme?

From: Ken Nussear <knussear_at_usgs.gov>
Date: Thu, 24 May 2007 14:23:15 -0700

On May 24, 2007, at 2:12 PM, Douglas Bates wrote:

> On 5/24/07, Ken Nussear <knussear_at_usgs.gov> wrote:

>> > Ken Nussear <knussear <at> mac.com> writes:
>> >
>> > > I'm running a series of mixed models using lme, and I wonder if
>> > there
>> > > is a way to sort them by AIC prior to testing using anova
>> > > (lme1,lme2,lme3,....lme7) other than by hand.
>> >
>> > You can try something like the following. However, also consider
>> > using dropterm or stepAIC in MASS.
>> >
>> > Dieter
>> >
>> > #---------------------
>> >
>> > library(nlme)
>> > fmlist = vector("list",2)
>> > fmlist[[1]] = lme(distance ~ age, data = Orthodont,method="ML")
>> > fmlist[[2]] = lme(distance ~ age + Sex, data = Orthodont, random
>> > = ~ 1,method="ML")
>> > aic0 = unlist(lapply(fmlist,AIC))
>> > aic0 # larger first
>> > fmlist1 = fmlist[order(aic0)]
>> > unlist(lapply(fmlist1,AIC)) # smaller first
>>
>> I looked at stepAIC, but wanted to specify my own models.
>>
>> I have come pretty close with this
>>
>> # grab all lme objects
>> tst1 <- ls(pat=".ml")
>> > tst1
>> [1] "lme.T972way.ml" "lme.T97FULL.ml" "lme.T97NOINT.ml"
>> "lme.T97NULL.ml" "lme.T97fc.ml" "lme.T97min.ml"
>> [7] "lme.T97ns.ml"
>>
>> #create array of AIC for all in tst1
>> tst2 <- lapply(lapply(tst1,get),AIC)
>> > tst2
>> [[1]]
>> [1] 507.0991
>>
>> [[2]]
>> [1] 508.7594
>>
>> [[3]]
>> [1] 564.8574
>>
>> [[4]]
>> [1] 624.3053
>>
>> [[5]]
>> [1] 502.5878
>>
>> [[6]]
>> [1] 569.8188
>>
>> [[7]]
>> [1] 504.8971
>>
>> #sort tst1 by order of tst2
>> tst3 <- tst1[order(unlist(tst2))]
>>
>> > tst3
>> [1] "lme.T97fc.ml" "lme.T97ns.ml" "lme.T972way.ml"
>> "lme.T97FULL.ml" "lme.T97NOINT.ml" "lme.T97min.ml"
>> [7] "lme.T97NULL.ml"
>>
>>
>> The problem comes with getting the final anova statement to run.
>>
>> >anova(tst3)
>> Error in anova(tst3) : no applicable method for "anova"
>>
>> I also tried
>>
>> tst4 <- paste(toString(tst3),collapse="")
>>
>> > tst4
>> [1] "lme.T97fc.ml, lme.T97ns.ml, lme.T972way.ml, lme.T97FULL.ml,
>> lme.T97NOINT.ml, lme.T97min.ml, lme.T97NULL.ml"
>> >
>> > anova(tst4)
>> Error in anova(tst4) : no applicable method for "anova"
>>
>> Any ideas on getting the last part to work?
>
>
> tst3 is a character vector of names so you would need to use
>
> do.call(anova, lapply(tst3, get))
>
> I think it is better to create a list of fitted models initially
> instead of working with names.  It would look something like this
> (this code is untested)
>
> tst2 <- lapply(tst1, get)
> names(tst2) <- tst1
> do.call(anova, tst2[order(unlist(lapply(tst2, AIC)))])



I get errors with each method that I'm not sure how to solve. Any Ideas?

Method 1

 > do.call(anova, lapply(tst3, get))
Error in `row.names<-.data.frame`(`*tmp*`, value = c("structure(list (modelStruct = structure(list(reStruct = structure(list(", :

        invalid 'row.names' length

Method 2
 > names(tst2) <- tst1
 > tst2
$lme.T972way.ml
[1] 507.0991

$lme.T97FULL.ml
[1] 508.7594

$lme.T97NOINT.ml
[1] 564.8574

$lme.T97NULL.ml
[1] 624.3053

$lme.T97fc.ml
[1] 502.5878

$lme.T97min.ml
[1] 569.8188

$lme.T97ns.ml
[1] 504.8971

 > do.call(anova, tst2[order(unlist(lapply(tst2, AIC)))]) Error in logLik(object) : no applicable method for "logLik"  >



R-help_at_stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Thu 24 May 2007 - 21:36:16 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Thu 24 May 2007 - 22:31:10 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.