Regularized MIMIC

Joshua Pritikin and Ross Jacobucci and Timothy R. Brick

2023-11-27

Regularized MIMIC model

This example uses the immortal Holzinger Swineford data set.

library(OpenMx)
data(HS.ability.data)

The OpenMx model looks like this:

HS.ability.data$ageym <- HS.ability.data$agey*12 + HS.ability.data$agem
HS.ability.data$male <- as.numeric(HS.ability.data$Gender == 'Male')

# Specify variables
indicators <- c('visual','cubes','paper','flags','paperrev','flagssub',
                'general','paragrap','sentence','wordc','wordm')
covariates <- c("male","ageym","grade")
latents = c("g", covariates)

# Build the model
mimicModel <- mxModel(
  "MIMIC", type="RAM",
  manifestVars = indicators, latentVars = latents,

  # Set up exogenous predictors
  mxPath("one", covariates, labels=paste0('data.',covariates), free=FALSE),

  # Fix factor variance
  mxPath('g', arrows=2, free=FALSE, values=1),

  # Error variances:
  mxPath(from=c(indicators), arrows=2, free=TRUE, values=10),

  # Means (saturated means model):
  mxPath(from="one", to=indicators, values=rep(5, length(indicators))),

  # Loadings:
  mxPath(from="g", to=indicators, values=.5),

  # Covariate paths
  mxPath(covariates, "g", labels=covariates),

  # Data
  mxData(observed = HS.ability.data, type = "raw"))

# Get some good starting values for regularization. This
# saves 2-3 minutes on my laptop.
mimicModel <- mxRun(mimicModel)
## Running MIMIC with 36 parameters

Add the penalty:

mimicModel <- mxModel(
  mimicModel,
  mxMatrix('Full',1,1,free=TRUE,values=0,labels="lambda",name="hparam"),
  # Set scale to ML estimates for adaptive lasso
  mxPenaltyLASSO(what=covariates, name="LASSO",
                    scale = coef(mimicModel)[covariates],
                    lambda =  0, lambda.max =2, lambda.step=.04)
)

Run the regularization. With only three covariates, the plot of results is not very exciting. We learn that sex is not a good predictor of this factor.

regMIMIC <- mxPenaltySearch(mimicModel)
## Running MIMIC with 37 parameters
## Warning: In model 'MIMIC' Optimizer returned a non-zero status code 6. The
## model does not satisfy the first-order optimality conditions to the required
## accuracy, and no improved point for the merit function could be found during
## the final linesearch (Mx status RED)
detail <- regMIMIC$compute$steps$PS$output$detail

library(reshape2)
library(ggplot2)

est <- detail[,c(covariates, 'lambda')]
ggplot(melt(est, id.vars = 'lambda')) +
  geom_line(aes(x=lambda, y=value, color=variable)) +
  geom_vline(aes(xintercept=coef(regMIMIC)['lambda']),
             linetype="dashed", alpha=.5)

The regularized factor loadings can be found here,

detail[detail$EBIC == min(detail$EBIC), covariates]
##            male       ageym    grade
## 35 3.551246e-07 -0.02792019 1.060184

The regularization causes a lot of bias. One way to deal with this is to fix zerod parameters to zero, discard the regularization penalty, and re-fit model.

