[R-sig-ME] Why lme4 fails on schizo dataset, but glmmTMB succeeds?

[Sun, John]

Dear All,

I am writing to ask why lme4 fails to work with random-slopes+intercepts on schizo dataset, but glmmTMB works. 
    
>## data from http://www.uic.edu/~hedeker/SCHIZX1.DAT.txt
    ## inspiration from http://www.uic.edu/~hedeker/long.html
    #may have to manually download from website.

    fn <- "schizx1.dat"
    if (!file.exists(fn)) {
        download.file("http://www.uic.edu/~hedeker/SCHIZX1.DAT.txt",dest=fn)
    }
    ## almost works: get extra ^Z, need to strip it ...
    dd0 <- read.table("schizx1.dat")
    names(dd0) <- c("id","imps79","imps79b","imps79o","int","tx",
                   "week","sweek","txswk")
    apply(dd0,2,function(x) sum(x==-9,na.rm=TRUE))
    
    ## dangerous to use na.strings==c("-9","NA") in general?
    ## what if there are legitimate -9 values in *some* columns?
    na.vals <- function(x,crit) {
        x[crit] <- NA
        return(x)
    }
    library(plyr)  ## for mutate()
    library(MASS)  ## for glmmPQL
    library(lme4)
    dd <- mutate(dd0,
                    imps79=na.vals(imps79,imps79<0),
                    imps79b=na.vals(imps79b,imps79b<0),
                    imps79b=factor(imps79b,labels=c("le mild","ge moderate")),
                    tx=factor(tx,labels=c("placebo","drug")))
    
   
    with(dd,table(tx,sweek))
    with(na.omit(dd),table(tx,sweek))
    
    ## reduce this to sequences: most data are observed at only 4 time periods
    dd$rweek <- round(dd$sweek^2)
    dd2 <- na.omit(subset(dd,rweek %in% c(0,1,3,6)))
    nrow(na.omit(dd))-nrow(dd2)  ## lose only 34 cases (out of 1600)
    dd3 <- ddply(dd2,"id",summarise,
                 tx=tx[1],
                 trans=paste(imps79b,collapse=""))
    ttab <- with(dd3,table(trans,tx))
    ttab[rowSums(ttab)>10,]
    
    m0 <- glm(imps79b~tx*sweek,dd,family=binomial)
    
    m3 <- glmer(imps79b~tx*sweek+(1|id),dd,family=binomial)
    ## extremely slow: lots of 'step-halving' action
    m4 <- glmer(imps79b~tx*sweek+(1+sweek|id),dd,family=binomial,
                verbose=100)
    ## eventually fails to converge with maxgrad=180.4 (!!)

    library(glmmTMB)
    
    m5<-glmmTMB(imps79b~tx*sweek+us(1+sweek|id),dd,family=binomial(link="logit"))
    
    summary(m5)
     Family: binomial  ( logit )
    Formula:          
    imps79b ~ tx * sweek + us(1 + sweek | id)
    Data: dd
    
         AIC      BIC   logLik deviance df.resid 
      1173.7   1211.4   -579.9   1159.7     1596 
    
    Random effects:
    
    Conditional model:
     Groups Name        Variance Std.Dev. Corr  
     id     (Intercept) 2184.9   46.74          
            sweek        365.9   19.13    -0.81 
    Number of obs: 1603, groups:  id, 437
    
    Conditional model:
                 Estimate Std. Error z value
    (Intercept)   45.2732     7.6836   5.892
    txdrug         0.8314     6.9096   0.120
    sweek        -14.6736     3.1210  -4.702
    txdrug:sweek  -7.3960     2.9027  -2.548
                 Pr(>|z|)    
    (Intercept)  3.81e-09 ***
    txdrug         0.9042    
    sweek        2.58e-06 ***
    txdrug:sweek   0.0108 *  
    ---
    Signif. codes:  
    0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

https://github.com/lme4/lme4/blob/master/misc/issues/schizo.R

Best regards,
John