[R-sig-ME] Error message
Ken Beath
ken.beath at mq.edu.au
Thu Nov 6 00:02:31 CET 2014
There is nothing I know of, which I find surprising. If you look
at Rabe-Hesketh, S., Skrondal, A., & Pickles, A. (2005). Maximum likelihood
estimation of limited and discrete dependent variable models with nested
random effects. Journal of Econometrics, 128(2), 301–323 they compare
different number of quadrature points starting at 5, and find that they
need more in situations where the clusters are small and random effects
variance high. In Rabe-Hesketh and Skrondal's book "Multilevel and
Longitudinal Modeling Using Stata" they recommend a minimum of 5 adaptive
quadrature points for binary data. Obviously not peer reviewed but this
entry from Andrew Gelmans blog is interesting
http://andrewgelman.com/2010/09/10/r_vs_stata_or_d/ The best strategy is
what I was taught with non-adaptive Gauss_Hermite and that is to increase
the quadrature points and see if it makes a difference. If it doesn't then
you have the right number of quadrature points.
One of my experience with this was fitting some unusual random effects
models to binary data using adaptive Gauss-Hermite, using routines I wrote.
I wondered why I needed 21 points, so I plotted the log-likelihood for a
cluster and one side was very steep (almost cliff-like) compared to the
other, so it just wasn't going to work.
On 6 November 2014 08:42, Ben Bolker <bbolker at gmail.com> wrote:
> On 14-11-05 04:37 PM, Ken Beath wrote:
> > I would try it using adaptive Gauss-Hermite, by setting nAgQ=3 or more
> and
> > seeing how that works. It really should be your first option when
> fitting a
> > GLMM, and something that should be checked anyway. In your case with
> binary
> > data and approx 2 per group the Laplace approximation is almost certainly
> > poor.
>
> Ken, can you point me to heuristic and/or anecdotal and/or
> (preferably) official or peer-reviewed discussions of when Laplace
> approximation is most likely to fail? (I know it fails when the
> sampling distribution of the conditional modes is non-Normal, it makes
> sense that that would occur esp. for binary data and small samples per
> group, but I'm trying to get a more precise handle on it ...)
>
> cheers
> Ben Bolker
>
> >
> > On 5 November 2014 22:55, Luciano La Sala <lucianolasala at yahoo.com.ar>
> > wrote:
> >
> >> Thank you Dan,
> >>
> >> According to the new version of lme4 I refited my model as follows:
> >>
> >> model <- glmer(Death ~ Year + Sex + Egg Volume + Hatch Order + (1|Nest
> >> ID), family = binomial, data = Data)
> >> summary(model)
> >>
> >> However, the same error message keeps showing up:
> >>
> >>
> >> Error: (maxstephalfit) PIRLS step-halvings failed to reduce deviance in
> >> pwrssUpdate
> >>
> >>
> >> Interestingly, if I reduce the model to contain only one main effect
> >> (whichever), say Hatch_Order, things look better:
> >>
> >> model2 <- glmer(Death 2 ~ Hatch Order + (1|Nest_ID), family = binomial,
> >> data = Data) summary(model2)
> >>
> >>
> >> Generalized linear mixed model fit by maximum likelihood (Laplace
> >> Approximation) ['glmerMod']
> >> Family: binomial ( logit )
> >> Formula: Death_2 ~ Hatch_Order + (1 | Nest_ID)
> >> Data: surv.2
> >>
> >> AIC BIC logLik deviance df.resid
> >> 118.5 131.8 -55.2 110.5 205
> >>
> >> Scaled residuals:
> >> Min 1Q Median 3Q Max
> >> -0.7390 -0.1714 -0.1682 -0.1506 3.7689
> >>
> >> Random effects:
> >> Groups Name Variance Std.Dev.
> >> Nest_ID (Intercept) 1.586 1.259
> >> Number of obs: 209, groups: Nest ID, 115
> >>
> >> Fixed effects:
> >> Estimate Std. Error z value Pr(>|z|)
> >> (Intercept) -3.4824 1.1274 -3.089 0.00201 **
> >> Hatch_OrderSecond -0.1266 0.7576 -0.167 0.86729
> >> Hatch_OrderThird 2.0486 0.7572 2.705 0.00682 **
> >> ---
> >> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> >>
> >> Correlation of Fixed Effects:
> >> (Intr) Htc_OS
> >> Htch_OrdrSc -0.111
> >> Htch_OrdrTh -0.709 0.276
> >>
> >>
> >> Any pointers please? Best. Luciano
> >>
> >>
> >>
> >> El 10/22/2014 6:35 PM, Daniel Wright escribió:
> >>> The lme4 package has changed some. Details are inhttp://
> >> arxiv.org/pdf/1406.5823.pdf
> >>>
> >>> For your problem, the first thing to note is glmer is now used instead
> >> of lmer for generalized linear models. Glancing at your model the other
> >> bits look like they should work.
> >>>
> >>> Dan
> >>>
> >>> Daniel B. Wright, Ph.D.
> >>> Statistical Research Division
> >>> 8701 N. MoPac Expressway, Suite 200, Austin, TX 78759
> >>> (preferred method of communication is email, use cell if urgent)
> >>> Office: 512.320.1827
> >>> Cell: 786 342 4656
> >>>
> >>>
> >>>
> >>>
> >>>
> >>>
> >>> This email message is intended only for the personal use of the
> >> recipient(s) named above. If you are not an intended recipient, you may
> not
> >> review, copy, or distribute this message. If you have received this
> >> communication in error, please notify the sender immediately by email
> and
> >> delete the original message.
> >>>
> >>>
> >>>
> >>> -----Original Message-----
> >>> From:r-sig-mixed-models-bounces at r-project.org [mailto:
> >> r-sig-mixed-models-bounces at r-project.org] On Behalf Of Luciano La Sala
> >>> Sent: Wednesday, October 22, 2014 4:20 PM
> >>> Cc:r-sig-mixed-models at r-project.org
> >>> Subject: [R-sig-ME] Error message
> >>>
> >>> Hello,
> >>>
> >>> A few years back I used to fit GLMM (binomial response) using lmer
> >> function in lme4. Back then I had to specify the family of response
> >> variable (dead /alive) as binomial. Now I have to refit those models
> using
> >> quite newer versions of both R (R x64 3.1.1) and lme4 (lme4_1.1-7), but
> >> things seem to have changed quite a bit.
> >>>
> >>> My response variable is death (yes/no), and independent variables are
> >> Year (2006 / 2007), Sex (M / F), Egg volume (continuous), and Hatching
> >> Order (ordered factor variable, namely first, second, third). I need to
> >> control autocorrelation among siblings, so I use "Nest ID" to fit random
> >> intercepts for different nests.
> >>>
> >>> My model is:
> >>>
> >>> model.1 <- lmer(Death_2 ~ Year + Sex + Egg_Volume + Hatch_Order +
> >> (1|Nest_ID), family = binomial, data = Data)
> >>> summary(model.1)
> >>>
> >>> But I get the error and warning messages below:
> >>>
> >>> Error in eval(expr, envir, enclos) :
> >>> (maxstephalfit) PIRLS step-halvings failed to reduce deviance in
> >> pwrssUpdate In addition:Warning message:
> >>> In lmer(Death_2 ~ Year + Sex + Egg_Volume + Hatch_Order + (1 |
> >> Nest_ID), :
> >>> calling lmer with 'family' is deprecated; please use glmer()
> instead
> >>>
> >>>
> >>> Question: how can I circumvent these two issues?
> >>>
> >>> Thanks in advance.
> >>>
> >>> Luciano
> >>>
> >>>
> >>> [[alternative HTML version deleted]]
> >>>
> >>> _______________________________________________
> >>> R-sig-mixed-models at r-project.org mailing listhttps://
> >> stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >>>
> >>
> >> --
> >> Luciano F. La Sala
> >> Consejo Nacional de Investigaciones Científicas y Técnicas (CONICET)
> >> Cátedra de Epidemiología
> >> Departamento de Biología, Bioquímica y Farmacia
> >> Universidad Nacional del Sur
> >> San Juan 670
> >> Bahía Blanca (8000)
> >> Argentina
> >>
> >>
> >> [[alternative HTML version deleted]]
> >>
> >>
> >> _______________________________________________
> >> R-sig-mixed-models at r-project.org mailing list
> >> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
> >>
> >>
> >
> >
>
> _______________________________________________
> R-sig-mixed-models at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
--
*Ken Beath*
Lecturer
Statistics Department
MACQUARIE UNIVERSITY NSW 2109, Australia
Phone: +61 (0)2 9850 8516
Building E4A, room 526
http://stat.mq.edu.au/our_staff/staff_-_alphabetical/staff/beath,_ken/
CRICOS Provider No 00002J
This message is intended for the addressee named and may...{{dropped:9}}
More information about the R-sig-mixed-models
mailing list