esem_invariance: Measurement Invariance Testing for ESEM Models
In bifactory: (Bifactor) ESEM with Continuous (MLR) or Ordinal (WLSMV) Data

esem_invariance

R Documentation

Measurement Invariance Testing for ESEM Models

Description

Tests configural, weak (metric), strong (scalar), and strict invariance for an ESEM model across groups, optionally followed by latent variance/covariance and latent mean invariance (through). Returns a formatted table of fit indices and chi-square difference tests, analogous to Mplus's multi-group output.

Usage

esem_invariance(
  spec,
  model = c("esem", "besem"),
  missing = NULL,
  verbose = TRUE,
  through = c("strict", "varcov", "means"),
  cores = 2L,
  ...
)

Arguments

`spec`	An `esem_spec` object from `specify_model` that includes a `group` variable. If `spec$group` is `NULL`, prints a message and returns `invisible(NULL)` (not an error).
`model`	Character. Which model type to test: `"esem"` (default) for standard ESEM, or `"besem"` for Bifactor ESEM. Ordered B-ESEM uses lavaan multi-group `efa()` with explicit syntax patches (orthogonal rotation, Theta identification, threshold labels); see Scope below.
`missing`	Character. Missing data handling passed to `lavaan::cfa()`. Default `NULL` uses `spec$missing` from `specify_model` (`"pairwise"` for ordered data, `"listwise"` for continuous).
`verbose`	Logical. Print progress messages. Default `TRUE`. The `⁠Fitting <level> ...⁠` line is shown in bold blue where the console supports ANSI colour (RStudio, or a terminal whose `TERM` is not "dumb"); `options(bifactory.ansi_colour = FALSE)` turns the colour off, `TRUE` forces it.
`through`	Character. Last level of the sequence to fit. `"strict"` (default) stops at the four measurement-invariance levels. `"varcov"` adds latent variance/covariance invariance (level 5: factor variances fixed to 1 in every group and factor covariances equal across groups, i.e. 0 in every group for B-ESEM). `"means"` adds latent mean invariance (level 6: factor means fixed to 0 in every group). Levels are cumulative, so `"means"` fits all six models.
`cores`	Integer. Worker processes for the two starting-value fits of each level (lavaan's default start and `start = "simple"`), which are independent and run at the same time when `cores >= 2` (default 2, the most that helps). The results are identical; a level then takes as long as its slower start instead of both together (`psych::bfi`, three groups, six levels: 13 instead of 18 minutes). `1` fits sequentially. The workers load this package and the calling session's library path.
`...`	Additional arguments passed to the underlying fit function. Do not pass `group`, `group_equal`, `ordered`, or `parameterization` here – these are managed internally.

Value

An object of class "esem_invariance" containing:

table: Data frame with fit indices and D-statistics. print() renders it formatted.
models: Named list of fit objects, one per fitted level: configural, weak, strong, strict, and, with through, varcov and means. NULL entries indicate a model that failed to fit.
lrt: Named list of lavTestLRT() outputs, one per fitted level after configural (each compared with the level before it).
spec: The original model specification.
model: Character: "esem" or "besem".
notes: Named list, one character vector per level that raised a lavaan warning while being fitted or compared (singular information matrix for the robust test statistic, non-positive-definite latent or residual covariance matrix, negative variance). The warnings are caught, reworded and shown by print() instead of being raised, so a run that concludes does not end with a lavaan warning. Empty when nothing was raised.
fallback_from, fallback_note: If B-ESEM configural fails to converge, the result may be from an ESEM fallback (no general factor); not comparable df-for-df to Mplus B-ESEM.

Scope and caveats

Ordered B-ESEM: Validated against Mplus on BFI for G \in \{2, 3, 4\} (max |\Delta\mathrm{CFI}| \le 0.001, |\Delta\mathrm{SRMR}| \le 0.002, df match). At G \ge 5 with \ge 18 items a warning is issued; verify externally before reporting.
WLSMV SRMR: The invariance table applies an Mplus-style SRMR denominator correction for ordered models where implemented.
B-ESEM configural failure: After retries, may fall back to ESEM (fallback_from = "besem"). The printed table is structurally different from B-ESEM; do not compare to Mplus B-ESEM output.
Empty categories: Ordered multi-group fits may recode categories absent in one group (see verbose output).

See system.file("VALIDATION.md", package = "bifactory") for a summary.

Constraint mapping

Level	Continuous (MLR)	Ordered (WLSMV/Theta)
Configural	free	free
Weak	loadings	loadings
Strong	loadings + intercepts	loadings + thresholds
Strict	+ residuals	+ residuals
Latent var/cov (`through = "varcov"`)	+ factor variances (1) and covariances	+ factor variances (1) and covariances
Latent means (`through = "means"`)	+ factor means (0)	+ factor means (0)

The last two levels follow Morin, Arens and Marsh (2016): the degrees of freedom grow by k(k+1)/2 (variances and covariances of k factors) and then by k (means) per non-reference group.

Scaled chi-square difference tests

Simple subtraction of scaled chi-square values is not valid for MLR or WLSMV. lavaan::lavTestLRT() is used, which applies the Satorra-Bentler (2001) correction for MLR and a mean-variance-adjusted difference test for WLSMV. WLSMV difference test results may differ numerically from Mplus's DIFFTEST procedure.

Identification across groups

At the configural level std.lv = TRUE fixes factor variances to 1 in all groups. Under weak invariance lavaan automatically frees factor variances in Group 2+ and keeps them at 1 in Group 1. Under strong invariance factor means are freed in Group 2+ and fixed to 0 in Group 1. The optional latent var/cov level fixes the Group 2+ variances back to 1 and equates the covariances; the latent mean level fixes the Group 2+ means back to 0.

With ordered indicators (theta parameterization) the configural model uses the standard identification in every group: residual variances 1, factor means 0, all thresholds free. From the weak level on, the non-reference residual variances and factor means are free and identified through the equal loadings plus partial threshold equalities (two thresholds of one referent item per factor, one threshold of every other item), which is the Mplus convention. The same scheme is not used at the configural level: with free loadings an item's single equated threshold is its only link across groups, and when that threshold lies near zero in a non-reference group the item's residual variance is empirically unidentified (it ran to several hundred on psych::bfi, with the chi-square drifting by 14 to 36 units between programs). The standard identification is the equivalent, well-conditioned model.

Group 1 (the reference group) is the first sorted level of the grouping variable, spec$group_levels, which is also the reference group in the Mplus syntax written by run_mplus_besem_invariance. The fits pass this order to lavaan as group.label; lavaan's own default would be the order in which the groups appear in the data, and a different reference group gives an equivalent model with the same fit but differently oriented loadings and factor covariances.

Estimator paths

Continuous ESEM / B-ESEM: esem and besem use lavaan::cfa() with a native efa() block (integrated MLR estimation and rotation).
Ordered ESEM: esem(ordered = ...) calls esem_ordered. Default method = "lavaan" (lavaan WLSMV, efa() block, post-hoc rotation). method = "rotation" uses a custom DWLS pipeline (polychoric correlations + GPArotation).
Ordered B-ESEM: besem(ordered = ...) calls besem_ordered. Default method = "rotation" (custom DWLS/WLSMV + orthogonal targetT; this is what run_comparison uses for ordered data). method = "set-esem" fits a lavaan WLSMV bifactor CFA with non-primary specific loadings fixed at zero; it does not match Mplus B-ESEM loadings.
Multi-group invariance: esem_invariance uses lavaan multi-group efa() models plus explicit syntax patches for ordered B-ESEM. See that help page for scope limits.

Examples

data("HolzingerSwineford1939", package = "lavaan")

spec <- specify_model(
  Visual  = c("x1", "x2", "x3"),
  Textual = c("x4", "x5", "x6"),
  Speed   = c("x7", "x8", "x9"),
  data  = HolzingerSwineford1939,
  group = "school",
  label = "Holzinger-Swineford"
)


# Continuous ESEM measurement invariance across schools
inv <- esem_invariance(spec)
print(inv)

# Access individual model fits
summary(inv$models$strong, fit.measures = TRUE, standardized = TRUE)
lavaan::lavTestScore(inv$models$strong$lavaan_fit)

# Ordered data (WLSMV) and B-ESEM invariance follow the same pattern:
# add `ordered = TRUE` to specify_model(), then optionally `model = "besem"`.
items <- paste0("x", 1:9)
ord <- as.data.frame(lapply(HolzingerSwineford1939[, items], function(v) {
  as.integer(cut(v, breaks = quantile(v, probs = seq(0, 1, 0.2)),
                 include.lowest = TRUE))
}))
names(ord) <- items
ord$school <- HolzingerSwineford1939$school

spec_ord <- specify_model(
  Visual = c("x1", "x2", "x3"), Textual = c("x4", "x5", "x6"),
  Speed = c("x7", "x8", "x9"),
  data = ord, group = "school", ordered = TRUE
)
inv_besem <- esem_invariance(spec_ord, model = "besem")
print(inv_besem)

bifactory documentation built on Oct. 28, 2026, 1:06 a.m.