Skip to contents
library(modsem)
#> This is modsem (1.0.12). Please report any bugs!

Accelerated EM and Adaptive Quadrature

By default (as of v1.0.9) the LMS approach uses an accelerated EM procedure ("EMA") that uses Quasi-Newton and Fisher Scoring optimization steps when needed. If desireable, this can be switched to the standard Expectation-Maximization (EM) algorithm, by setting algorithm = "EM".

By default the LMS approach also uses a fixed Gauss-Hermite quadrature, to estimate a numerical approximation of the log-likelihood. Instead of a fixed quadrature, it is possible to use a quasi-adaptive quadrature instead. Due to performance reasons, the adaptive quadrature does not fit an individual quadrature to each participant, but instead one for the entire sample (at each EM iteration), based on the whole sample densities of the likelihood function. It essentially works by removing irrelevant nodes which don’t contribute to the integral, and increasing the number of nodes which actually contribute to the integral. This usually means that more nodes are placed towards the center of the distribution, compared to a standard fixed Gauss-Hermite quadrature. Using the EMA and adaptive quadrature might yield estimates that are closer to results from Mplus.

If the model struggles to converge, you can try changing the EM procedure by setting algorithm = "EMA", or algorithm = "EM", and adaptive.quad = TRUE in the modsem() function. Additionally it is possible to tweak these parameters:

  • max.iter: Maximum number of iterations for the EM algorithm (default is 500).
  • max.step: Maximum number of steps used in the Maximization step of the EM algorithm (default is 1).
  • convergence.rel: Relative convergence criterion for the EM algorithm.
  • convergence.abs: Absolute convergence criterion for the EM algorithm.
  • nodes: Number of nodes for numerical integration (default is 24). Increasing this number can improve the accuracy of the estimates, especially for complex models.
  • quad.range: Integration range for quadrature. Smaller ranges means that the integral is more focused. Applies to only when using a quasi-adaptive quadrature.
  • adaptive.frequency: How often should the quasi-adaptive quadrature be calculated? Defaults to every third EM iteration.
  • adaptive.quad.tol: Relative error tolerance when calculating the quasi-adaptive quadrature.

Here we can see an example using the TPB_UK dataaset, which is more troublesome to estimate than the simulated TPB dataset.

tpb_uk <- "
# Outer Model (Based on Hagger et al., 2007)
  ATT =~ att3 + att2 + att1 + att4
  SN =~ sn4 + sn2 + sn3 + sn1
  PBC =~ pbc2 + pbc1 + pbc3 + pbc4
  INT =~ int2 + int1 + int3 + int4
  BEH =~ beh3 + beh2 + beh1 + beh4

# Inner Model (Based on Steinmetz et al., 2011)
  INT ~ ATT + SN + PBC
  BEH ~ INT + PBC
  BEH ~ INT:PBC
"

fit <- modsem(tpb_uk,
              data = TPB_UK,
              method = "lms",
              nodes = 32, # Number of nodes for numerical integration
              adaptive.quad = TRUE, # Use quasi-adaptive quadrature
              adaptive.frequency = 3, # Update the quasi-adaptive quadrature every third EM-iteration
              adaptive.quad.tol = 1e-12, # Relative error tolerance for quasi-adaptive quad
              algorithm ="EMA", # Use accelerated EM algorithm (Default)
              convergence.abs = 1e-4, # Relative convergence criterion
              convergence.rel = 1e-10, # Relative convergence criterion
              max.iter = 500, # Maximum number of iterations
              max.step = 1) # Maximum number of steps in the maximization step
