Dear Doug, I assume that in the absence of a suitable anova() method, there's nothing wrong with testing the difference between two GLMM objects nested in their random effects by manually comparing the log-likelihoods as printed or returned by logLik(). Is that correct?
Regards, John -------------------------------- John Fox Department of Sociology McMaster University Hamilton, Ontario Canada L8S 4M4 905-525-9140x23604 http://socserv.mcmaster.ca/jfox -------------------------------- > -----Original Message----- > From: [EMAIL PROTECTED] > [mailto:[EMAIL PROTECTED] On Behalf Of Douglas Bates > Sent: Friday, November 26, 2004 3:13 PM > To: Spencer Graves > Cc: R-help; [EMAIL PROTECTED]; [EMAIL PROTECTED] > Subject: Re: [R] help with glmmPQL > > Spencer Graves wrote: > > HI, DOUG & JOSE: > > Is there some reason that "anova.lme" should NOT > accept an object > > of class "glmmPQL" in the example below? If you don't see > one either, > > then I suggest you consider modifying the code as described below. > > I don't think that would be appropriate. AFAIK the logLik > generic applied to a glmmPQL object does not return the > likelihood of the fitted model or an approximation to that likelihood. > > That's why there is a difference between the values of the > likelihood for the same model fit by glmmPQL and fit using > the PQL method in GLMM from the lme4 package. > > I will add anova methods for GLMM objects whenever I manage > to dig myself out of the current rewrite of the > representation of linear mixed-effects models. > > > HI, ANDREW: > > I couldn't find your data "learning" in my Windows > installation > > of R 2.0.0pat, which meant that I had to take the time to > find another > > example before I could get the error message you described. > I got it > > from modifying the example in the documentation for > "glmmPQL" as follows: > > fit1 <- glmmPQL(y ~ trt, random = ~ 1 | ID, > > family = binomial, data = bacteria) > > fit2 <- glmmPQL(y ~ trt + I(week > 2), random = ~ 1 | ID, > > family = binomial, data = bacteria) anova(fit1, > > fit2) > > Error in anova.lme(fit1, fit2) : Objects must inherit > from classes > > "gls", "gnls" "lm","lmList", "lme","nlme","nlsList", or "nls" > > > > I then checked the class of fit1 and fit2: > class(fit1) > > [1] "glmmPQL" "lme" > class(fit2) > > [1] "glmmPQL" "lme" > > As an experiment, I changed the class of fit1 and fit2: > > > class(fit1) <- "lme" > > > class(fit2) <- "lme" > > > anova(fit1, fit2) > > Model df AIC BIC logLik Test L.Ratio p-value > > fit1 1 5 1054.623 1071.592 -522.3117 > fit2 > > 2 6 1113.622 1133.984 -550.8111 1 vs 2 56.99884 <.0001 > > > > Unless someone like Doug or Jose tells us, "Don't do that", I > > would use these answers. > > I looked briefly at the code for "anova.lme", and it > looks to me > > like "glmmPQL" could be added to the list of allowed options in the > > following test roughly half way through the code: > > if (!all(match(termsClass, c("gls", "gnls", "lm", "lmList", > > "lme", "nlme", "nlsList", "nls"), 0))) { > > stop(paste("Objects must inherit from classes \"gls\", > > \"gnls\"", > > "\"lm\",\"lmList\", \"lme\",\"nlme\",\"nlsList\", or > > \"nls\"")) > > } > > > > hope this helps. spencer graves > > > > A.J. Rossini wrote: > > > >> For better or worse, it's holidays in the states. Very > amusing for > >> me being in a non-Thanksgiving celebrating country. > >> > >> In addition, it's not a problem. The complaint is valid. > Probably > >> no one has coded up the right solution yet for comparison. > >> > >> I can't recall if one would want those statistics for a binomial > >> random effects model, but I do recall some issues with model > >> comparison in that setting, though they are a bit dated > (say, 2 years > >> or so). > >> > >> On Fri, 26 Nov 2004 09:31:40 +0700, Andrew Criswell > >> <[EMAIL PROTECTED]> wrote: > >> > >> > >>> Hello: > >>> > >>> Will someone PLEASE help me with this problem. This is the third > >>> time I've posted it. > >>> > >>> When I appply anova() to two equations estimated using glmmPQL, I > >>> get a complaint, > >>> > >>> > >>> > >>>> anova(fm1, fm2) > >>>> > >>> > >>> Error in anova.lme(fm1, fm2) : Objects must inherit from classes > >>> "gls", "gnls" "lm","lmList", "lme","nlme","nlsList", or "nls" > >>> > >>> > >>> The two equations I estimated are these: > >>> > >>> > >>> > >>>> fm1 <- glmmPQL(choice ~ day + stereotypy, > >>>> > >>> > >>> + random = ~ 1 | bear, data = learning, family = > >>> binomial) > >>> > >>> > >>> > >>>> fm2 <- glmmPQL(choice ~ day + envir + stereotypy, > >>>> > >>> > >>> + random = ~ 1 | bear, data = learning, family = > >>> binomial) > >>> > >>> Individually, I get results from anova(): > >>> > >>> > >>> > >>>> anova(fm1) > >>>> > >>> > >>> numDF denDF F-value p-value > >>> (Intercept) 1 2032 7.95709 0.0048 > >>> day 1 2032 213.98391 <.0001 > >>> stereotypy 1 2032 0.42810 0.5130 > >>> > >>> > >>> > >>>> anova(fm2) > >>>> > >>> > >>> numDF denDF F-value p-value > >>> (Intercept) 1 2031 5.70343 0.0170 > >>> day 1 2031 213.21673 <.0001 > >>> envir 1 2031 12.50388 0.0004 > >>> stereotypy 1 2031 0.27256 0.6017 > >>> > >>> > >>> I did look through the archives but didn't finding > anything relevant > >>> to my problem. > >>> > >>> Hope someone can help. > >>> > >>> ANDREW > >>> ____________________________ > >>> _ > >>> platform i586-mandrake-linux-gnu > >>> arch i586 > >>> os linux-gnu > >>> system i586, linux-gnu > >>> status > >>> major 2 > >>> minor 0.0 > >>> year 2004 > >>> month 10 > >>> day 04 > >>> language R > >>> > >>> -- > >>> Andrew R. Criswell, Ph.D. > >>> Graduate School, Bangkok University > >>> > >>> ______________________________________________ > >>> [EMAIL PROTECTED] mailing list > >>> https://stat.ethz.ch/mailman/listinfo/r-help > >>> PLEASE do read the posting guide! > >>> http://www.R-project.org/posting-guide.html > >>> > >>> > >> > >> > >> > >> > >> > > > > ______________________________________________ > [EMAIL PROTECTED] mailing list > https://stat.ethz.ch/mailman/listinfo/r-help > PLEASE do read the posting guide! > http://www.R-project.org/posting-guide.html ______________________________________________ [EMAIL PROTECTED] mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
