rlme4

How to ensure that modular lme4 takes the control argument?


I have been using the modular lme4 functions to customize some of our problems and I noticed that the control argument seems to be ignored in the modular functions or I am just not placing it in the right place. Here a minimal example based on the rpubs from https://rpubs.com/bbolker/groupmembers

Can you help me find what am I missing that I can't change the optimizer using the modular lme4 functions?

nm <- 20
nobs <- 500
set.seed(101)
## choose items for observations
pres <- matrix(rbinom(nobs*nm, prob=0.25, size=1), nrow=nobs, ncol=nm)
dimnames(pres) <- list(NULL, LETTERS[seq(nm)])
pres[1:5, ]

b <- rnorm(nm)  ## item-level effects
## n.b. get in trouble if we don't add residual error
## (theta is scaled relative to residual error)
## here, use theta=sigma=1
y <- c(pres %*% b) + rnorm(nobs, sd=1)

## helpful to specify a factor with the right levels:
## actual values are unimportant since we will specify Zt/Ztlist directly
fake <- rep(LETTERS[seq(nm)], length.out=nobs)
lmod <- lFormula(y ~ 1 + (1 | fake), control=lmerControl(
  optimizer="bobyqa", 
  optCtrl=list(maxeval=1)
))
lmod$reTrms$Zt <- lmod$reTrms$Ztlist[[1]] <- Matrix(t(pres))
devfun <- do.call(mkLmerDevfun, lmod)
opt <- optimizeLmer(devfun)
m1 <- mkMerMod(environment(devfun), opt, lmod$reTrms, fr=lmod$fr)
summary(m1)
m1@optinfo$optimizer # does not coincide with requested optimizer
# [1] "nloptwrap"

Solution

  • tl;dr make sure to include the control information in the list of arguments passed to mkLmerDevfun (and optimizeLmer). If we look at the arguments of mkLmerDevfun:

    args(mkLmerDevfun)
    function (fr, X, reTrms, REML = TRUE, start = NULL, verbose = 0, 
        control = lmerControl(), ...) 
    

    and if we look at how lmer itself calls mkLmerDevfun (here):

    devfun <- do.call(mkLmerDevfun,
                      c(lmod,
                        list(start=start, verbose=verbose, control=control)))
    

    we can see that start, verbose, and control are not carried over from lFormula. Logically, this is because these arguments control the way the optimization is done, and not the construction of the data structures needed to fit the model.

    So you need to include the start, verbose, and control components (if different from the default) at the mkLmerDevfun() step of modular fitting — and in fact it's slightly worse than that, some arguments have to be supplied at the optimizeLmer step as well1:

    ## we have to use 'maxfun' and not 'maxeval' for bobyqa()
    ctrl <-lmerControl(optimizer="bobyqa", optCtrl=list(maxfun=1))
    devfun <- do.call(mkLmerDevfun, c(lmod, list(control = ctrl)))
    opt <- optimizeLmer(devfun, optimizer = ctrl$optimizer,
                        control = ctrl$optCtrl)
    

    1: this is a little messy for historical reasons