summary(fit)
#> 
#> modsem (1.0.12) ended normally after 132 iterations
#> 
#>   Estimator                                        LMS
#>   Optimization method                       EMA-NLMINB
#>   Number of model parameters                        69
#>                                                       
#>   Number of observations                          1169
#>  
#> Loglikelihood and Information Criteria:
#>   Loglikelihood                              -35375.33
#>   Akaike (AIC)                                70888.66
#>   Bayesian (BIC)                              71238.06
#>  
#> Numerical Integration:
#>   Points of integration (per dim)                   32
#>   Dimensions                                         1
#>   Total points of integration                       32
#>  
#> Fit Measures for Baseline Model (H0):
#>   Loglikelihood                              -35522.87
#>   Akaike (AIC)                                71181.74
#>   Bayesian (BIC)                              71526.09
#>   Chi-square                                   5519.01
#>   Degrees of Freedom (Chi-square)                  162
#>   P-value (Chi-square)                           0.000
#>   RMSEA                                          0.168
#>  
#> Comparative Fit to H0 (LRT test):
#>   Loglikelihood change                          147.54
#>   Difference test (D)                           295.09
#>   Degrees of freedom (D)                             1
#>   P-value (D)                                    0.000
#>  
#> R-Squared Interaction Model (H1):
#>   INT                                            0.898
#>   BEH                                            0.922
#> R-Squared Baseline Model (H0):
#>   INT                                            0.896
#>   BEH                                            0.867
#> R-Squared Change (H1 - H0):
#>   INT                                            0.002
#>   BEH                                            0.055
#> 
#> Parameter Estimates:
#>   Coefficients                          unstandardized
#>   Information                                 observed
#>   Standard errors                             standard
#>  
#> Latent Variables:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>   ATT =~        
#>     att3            1.000                             
#>     att2            0.965      0.011   86.348    0.000
#>     att1            0.812      0.017   47.181    0.000
#>     att4            0.870      0.019   45.460    0.000
#>   SN =~         
#>     sn4             1.000                             
#>     sn2             1.313      0.041   32.302    0.000
#>     sn3             1.350      0.041   32.721    0.000
#>     sn1             1.000      0.038   26.608    0.000
#>   PBC =~        
#>     pbc2            1.000                             
#>     pbc1            0.859      0.021   41.277    0.000
#>     pbc3            0.935      0.017   54.895    0.000
#>     pbc4            0.819      0.021   39.796    0.000
#>   INT =~        
#>     int2            1.000                             
#>     int1            0.970      0.011   92.124    0.000
#>     int3            0.984      0.010   98.419    0.000
#>     int4            0.992      0.009  104.574    0.000
#>   BEH =~        
#>     beh3            1.000                             
#>     beh2            0.986      0.013   77.653    0.000
#>     beh1            0.814      0.019   42.777    0.000
#>     beh4            0.803      0.019   41.556    0.000
#> 
#> Regressions:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>   INT ~         
#>     ATT            -0.061      0.030   -2.067    0.039
#>     SN              0.051      0.033    1.557    0.120
#>     PBC             1.038      0.037   28.298    0.000
#>   BEH ~         
#>     PBC             0.400      0.052    7.653    0.000
#>     INT             0.589      0.049   12.114    0.000
#>     INT:PBC         0.141      0.008   17.661    0.000
#> 
#> Intercepts:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>    .pbc2            4.012      0.066   60.988    0.000
#>    .pbc1            3.981      0.063   62.978    0.000
#>    .pbc3            3.749      0.063   59.694    0.000
#>    .pbc4            3.785      0.061   61.861    0.000
#>    .att3            3.718      0.064   58.337    0.000
#>    .att2            3.833      0.062   62.119    0.000
#>    .att1            4.206      0.060   69.989    0.000
#>    .att4            3.685      0.065   56.616    0.000
#>    .sn4             4.498      0.051   87.628    0.000
#>    .sn2             4.344      0.054   79.943    0.000
#>    .sn3             4.384      0.055   80.052    0.000
#>    .sn1             4.467      0.052   85.867    0.000
#>    .int2            3.715      0.067   55.794    0.000
#>    .int1            3.860      0.066   58.693    0.000
#>    .int3            3.732      0.066   56.373    0.000
#>    .int4            3.776      0.066   56.951    0.000
#>    .beh3            2.650      0.076   34.881    0.000
#>    .beh2            2.569      0.075   34.067    0.000
#>    .beh1            2.532      0.072   34.961    0.000
#>    .beh4            2.674      0.072   37.011    0.000
#> 
#> Covariances:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>   ATT ~~        
#>     SN              1.679      0.110   15.299    0.000
#>   PBC ~~        
#>     ATT             3.676      0.177   20.759    0.000
#>     SN              1.934      0.117   16.571    0.000
#> 
#> Variances:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>    .pbc2            0.702      0.042   16.688    0.000
#>    .pbc1            1.453      0.069   21.009    0.000
#>    .pbc3            0.803      0.043   18.615    0.000
#>    .pbc4            1.456      0.068   21.285    0.000
#>    .att3            0.296      0.023   12.807    0.000
#>    .att2            0.306      0.022   13.809    0.000
#>    .att1            1.286      0.057   22.576    0.000
#>    .att4            1.584      0.071   22.425    0.000
#>    .sn4             1.362      0.065   21.030    0.000
#>    .sn2             0.491      0.032   15.393    0.000
#>    .sn3             0.377      0.032   11.889    0.000
#>    .sn1             1.445      0.068   21.131    0.000
#>    .int2            0.237      0.014   16.970    0.000
#>    .int1            0.404      0.020   20.244    0.000
#>    .int3            0.334      0.017   19.346    0.000
#>    .int4            0.271      0.015   17.899    0.000
#>    .beh3            0.457      0.030   15.262    0.000
#>    .beh2            0.514      0.031   16.555    0.000
#>    .beh1            1.830      0.082   22.356    0.000
#>    .beh4            1.909      0.085   22.458    0.000
#>     ATT             4.450      0.197   22.605    0.000
#>     SN              1.717      0.118   14.490    0.000
#>     PBC             4.356      0.210   20.702    0.000
#>    .INT             0.506      0.038   13.263    0.000
#>    .BEH             0.448      0.034   13.217    0.000