Search code examples
rmixed-modelsrandom-effects

How to input model specification for glmmLasso to get expected iterations?


I'm trying to run glmmLasso on some data with the subject IDs established as the groups for random effects. This is seen in several posted examples of glmmLasso and the example soccer data set included in the CRAN package.

However, whenever I try to run my code, I end up with the following warning:

Warning message:
In est.glmmLasso.RE(fix = fix, rnd = rnd, data = data, lambda = lambda,  :
Cluster variable should be specified as a factor variable!

This is despite no variable I'm using being categorical aside from the subject IDs.

Below is an example of the code I'm trying to run:

library(lme4)
library(tidyverse)
library(tools)
library(glmmLasso)

ChoiceData <- read_csv("FatigueData.csv"))

dfs <- ChoiceData
for(i in c(2,3,5)){
  dfs[,i] <- scale(ChoiceData[,i])
}

LassodModel = glmmLasso(C ~ E_sq + R + trialNum, 
                        rnd = list(subj =~ 1),
                        data = na.omit(dfs), lambda=5, family = binomial())

summary(LassodModel)

The data I'm using can be found at this link: https://drive.google.com/file/d/1xkJulctrH1f8W6tTWSRON5x4DOL0Thkm/view?usp=sharing

I tried using as.factor(.) on my subj variable in the random effect specification, and got a separate error. So labeling it as a factor in the random effects doesn't seem to be the solution. See the error below.

LassodModel = glmmLasso(C ~ E_sq + R + trialNum, rnd = list(as.factor(subj) =~ 1),
Error: unexpected '=' in:
"LassodModel = glmmLasso(C ~ E_sq + R + trialNum, 
                        rnd = list(as.factor(subj) =" 

I've also tried a couple lambda values and haven't seen any difference.

Running the code with rnd = NULL allows it to continue, but the program doesn't seem to run through any iterations, even though trialNum should probably be removed, and the StdErr, z.value, and p.value all end up as NA, see below:

Fixed Effects:

Coefficients:
             Estimate StdErr z.value p.value
(Intercept) -0.338451     NA      NA      NA
E_sq        -1.173294     NA      NA      NA
R            1.291759     NA      NA      NA
trialNum    -0.043156     NA      NA      NA

No random effects included!

I've run the demo from the glmmLasso package, and I get a very similar outcome to above. I've reinstalled and updated packages, so theoretically there shouldn't be an issue there.

I'm really not sure what's going on.


Solution

  • I don't think you can do the conversion on the fly. Running

    dfs$subj <- factor(dfs$subj)
    

    before fitting the model takes care of the warning. As for the standard errors etc., I think you need to include final.re = TRUE ("Should the final Fisher scoring re-estimation be performed?"; perhaps the default is FALSE because this could be time-consuming, but I really don't know ...)

    summary(LassodModel)
    Call:
    glmmLasso(fix = C ~ E_sq + R + trialNum, rnd = list(subj = ~1), 
        data = na.omit(dfs), lambda = 5, family = binomial(), final.re = TRUE)
    
    
    Fixed Effects:
    
    Coefficients:
                 Estimate    StdErr  z.value p.value    
    (Intercept) -0.461363  0.308141  -1.4972  0.1343    
    E_sq        -1.627749  0.083469 -19.5012  <2e-16 ***
    R            1.790404  0.085668  20.8994  <2e-16 ***
    trialNum    -0.047837  0.062996  -0.7594  0.4476    
    ---
    Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
    
    Random Effects:
    
    StdDev:
             subj
    subj 1.593785