[R-sig-ME] discrepancies in lme4
Ben Bolker
bbolker at gmail.com
Thu Jun 18 17:24:55 CEST 2015
-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1
On 15-06-17 12:00 PM, Julia Hoffmann wrote:
> The change in lme4 version 1.1-4 does explain the changes in the
> standard errors of fixed effects in our models, but there are some
> exceptions where the output suggests, that these values are not
> computed correctly.
>
> For example, I have a model where the output seems reasonably close
> to the older computations:
>
> Generalized linear mixed model fit by maximum likelihood (Laplace
> Approximation) [glmerMod] Family: gaussian ( inverse ) Formula:
> Area1_ha ~ Light * tel_round + SEX + (1 | Enclosure/animal) Data:
> kern
>
> Fixed effects: Estimate Std. Error t value Pr(>|z|) (Intercept)
> 16.745 3.611 4.638 3.52e-06 *** Light1 -1.010
> 4.023 -0.251 0.8017 tel_round2 0.442 1.695
> 0.261 0.7943 SEX1 1.818 2.944 0.618
> 0.5368
>
> Light1:tel_round2 -6.996 2.860 -2.446 0.0144 *
>
>
> But after model simplification, where only one fixed factor is
> omitted, the output is extremely different. All the calculated
> standard errors are very small and similar leading to
> unrealistically high p-values which do not correspond at all to
> older computations:
>
> Generalized linear mixed model fit by maximum likelihood (Laplace
> Approximation) [glmerMod] Family: gaussian ( inverse ) Formula:
> Area1_ha ~ Light * tel_round + (1 | Enclosure/animal) Data: kern
>
> Fixed effects: Estimate Std. Error t value Pr(>|z|) (Intercept)
> 17.748296 0.006443 2754.5 <2e-16 *** Light1
> -0.345857 0.006444 -53.7 <2e-16 *** tel_round2
> 0.455012 0.006445 70.6 <2e-16 *** Light1:tel_round2
> -7.368427 0.006445 -1143.3 <2e-16 ***
>
> Do you know what could be the problem in the calculation of the
> standard errors in the second model although it is so similar to
> the first one?
>
>
> Thank you for your help, Julia Hoffmann
Without looking at the data, I'm not sure -- I only have one guess.
For large data sets (number of obs > 10^4) we have noticed that the
relatively crude finite-difference calculation we do to estimate the
variance-covariance matrix/standard errors is sometimes unreliable.
Here's a comparison of three methods of computing standard errors;
(1) using numDeriv::hessian (the most accurate general method we
currently have available); (2) using an internal finite-difference
calculation; (3) using internal stored information (the old method).
In these cases the answers are quite similar (within a few percent),
but maybe they diverge in your case ...
library(lme4)
gm1 <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd),
data = cbpp, family = binomial)
library(numDeriv)
dd <- update(gm1,devFunOnly=TRUE)
pp <- unlist(getME(gm1,c("theta","beta")))
H <- hessian(dd,pp)
all.equal(H,H0 <- gm1 at optinfo$derivs$Hessian)
sqrt(diag(solve(H/2))[-1]) ## use numDeriv::hessian
sqrt(diag(vcov(gm1))) ## use internal finite-diff calc
## based on internal (obsolete) stored information
sqrt(diag(vcov(gm1,use.hessian=FALSE)))
>
>
> On Mon, 15 Jun 2015 09:47:08 -0400 Ben Bolker <bbolker at gmail.com>
> wrote:
>> I strongly suspect this is related to the following change in
>> lme4 version 1.1-4:
>>
>> \item Standard errors of fixed effects are now computed from the
>> approximate Hessian by default (see the \code{use.hessian}
>> argument in \code{vcov.merMod}); this gives better (correct)
>> answers when the estimates of the random- and fixed-effect
>> parameters are correlated (Github #47)
>>
>> (see https://github.com/lme4/lme4/blob/master/inst/NEWS.Rd )
>>
>> This should apply only to your glmer results, not to lmer
>> results.
>>
>> To check this, try summary(fitted_model,use.hessian=FALSE); it
>> should give you the old results.
>>
>> A further check would be to run
>> confint(fitted_model,parm="beta_"), which will give you more
>> reliable likelihood profile confidence intervals (rather than
>> relying on Wald approximations).
>>
>> Please follow up on r-sig-mixed-models at r-project.org if you have
>> further questions ...
>>
>>
>> On Mon, Jun 15, 2015 at 9:10 AM, Annika Schirmer
>> <aschirme at uni-potsdam.de> wrote:
>>> Dear Mr. Bolker, we experienced some discrepancies/problems
>>> with the lme4 package in R and would like to get your expertise
>>> on the subject. We´ve run mixed models with the functions lmer
>>> and glmer and as of late the outputs given from the summary
>>> function changed. The same models that were run a couple of
>>> month ago produce now a different output, specifically
>>> different standard errors, t values/ z values and p values and
>>> now scaled residuals are stated. Changes in the specifications
>>> of the models were not made. In most cases those differences
>>> are minimal and do not change the overall results of the
>>> models, but in some extreme cases the standard errors
>>> experienced an extreme change that led them to become very
>>> small and of the same value for all the fixed factors in the
>>> model. As a result these models now produce highly significant
>>> p values which was not the case a few month ago. Therefore the
>>> new results seem to us highly unlikely and untrustworthy. The
>>> same discrepancy in the models happens on two independent
>>> computers, therefore we exclude a general software problem. The
>>> R and lme4 versions we`re both working with are:
>>>
>>> R version 3.2.0 (2015-04-16) Platform: x86_64-w64-mingw32/x64
>>> (64-bit) Running under: Windows 8 x64 (build 9200
>>>
>>> locale: [1] LC_COLLATE=German_Germany.1252
>>> LC_CTYPE=German_Germany.1252 [3]
>>> LC_MONETARY=German_Germany.1252 LC_NUMERIC=C [5]
>>> LC_TIME=German_Germany.1252
>>>
>>> attached base packages: [1] stats graphics grDevices utils
>>> datasets methods base
>>>
>>> other attached packages: [1] lme4_1.1-7 Rcpp_0.11.6
>>> Matrix_1.2-0
>>>
>>> loaded via a namespace (and not attached): [1] minqa_1.2.4
>>> MASS_7.3-40 splines_3.2.0 nlme_3.1-120 [5] grid_3.2.0
>>> nloptr_1.0.4 lattice_0.20-31
>>>
>>> Was there a recent change in the lme4 package that could have
>>> lead to the new output? Or could there maybe be a compatibility
>>> problem between the newest version of R and lme4?
>>>
>>> We couldn`t find anything about similar problems on the
>>> internet thats why we turn directly to you, in case it might be
>>> a general problem that needs fixing or you have any idea what
>>> might be the problem.
>>>
>>> We would be very thankful for any kind of tip you could give
>>> us. Thanks in advance. Kind regards,
>>>
>>> Julia Hoffmann & Annika Schirmer
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.11 (GNU/Linux)
iQEcBAEBAgAGBQJVguLHAAoJEOCV5YRblxUHxjwIAMY5bnDvAXTSeVAKPxezvmpT
96L31nxhxePuDJsDj/JXxbCQfj+itU2dclY9WA6jWqzo42rXnQX+Ek64Oyu+ursi
eqkXWHtVJ6Fn54BKYT8czeVP+NripsFBy2VVf1awEZnsSbKGMXsu/vSwSNcykOKY
w+oG0reOsQxYk+F0i6vio755kaYNLo519m9ciIU6uv5IqncYyQvFFdpaCgPUF1yV
IeQveVx1O23EZGJ3aRD+UzaaDsr1mK0YHxycm4hW0j+DYI06pbDyw06h96oM4dE0
PGXzGMzGMWhnMeQvvd883AUub63ksh8tMgWq66vmtWK3vCGHMW9BYmzQ7nT1si8=
=boUJ
-----END PGP SIGNATURE-----
More information about the R-sig-mixed-models
mailing list