[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:30:04 CET 2008
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