[R-sig-ME] Binomial problem for lme4_0.999375-1: 0,1 vs xbar n

Douglas Bates bates at stat.wisc.edu
Wed Jan 23 21:38:54 CET 2008


It appears that the problem is more with the handling of the weights
argument than with the binomial per se.  Look at the examples in

?cbpp

The first example runs properly but the second, which is inside a
\dontrun block, doesn't.

On Jan 23, 2008 2:30 PM, Douglas Bates <bates at stat.wisc.edu> wrote:
> Thanks for the report, Van.  I have entered it in the bug tracker list.
>
> I had thought that I had checked such results for consistency but
> apparently I didn't.
>
>
> On Jan 23, 2008 2:14 PM, Parsons, Van L. (CDC/CCHIS/NCHS) <vlp1 at cdc.gov> wrote:
> > Hello,
> >   The lme4_0.999375-1 lmer run on a  full response (0,1) dataset
> >   is not consistent with the run on an equivalent, but condensed,
> >   (pbar,n) dataset for the binomial.
> >
> >   For lme4_0.99875-9, the outputs appear to be order-of-magnitude
> >   consistent, so the alpha version seems to have the problem.
> >   ( I am focusing on the random and fixed effects output )
> >
> >   Here are some sample code and outputs generated from sleepstudy
> >   in lme4 using 2 datasets and both new and old lme4.
> >
> >   Thanks,
> >   Van
> >
> > #=======================================================
> > # R code
> > library(lme4)
> > library(doBy)
> > sessionInfo()
> >
> >
> >  sleep1 = sleepstudy   # from lme4 pkg
> >
> >      # add binary variable
> >  sleep1$p1 = ifelse(sleep1$Reaction < mean(sleep1$Reaction),0,1)
> >
> >    # condense binary to mean and count within Subject
> >    # use doBy package
> > sleep1c= summaryBy(p1~Subject, data = sleep1,FUN = c(mean,length) )
> >    names(sleep1c)[ c(2,3)] = c("p1", "n1")
> >
> >
> >
> >  #--------------------------------------------
> >   # simple model
> >
> > mod1  =p1~ 1 + (1| Subject)
> >
> > # lmer run   on full   data
> > mfull   =  lmer(mod1, data=sleep1,  family=binomial )
> > summary(mfull)
> >
> > # lmer run   on   condensed data
> > mcond   =  lmer(mod1, data=sleep1c, weights= n1,family=binomial )
> > summary(mcond)
> >
> > #--------------------------------------------
> > #====================================================================
> >
> >
> >
> > ## 4  partial outputs
> >
> > #====================================================================
> >
> > lme4_0.999375-1  full data
> > mfull   =  lmer(mod1, data=sleep1,  family=binomial )
> >    Data: sleep1
> >    AIC   BIC logLik deviance
> >  238.5 244.9 -117.2    234.5
> > Random effects:
> >  Groups  Name        Variance Std.Dev.
> >  Subject (Intercept) 1.0231   1.0115
> > Number of obs: 180, groups: Subject, 18
> >
> > Fixed effects:
> >             Estimate Std. Error z value Pr(>|z|)
> > (Intercept)  -0.3213     0.2887  -1.113    0.266
> > #====================================================================
> > # this next  run is not consistent with the other 3
> >
> > lme4_0.999375-1  condensed  data
> >
> > mcond   =  lmer(mod1, data=sleep1c, weights= n1,family=binomial )
> >    Data: sleep1c
> >    AIC   BIC logLik deviance
> >  9.622 11.40 -2.811    5.622
> > Random effects:
> >  Groups  Name        Variance Std.Dev.
> >  Subject (Intercept)  0        0
> > Number of obs: 18, groups: Subject, 18
> >
> > Fixed effects:
> >             Estimate Std. Error z value Pr(>|z|)
> > (Intercept)  -0.2457     0.4750 -0.5173    0.605
> > #====================================================================
> > lme4_0.99875-9 full
> >
> > > mfull   =  lmer(mod1, data=sleep1,  family=binomial )
> >    AIC   BIC logLik deviance
> >  238.5 244.9 -117.2    234.5
> > Random effects:
> >  Groups  Name        Variance Std.Dev.
> >  Subject (Intercept) 1.0420   1.0208
> > number of obs: 180, groups: Subject, 18
> >
> > Estimated scale (compare to  1 )  0.9369863
> >
> > Fixed effects:
> >             Estimate Std. Error z value Pr(>|z|)
> > (Intercept)  -0.3309     0.2906  -1.139    0.255
> >
> >
> >
> > #====================================================================
> > lme4_0.99875-9 condensed
> >  mcond   =  lmer(mod1, data=sleep1c, weights= n1,family=binomial )
> >
> >    AIC   BIC logLik deviance
> >  47.88 49.66 -21.94    43.88
> > Random effects:
> >  Groups  Name        Variance Std.Dev.
> >  Subject (Intercept) 1.0241   1.012
> > number of obs: 18, groups: Subject, 18
> >
> > Estimated scale (compare to  1 )  1.012022
> >
> > Fixed effects:
> >             Estimate Std. Error z value Pr(>|z|)
> > (Intercept)  -0.3183     0.2888  -1.102    0.270
> >
> >
> >
> > #====================================================================
> > > sessionInfo()
> > R version 2.6.1 (2007-11-26)
> > i386-pc-mingw32
> >
> > locale:
> > LC_COLLATE=English_United States.1252;
> > LC_CTYPE=English_United States.1252;
> > LC_MONETARY=English_United States.1252;
> > LC_NUMERIC=C;LC_TIME=English_United States.1252
> >
> > attached base packages:
> > [1] stats     graphics  grDevices utils   datasets  methods   base
> >
> > other attached packages:
> > [1] doBy_2.1    lme4_0.999375-1   Matrix_0.999375-4 lattice_0.17-4
> >
> > loaded via a namespace (and not attached):
> > [1] cluster_1.11.9 grid_2.6.1     Hmisc_3.4-3
> >
> >
> >
> > _______________________________________________
> > R-sig-mixed-models at r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >
>




More information about the R-sig-mixed-models mailing list