estimation alternatives using the LMS approach
Source:vignettes/estimation_lms.Rmd
estimation_lms.RmdAccelerated 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 desirable, this can be
switched to the standard Expectation-Maximization (EM) algorithm, by
setting algorithm = "EM".
By default the LMS approach also uses a quasi-adaptive Gauss-Hermite
quadrature, to estimate a numerical approximation of the log-likelihood.
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. A fixed Gauss-Hermite
quadrature can be used instead, by setting
adaptive.quad = FALSE.
If the model struggles to converge, you can try changing the EM
procedure by setting algorithm = "EMA" or
algorithm = "EM", switching between the adaptive and fixed
quadrature (adaptive.quad), or increasing the number of
nodes (nodes) 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 (default is1e-10). -
convergence.abs: Absolute convergence criterion for the EM algorithm (default is1e-4). -
convergence.grad: Gradient convergence criterion (default is1e-3). When the absolute or relative criterion is reached, estimation only stops if the largest absolute gradient of the log-likelihood is belowconvergence.grad. Otherwise a Fisher scoring step is taken and the EM algorithm continues, for at mostem.control$grad.max.steps(default is 25) such steps. Set toInfto disable the check. -
nodes: Number of nodes for numerical integration (default is 24). Increasing this number can improve the accuracy of the estimates, especially for complex models. -
adaptive.quad: Should a quasi-adaptive quadrature be used (default isTRUE)? IfFALSE, a fixed Gauss-Hermite quadrature is used. -
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 (default is1e-12). -
em.control: A list of control parameters for the EM and EMA algorithms (e.g.,grad.max.steps), for fine-tuning the estimation procedure. -
n.threads: Number of threads used for computing the densities of the likelihood function. By default (NULL), at most two threads are used, while"max"uses all available threads.
Here we can see an example using the TPB_UK dataset,
which is more troublesome to estimate than the simulated
TPB dataset.
tpb_uk <- "
# Outer Model (Based on Hagger et al., 2007)
ATT =~ att1 + att2 + att3 + att4 + att5
SN =~ sn1 + sn2 + sn3
PBC =~ pbc1 + pbc2 + pbc3 + pbc4
INT =~ int1 + int2 + int3 + int4
BEH =~ beh1 + beh2 + beh3
# 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, # Absolute convergence criterion
convergence.rel = 1e-10, # Relative convergence criterion
convergence.grad = 1e-3, # Gradient 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.24) ended normally after 84 iterations
#>
#> Estimator LMS
#> Optimization method EMA-NLMINB
#> Number of model parameters 66
#>
#> Number of observations 1894
#>
#> Loglikelihood and Information Criteria:
#> Loglikelihood -53953.56
#> Akaike (AIC) 108039.12
#> Bayesian (BIC) 108405.19
#>
#> Numerical Integration:
#> Standard Effective
#> Points of integration (per dim) 32 339
#> Dimensions 1 1
#> Total points of integration 32 339
#>
#> Fit Measures for Baseline Model (H0):
#> Standard
#> Chi-square 3197.68
#> Degrees of Freedom (Chi-square) 144
#> P-value (Chi-square) 0.000
#> RMSEA 0.106
#>
#> Loglikelihood -54095.04
#> Akaike (AIC) 108320.07
#> Bayesian (BIC) 108680.59
#>
#> Comparative Fit to H0 (LRT test):
#> Loglikelihood change 141.48
#> Difference test (D) 282.95
#> Degrees of freedom (D) 1
#> P-value (D) 0.000
#>
#> R-Squared Interaction Model (H1):
#> INT 0.914
#> BEH 0.874
#> R-Squared Baseline Model (H0):
#> INT 0.914
#> BEH 0.829
#> R-Squared Change (H1 - H0):
#> INT 0.000
#> BEH 0.046
#>
#> Parameter Estimates:
#> Coefficients unstandardized
#> Information observed
#> Standard errors standard
#>
#> Latent Variables:
#> Estimate Std.Error z.value P(>|z|)
#> ATT =~
#> att1 1.000
#> att2 1.327 0.028 47.080 0.000
#> att3 1.394 0.029 47.695 0.000
#> att4 0.920 0.028 33.186 0.000
#> att5 0.992 0.029 33.862 0.000
#> SN =~
#> sn1 1.000
#> sn2 1.243 0.026 46.973 0.000
#> sn3 1.288 0.027 48.203 0.000
#> PBC =~
#> pbc1 1.000
#> pbc2 3.958 0.282 14.048 0.000
#> pbc3 3.661 0.262 13.990 0.000
#> pbc4 1.313 0.115 11.416 0.000
#> INT =~
#> int1 1.000
#> int2 1.040 0.010 107.633 0.000
#> int3 1.020 0.010 100.967 0.000
#> int4 1.028 0.010 104.085 0.000
#> BEH =~
#> beh1 1.000
#> beh2 1.203 0.016 73.486 0.000
#> beh3 1.226 0.016 78.279 0.000
#>
#> Regressions:
#> Estimate Std.Error z.value P(>|z|)
#> INT ~
#> ATT 0.283 0.032 8.944 0.000
#> SN -0.023 0.021 -1.102 0.271
#> PBC 3.069 0.239 12.829 0.000
#> BEH ~
#> PBC 0.675 0.168 4.019 0.000
#> INT 0.658 0.041 16.193 0.000
#> INT:PBC 0.390 0.036 10.818 0.000
#>
#> Intercepts:
#> Estimate Std.Error z.value P(>|z|)
#> .pbc1 5.783 0.038 151.255 0.000
#> .pbc2 4.084 0.050 81.490 0.000
#> .pbc3 3.837 0.048 80.077 0.000
#> .pbc4 5.378 0.039 137.171 0.000
#> .att1 4.359 0.046 95.460 0.000
#> .att2 3.835 0.048 79.202 0.000
#> .att3 3.797 0.050 75.442 0.000
#> .att4 4.963 0.045 111.421 0.000
#> .att5 4.986 0.047 105.961 0.000
#> .sn1 5.012 0.040 126.687 0.000
#> .sn2 4.397 0.041 106.010 0.000
#> .sn3 4.427 0.042 106.443 0.000
#> .int1 3.948 0.050 79.120 0.000
#> .int2 3.848 0.051 76.048 0.000
#> .int3 3.845 0.050 76.706 0.000
#> .int4 3.885 0.050 77.256 0.000
#> .beh1 2.380 0.052 45.884 0.000
#> .beh2 2.859 0.060 47.364 0.000
#> .beh3 2.910 0.060 48.117 0.000
#>
#> Covariances:
#> Estimate Std.Error z.value P(>|z|)
#> ATT ~~
#> SN 1.396 0.070 20.086 0.000
#> PBC ~~
#> ATT 0.695 0.057 12.291 0.000
#> SN 0.496 0.042 11.895 0.000
#>
#> Variances:
#> Estimate Std.Error z.value P(>|z|)
#> .pbc1 2.494 0.081 30.627 0.000
#> .pbc2 0.459 0.026 17.788 0.000
#> .pbc3 0.672 0.029 23.392 0.000
#> .pbc4 2.438 0.080 30.500 0.000
#> .att1 1.656 0.056 29.358 0.000
#> .att2 0.400 0.020 20.171 0.000
#> .att3 0.341 0.020 16.981 0.000
#> .att4 1.815 0.062 29.337 0.000
#> .att5 1.937 0.066 29.205 0.000
#> .sn1 1.151 0.041 28.031 0.000
#> .sn2 0.459 0.025 18.024 0.000
#> .sn3 0.269 0.024 11.251 0.000
#> .int1 0.457 0.017 26.405 0.000
#> .int2 0.240 0.011 21.835 0.000
#> .int3 0.327 0.013 24.517 0.000
#> .int4 0.287 0.012 23.415 0.000
#> .beh1 0.722 0.028 25.886 0.000
#> .beh2 0.599 0.027 21.871 0.000
#> .beh3 0.408 0.023 17.447 0.000
#> ATT 2.293 0.118 19.481 0.000
#> SN 1.813 0.090 20.065 0.000
#> PBC 0.274 0.040 6.886 0.000
#> .INT 0.368 0.022 17.119 0.000
#> .BEH 0.462 0.023 19.779 0.000When using the quasi-adaptive quadrature (default), the
Numerical Integration section of the summary reports two
sets of node counts. Standard is the requested number of
nodes (nodes), which is roughly the number of nodes kept
after pruning. Effective is the number of nodes you would
need with a fixed Gauss-Hermite quadrature
(adaptive.quad = FALSE) to get effectively the same
results. The same information can be retrieved using
modsem_inspect(fit, "quad").