regMIMIC <- mxPenaltyZap(regMIMIC)
## Zapping 'male'
## Fixing 'lambda'
## Tip: Use
##   model = mxRun(model)
## to re-estimate the model without any penalty terms.
regMIMIC <- mxRun(regMIMIC)
## Running MIMIC with 35 parameters
summary(regMIMIC)
## Summary of MIMIC 
##  
## free parameters:
##              name matrix      row      col    Estimate   Std.Error A
## 1   MIMIC.A[1,12]      A   visual        g  2.61817539 0.359734332  
## 2   MIMIC.A[2,12]      A    cubes        g  0.93495892 0.249140391  
## 3   MIMIC.A[3,12]      A    paper        g  0.70084647 0.148839958  
## 4   MIMIC.A[4,12]      A    flags        g  1.57351355 0.485114112  
## 5   MIMIC.A[5,12]      A paperrev        g  0.99010508 0.241840086  
## 6   MIMIC.A[6,12]      A flagssub        g  3.33338049 0.638126581  
## 7   MIMIC.A[7,12]      A  general        g  9.23700633 0.532147467  
## 8   MIMIC.A[8,12]      A paragrap        g  2.53491888 0.154443963  
## 9   MIMIC.A[9,12]      A sentence        g  3.96091929 0.218438121  
## 10 MIMIC.A[10,12]      A    wordc        g  3.80400196 0.254970074  
## 11 MIMIC.A[11,12]      A    wordm        g  5.73048429 0.332716170  
## 12          ageym      A        g    ageym -0.02870552 0.005178867  
## 13          grade      A        g    grade  1.08144998 0.144108407  
## 14   MIMIC.S[1,1]      S   visual   visual 40.27135796 3.352951308 !
## 15   MIMIC.S[2,2]      S    cubes    cubes 21.00804335 1.720709660  
## 16   MIMIC.S[3,3]      S    paper    paper  7.36561698 0.605079121  
## 17   MIMIC.S[4,4]      S    flags    flags 78.47431024 6.427213845  
## 18   MIMIC.S[5,5]      S paperrev paperrev  8.35238391 0.994207854 !
## 19   MIMIC.S[6,6]      S flagssub flagssub 56.56060067 6.793459063  
## 20   MIMIC.S[7,7]      S  general  general 45.64552429 4.717436705 !
## 21   MIMIC.S[8,8]      S paragrap paragrap  4.06597962 0.402767938  
## 22   MIMIC.S[9,9]      S sentence sentence  6.80450465 0.749101684  
## 23 MIMIC.S[10,10]      S    wordc    wordc 13.88565834 1.285595220 !
## 24 MIMIC.S[11,11]      S    wordm    wordm 17.27850943 1.798287305  
## 25   MIMIC.M[1,1]      M        1   visual 21.29314806 2.978783300  
## 26   MIMIC.M[1,2]      M        1    cubes 21.38053062 1.281420562  
## 27   MIMIC.M[1,3]      M        1    paper 12.00170261 0.885460579  
## 28   MIMIC.M[1,4]      M        1    flags 13.00215459 2.293226807  
## 29   MIMIC.M[1,5]      M        1 paperrev 12.12972056 1.350150836  
## 30   MIMIC.M[1,6]      M        1 flagssub 24.45758218 4.270239961  
## 31   MIMIC.M[1,7]      M        1  general 11.26617876 9.736136751  
## 32   MIMIC.M[1,8]      M        1 paragrap  1.12587264 2.698409922  
## 33   MIMIC.M[1,9]      M        1 sentence  4.77295161 4.165463544  
## 34  MIMIC.M[1,10]      M        1    wordc 14.03580671 4.049673917  
## 35  MIMIC.M[1,11]      M        1    wordm -2.91446624 6.048415610  
## 
## Model Statistics: 
##                |  Parameters  |  Degrees of Freedom  |  Fit (-2lnL units)
##        Model:             35                   2964              17843.68
##    Saturated:             77                   2922                    NA
## Independence:             22                   2977                    NA
## Number of observations/statistics: 301/2999
## 
## Information Criteria: 
##       |  df Penalty  |  Parameters Penalty  |  Sample-Size Adjusted
## AIC:     11915.6773               17913.68                 17923.19
## BIC:       927.8025               18043.43                 17932.43
## To get additional fit indices, see help(mxRefModels)
## timestamp: 2023-11-27 13:51:23 
## Wall clock time: 3.601799 secs 
## optimizer:  SLSQP 
## OpenMx version number: 2.21.11 
## Need help?  See help(mxSummary)