[R-sig-ME] glmer warning messages using version 1.1.7

[Sophia Kyriakou]

The bad news is that when I run the profList function I get an error
message instead

> profList <- lapply(list(default=g0.bobyqa,
+                    bobyqa=g0.bobyqa,
+                    NM=g0.NM,
+                    nlminb=g0.nlminb,
+                    LBFGSB=g0.LBFGSB),profile,dev.tol=1e-3,
+                    delta=0.1)
There were 14 warnings (use warnings() to see them)
Error in diag(vcov(object, use.hessian = use.hessian)) :
  error in evaluating the argument 'x' in selecting a method for function
'diag': Error in solve.default(h) : system is computationally singular:
reciprocal condition number = 3.04588e-17

warnings()
Warning messages:
1: In zeta(shiftpar, start = opt[seqpar1][-w]) :
  slightly lower deviances (diff=-3.55271e-15) detected
2: In nextpar(mat, cc, i, delta, lowcut, upcut) :
  Last two rows have identical or NA .zeta values: using minstep
3: In nextpar(mat, cc, i, delta, lowcut, upcut) :
  Last two rows have identical or NA .zeta values: using minstep
4: In profile.merMod(X[[1L]], ...) : non-monotonic profile
5: In zeta(shiftpar, start = opt[seqpar1][-w]) :
  slightly lower deviances (diff=-3.55271e-15) detected
6: In nextpar(mat, cc, i, delta, lowcut, upcut) :
  Last two rows have identical or NA .zeta values: using minstep
7: In nextpar(mat, cc, i, delta, lowcut, upcut) :
  Last two rows have identical or NA .zeta values: using minstep
8: In profile.merMod(X[[2L]], ...) : non-monotonic profile
9: In nextpar(mat, cc, i, delta, lowcut, upcut) :
  Last two rows have identical or NA .zeta values: using minstep
10: In profile.merMod(X[[3L]], ...) : non-monotonic profile
11: In zetafun(np, ns) : slightly lower deviances (diff=-3.55271e-15)
detected
12: In zetafun(np, ns) : slightly lower deviances (diff=-3.55271e-15)
detected
13: In nextpar(mat, cc, i, delta, lowcut, upcut) :
  Last two rows have identical or NA .zeta values: using minstep
14: In profile.merMod(X[[4L]], ...) : non-monotonic profile


 Many thanks,
Sophia

Message: 5
> Date: Thu, 12 Mar 2015 14:45:05 -0400
> From: Ben Bolker <bbolker at gmail.com>
> To: r-sig-mixed-models at r-project.org
> Subject: Re: [R-sig-ME] glmer warning messages using version 1.1.7
> Message-ID: <5501DEB1.30100 at gmail.com>
> Content-Type: text/plain; charset=windows-1252
>
> -----BEGIN PGP SIGNED MESSAGE-----
> Hash: SHA1
>
>
>    I think these are ignorable.  If you want to turn them off,
> check out ?glmerControl.
>
> Here's what I looked at:
>
> H <- g0.bobyqa at optinfo$derivs$Hessian
> V <- solve(H)
> diag(V)
> ## [1] 2.674870e+05 5.314489e-02
> sqrt(diag(V))
> ## [1] 517.1914636   0.2305318
>
> So lme4 is detecting that the surface is much flatter
> in one direction than the other (the Wald std dev is about 2500
> times greater for the first parameter (the RE std. dev) than for
> the second (the intercept).
>
> I have been working on the profiling stuff a bit in the development
> version, so your results may not be identical, but mine seem to show
> that only with Nelder-Mead (which is a bit funny since it is otherwise
> often the worst optimizer out of the bunch) do we successfully profile
> the RE standard dev.  I think I need to work on the robustness of
> the profiler a little bit more.  The good news is that we do manage
> to fit the profile, we see what we expect -- the sqrt(zeta) line is
> very flat near the origin (corresponding to low curvature/large
> uncertainty),
> but steepens as we get away from it so that the actual profile
> confidence intervals are not too big.
>
> profList <- lapply(list(default=g0.bobyqa,
>                    bobyqa=g0.bobyqa,
>                    NM=g0.NM,
>                    nlminb=g0.nlminb,
>                    LBFGSB=g0.LBFGSB),profile,dev.tol=1e-3,
>                    delta=0.1)
> profdfList <-   Map(function(x,n) {
>                           data.frame(as.data.frame(x),opt=n)
>                        },
>                        profList,names(profList))
> profdf <- do.call(rbind,profdfList)
> library("ggplot2"); theme_set(theme_bw())
> ggplot(profdf,aes(.focal,abs(.zeta),colour=opt))+
>          geom_point()+geom_line()+facet_grid(opt~.par,scale="free")+
>              geom_hline(yintercept=1.96,col="gray")
>
>
> On 15-03-12 01:56 PM, Sophia Kyriakou wrote:
> > Hello,
> >
> > I run a simulation study fitting a GLMM random intercept model and
> > I get warning messages. I tried various optimizers but they all
> > give warning messages. Please find below one of the simulated
> > samples. I am using the lme4 package version 1.1.7 on R 3.1.3. Do I
> > ignore these warning messages or is there a way of dealing with
> > them?
> >
> > Thanks in advance.
> >
> library(lme4)
> library(optimx)
> m <- 6
> q <- 20
> beta  <- 2
> sigma2 <- 0.1
> set.seed(206)
> alpha <- rnorm(q, 0, sqrt(sigma2))
> Y <- rbinom(q,m,plogis(alpha+beta))
> cl <- seq.int(q)
> tot <- rep(m,q)
> dat <- data.frame(y = Y, tot = tot, cl = cl)
> g0.default <- glmer(cbind(y, tot - y) ~ 1 + (1 | cl), data = dat,family =
>   binomial(),nAGQ = 20)
> g0.bobyqa <- glmer(cbind(y, tot - y) ~ 1 + (1 | cl), data = dat,family =
>   binomial(),nAGQ = 20,control=glmerControl(optimizer="bobyqa"))
> g0.NM <- update(g0.bobyqa,control=glmerControl(optimizer="Nelder_Mead"))
> g0.nlminb <-
>
>  update(g0.bobyqa,control=glmerControl(optimizer="optimx",optCtrl=list(method="nlminb")))
> g0.LBFGSB <-
>
> update(g0.bobyqa,control=glmerControl(optimizer="optimx",optCtrl=list(method="L-BFGS-B")))
> >
> > library(nloptr) ## from https://github.com/lme4/lme4/issues/98:
> > defaultControl <-
> > list(algorithm="NLOPT_LN_BOBYQA",xtol_rel=1e-6,maxeval=1e5)
> > nloptwrap2 <- function(fn,par,lower,upper,control=list(),...) { for
> > (n in names(defaultControl)) if (is.null(control[[n]]))
> > control[[n]] <- defaultControl[[n]] res <-
> > nloptr(x0=par,eval_f=fn,lb=lower,ub=upper,opts=control,...)
> > with(res,list(par=solution, fval=objective, feval=iterations,
> > conv=if (status>0) 0 else status, message=message)) } g0.bobyqa2 <-
> > update(g0.bobyqa,control=glmerControl(optimizer=nloptwrap2)) g0.NM2
> > <-
> >
> update(g0.bobyqa,control=glmerControl(optimizer=nloptwrap2,optCtrl=list(algorithm="NLOPT_LN_NELDERMEAD")))
> >
> >
> >
> > dd <- g0.bobyqa at optinfo$derivs
> > with(dd,max(abs(solve(Hessian,gradient)))<2e-3)
>
>

	[[alternative HTML version deleted]]