Hi all,

I may be doing something very silly - but I'm stuck on something that has always worked fine for me!

Is there a problem with the 0.9993725-32 version of lme4 and the update() function?

I’ve been trying to use update() to drop a single fixed factor from a previous model fit, and then compare the two fits with anova(). I’ve done this successfully previously, but I can’t seem to get this to work now and am at a loss why.

Anova() returns “Error in anova(mod1, mod2) : all models must be fit to the same data object”, suggesting that the two models have been fit on different datasets which is presumably not the case, although if I re-fit the model manually, it does work?

See code and output to follow - any ideas?

######################################################
> library(lme4)
### The following dataset is a snippet of the full set for convenience
###  (the full dataset produces the same issue:)

> all_test
   fledge   mass nmass NEST
1       2 61.500     3   17
2       2 61.500     3   17
3       2 61.500     3   17
4       2 61.500     3   17
5       2 61.500     3   17
6       2 61.500     3   17
7       2 61.500     3   17
8       2 62.500     2   20
9       2 62.500     2   20
10      2 62.500     2   20
11      3 63.500     4   14
12      5 64.000     5   15
13      5 64.000     5   15
14      2 62.500     2   20
15      3 63.500     4   14
16      3 63.500     4   14
17      3 63.500     4   14
18      2 62.500     2   20
19      5 64.000     5   15
20      2 62.500     2   20
21      2 62.500     2   20
22      2 62.500     2   20
23      3 63.500     4   14
24      5 64.000     5   15
25      1 58.375     4   23
26      1 58.375     4   23
27      1 58.375     4   23
28      1 58.375     4   23
29      1 58.375     4   23
30      1 58.375     4   23

> mod1 <- lmer(fledge ~ mass + nmass + (1|NEST), family="poisson", data=all_test)
> summary(mod1)
Generalized linear mixed model fit by the Laplace approximation
Formula: fledge ~ mass + nmass + (1 | NEST)
   Data: all_test
   AIC   BIC  logLik deviance
 8.574 14.18 -0.2869   0.5737
Random effects:
 Groups Name        Variance Std.Dev.
 NEST   (Intercept)  0        0
Number of obs: 30, groups: NEST, 5

Fixed effects:
             Estimate Std. Error z value Pr(>|z|)
(Intercept) -14.26458    4.86333  -2.933  0.00336 **
mass          0.23377    0.07977   2.931  0.00338 **
nmass         0.16731    0.11564   1.447  0.14793
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
      (Intr) mass
mass  -0.996
nmass  0.271 -0.349

### .... attempting to drop 'nmass' using update() ###################
### .. Note that this only seems to work if specifying 'data=all_test'
### ..  again, I don't remember having to do this before!
### .. Also, something weird seems to be going on in the 'Data:..'
### ..  bit in the model summary?

> mod2 <- update(mod1, .~. -nmass, data=all_test)
> summary(mod2)
Generalized linear mixed model fit by the Laplace approximation
Formula: fledge ~ mass + (1 | NEST)
   Data: ..2
   AIC   BIC logLik deviance
 8.708 12.91 -1.354    2.708
Random effects:
 Groups Name        Variance   Std.Dev.
 NEST   (Intercept) 3.1459e-21 5.6088e-11
Number of obs: 30, groups: NEST, 5

Fixed effects:
             Estimate Std. Error z value Pr(>|z|)
(Intercept) -16.34354    5.15302  -3.172 0.001516 **
mass          0.27646    0.08231   3.359 0.000783 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
     (Intr)
mass -1.000

### Now, attempting to compare these two fits using anova(): #############

> anova(mod1, mod2)
Error in anova(mod1, mod2) :
  all models must be fit to the same data object

### However, if I 'refit' mod2 manually...:

> mod2 <- lmer(fledge ~ mass + (1|NEST), family="poisson", data=all_test)

### ... the two models compare fine using anova(), except for a couple of  #####
###  warning messages that I find puzzling and alarming.

> anova(mod1, mod2)
Data: all_test
Models:
mod2: fledge ~ mass + (1 | NEST)
mod1: fledge ~ mass + nmass + (1 | NEST)
     Df    AIC    BIC   logLik  Chisq Chi Df Pr(>Chisq)
mod2  3 8.7075 12.911 -1.35376
mod1  4 8.5737 14.178 -0.28686 2.1338      1     0.1441
Warning messages:
1: In deparse(expr) :
  it is not known that wchar_t is Unicode on this platform
2: In deparse(expr) :
  it is not known that wchar_t is Unicode on this platform

############################################################

Any suggestions would be very helpful, this has had me mystified for the past few days??

I’m using R version 2.10.1 (2009-12-14), on Mac OSX 10.6.2 (Snow Leopard), x86_64-apple-darwin9.8.0.

Thanks and best wishes,

Jeroen Minderman


-- 
The Sunday Times Scottish University of the Year 2009/2010
The University of Stirling is a charity registered in Scotland, 
 number SC 011159.


	[[alternative HTML version deleted]]


