search

singmann / afex

Analysis of Factorial EXperiments (R package)
120 stars 32 forks source link

Cannot increase number of iterations from where the GLMM models that will not converge left off during likelihood ratio testing #80

Open TJWukitsch opened 4 years ago

TJWukitsch commented 4 years ago

Hey Dr. Singmann,

I really like afex. It has been helpful so far with getting overall tests for categorical variables with more than 2 levels except when I have a GLMM that I need to use more iterations on to get it to converge. For example, I usually run the following code on my GLMER objects that have trouble converging and may have to run it a second time (I have yet to need a third to reach convergence).

 Savers <-glmer(Total.Aversive ~ c.conc*Sex*Condition*Environment
               + (c.conc|RatID), data=mydsuc, family=poisson)

 # Model did not converge, used code below to extend # of iterations and start from where the previous model left off.
ss2 <- getME(Savers,c("theta","fixef"))
Savers <- update(Savers,start=ss2,control=glmerControl(optCtrl=list(maxfun=2e9)))

#Model did not converge again, Rerun above code for further iterations of the function picking up where the previous model left off and overwriting ss2 with updated variables.
ss2 <- getME(Savers,c("theta","fixef"))
Savers <- update(Savers,start=ss2,control=glmerControl(optCtrl=list(maxfun=2e9)))

 summary(Savers)

No problem there. The model converges after the second iteration extension from where the last model left off. However, when I run afex::mixed to get the overall estimates for my GLMM effects with 3+ level categorical variables involved I run into problems with convergence that I cannot seem to correct.

I am running an LRT because I want to use REML due to its help with unequal cell sizes and not fully completed repeated measures by all subjects that are intrinsic to the paradigm I am working with. I have tried several things and this post would get very long if I explained fully. Briefly:

Is there a way to start from where the previous models left off (theta & fixef) and extend the number of iterations of all the models run during the LRT and rerun the LRT with those updated? If so, do you have some example code? If not, do you have any suggestions? I am not a newbie but not a super advanced user in R. Any help is greatly appreciated here as I imagine others will encounter this issue.

singmann commented 4 years ago

Can you provide a working example? Otherwise it will be difficult for me to provide a fix.

SDAMcIntyre commented 4 years ago

I had the same/similar problem, so here is a reprex.

# generate data
set.seed(2)
n <- 10
randomeffect <- runif(n, min = 0.2, max = 0.8)
fixedandrandom <- c(rowMeans(cbind(rep(0.5,n), randomeffect)), 
                    rowMeans(cbind(rep(0.8,n), randomeffect)))
df <- data.frame(ID = formatC(rep(1:n,2), width = 2, flag = '0'),
                 time = c(rep('early',n), rep('late',n)),
                 pCorrect = rbinom(n*2, 100, fixedandrandom)/100,
                 nTotal = rep(100, n*2))
# fit mixed model
library(afex)
set_sum_contrasts()

m1 <- mixed(pCorrect ~ time + (1 + time|ID),
            data = df, method = 'LRT',
            family = binomial, weights = df$nTotal)

Fitting 2 (g)lmer() models: [.boundary (singular) fit: see ?isSingular .] Warning message: In checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model failed to converge with max|grad| = 0.00306567 (tol = 0.001, component 1)

# get estimates
ss <- getME(m1$full_model, c('theta', 'fixef'))

# update, try again, starting at previous estimates
m2 <- update(m1, start = ss)

Fitting 2 (g)lmer() models: [boundary (singular) fit: see ?isSingular .Error in getStart(start, lower = rho$lower, pred = rho$pp, "theta") : incorrect number of fixef components (!=1)

length(ss$fixef)

[1] 2

# try to give it fixef length 1 it asked for
# not sensible, but to illustrate something weird is happening with the length
ss$fixef <- ss$fixef[-1]
length(ss$fixef)

[1] 1

m2 <- update(m1, start = ss)

Fitting 2 (g)lmer() models: [Error in getStart(start, lower = rho$lower, pred = rho$pp, "theta") : incorrect number of fixef components (!=2)

SDAMcIntyre commented 4 years ago

This is not limited to method = 'LRT' by the way, also happens with 'PB'.

TJWukitsch commented 4 years ago

This is not limited to method = 'LRT' by the way, also happens with 'PB'.

@SDAMcIntyre I found that update attempts with 'PB' also did not work. Thank you for the working example. I am working with proprietary data and thus cannot provide my working example and have had limited success generating data sets with such issues mostly due to my current skill level.