Structural equation modeling and confirmatory factor analysis – method notes #

The decisions behind the structural equation modeling and confirmatory factor analysis module – which construction, estimator, default, threshold, refusal or fallback stands behind each card, and the numbers that validate them – on simulated data where the note says so, and otherwise on the example dataset (cafe.csv, its spec beside it). Sections follow the manual’s.

Model specification #

The matrices are a projection of the syntax box, never the other way round. The parsed text is the model: every widget edit is written into the text and read back from lavaan’s own parameter table, so the matrices cannot show a model the text does not say. That is why they freeze on a syntax error instead of clearing – they still describe the last text that parsed, and an edit made to a stale projection would be discarded by the next successful parse – and why the module checks on load that lavaan reads its generated syntax back into the parameter table it was written from. A mismatch means the projection lies, so the widgets are hidden while the text, which is what lavaan fits, keeps running; the same round trip runs after every widget write, so a drift between the formatter and lavaan surfaces as a banner rather than as a fitted model that differs from the one on screen.

Column names that aren’t lavaan identifiers are aliased, and the alias is what the buffer holds. lavaan has no quoting mechanism: a backtick is a parse error, a space splits a name and a-b parses as b. So a header that isn’t a bare identifier is rewritten – every character outside letters, digits, _ and . to _, an X prefixed when the result doesn’t open with a letter (or with a dot no digit follows), a . appended to an R reserved word, and _2, _3, … appended until the alias is unique – and the syntax box, the completions and the diagram all speak the alias, while the original header still resolves when typed. Internally the fit runs on positional names, which is what lets a rename in the data view be carried into the text without the model being re-read.

Every generated name is minted against one namespace. In lavaan a label is a parameter’s name, and a label written twice equates the two parameters; a label that is also a column name reads as a column in a ~ term and as a parameter in a := expression. So the helper that labels mediation paths (a, b, c, …, then m1, m2, …), the one that names defined parameters (indirect, indirect2, …) and the one that names an interaction (XZ, XZ2, …) all draw on the union of the variable names, their aliases, the factor names, every label on a loading, path, covariance or intercept, and every := name. A suggested name is bumped past a collision; a typed one is refused, since the user meant that name and silently equating two parameters is the one outcome worse than a toast.

The mediation helper always writes the direct path, because a*b is the indirect effect only while it is free. Without an X→Y regression lavaan fits c′ = 0, and under that imposed full mediation a*b is the total effect, not the indirect one. The helper therefore adds X→Y when it is missing and emits direct and total beside the indirect line. When X→Y is fixed by the user, or differs across groups so that no single label names it, the indirect line is emitted alone and the toast says so – folding a fixed coefficient into a := term would silently drop the constraint – and a fixed or group-wise path on the indirect route stops the helper for the same reason.

The serial total sums every route, and is withheld past six mediators. A saturated forward chain reaches Y through every ordered subset of the mediators, not only the full one, so totalIndirect sums all 2^k − 1 products and total adds the direct path; summing the named chain alone would understate the total effect. At six mediators that expression already has 63 terms, and past six it stops being readable, so the form names the chain alone and points at the effect decomposition, which prints every route as a row of its own. The forward set is saturated in the first place so that the skip routes are estimated rather than fixed at zero by omission – which is what a := sum over them presupposes.

Product indicators, not LMS. A latent interaction is fitted the way lavaan can fit one: as an ordinary latent whose indicators are products of the two components’ indicators (Kenny & Judd, 1984; Marsh, Wen & Hau, 2004), the unconstrained variant that leaves the product loadings and variances free. Latent moderated structural equations (LMS) and QML estimate the interaction without forming products but need an engine outside lavaan, so they are not offered. The construction is kept latent rather than reduced to a product of scale scores because the reliability of a product is roughly the product of its components’ reliabilities: an interaction is the estimand measurement error damages hardest, and the one place where collapsing constructs to scores costs most.

Matched pairs follow the measurement order, and are the default. Marsh, Wen & Hau pair indicators by reliability rank; nothing knows a reliability before the fit, so the order the measurement matrix lists each factor’s indicators stands in for it, and reordering the indicators is how another pairing is imposed. Matched pairing uses each indicator once, so no two products share a component and no residual covariance is needed – the reason it is the default. Under all pairs two products sharing a component covary through it whatever the interaction does, and leaving that in their residuals would load it onto the factor, so the ~~ row between every such pair is written with the definition. An interaction needs two product indicators to be identified, like any latent, so a plan that yields fewer is refused.

Double mean centering is the default; orthogonalizing is refused on collinear components. Double mean centering (Lin, Wen, Marsh & Lin, 2010) centers the indicators, multiplies, and centers the products again, so the interaction needs no mean structure and its indicators have mean zero; plain mean centering leaves the product with the components’ covariance as its mean, which the model has to absorb. Orthogonalizing (Little, Bovaird & Widaman, 2006) regresses each raw product on its two raw components and keeps the residual, so the product is uncorrelated with both by construction; when a pair’s two component columns are perfectly collinear the normal equations are singular and there is nothing to separate, so the run is refused with a message naming double mean centering as the alternative.

The product columns are built when the model runs and never enter the dataset. Centering has to follow the sample that reaches the fit: a data-transformation rule would center over the whole column, rows that listwise deletion then drops included, and the products would be centered on a mean the fit never sees. So the columns are computed into the fit’s data from the rows where both components are present, and a row missing either carries no product. What is auditable is the definition, kept in the buffer as the # @interaction comment beside the =~ line: nothing in a parameter table distinguishes XZ =~ XZ.1 + XZ.2 from a factor measured by two columns that happen to be absent, so the definition rides in the one channel the parser leaves untouched and round-trips with the comments.

A term written in only some group: blocks is rewritten into c(...) form. A grouped parameter table states every term once per group, and the c(...) spelling the widgets write back has no way to say “absent in this group” – so a term one block omits comes back present in every group, which lavaan may fit as a different model than the blocks described. The projection cannot refuse it without refusing the paste, so it rewrites and warns: a dismissible banner names the statements it did that to and asks for the free-parameter count to be checked before running.

Options #

Ordinal indicators substitute WLSMV rather than letting lavaan error. lavaan refuses both halves of a mismatch: the ML family and GLS error on ordered=, and the five least-squares estimators (WLSMV, ULSMV, DWLS, ULS, WLS) error on continuous data unless ordered = FALSE is passed. So the pair is resolved before the call – WLSMV for any other estimator when the model names a column typed Ordinal, ordered = FALSE for a least-squares estimator on continuous columns – and the summary reports the substitution beside the estimator, since lavaan itself never upgrades. The scope is the model text, not the selection: one ordinal column is enough to put the whole fit into categorical mode, so a selected column the model never mentions must not reach ordered=, and the ordered set is read from lavaan’s own parameter table once a parse has landed, from the text before that.

FIML follows the resolved estimator. missing = "ml" is an ML-family option: lavaan errors on it under any other estimator (missing=“ml” is not allowed for estimator “MLM”, and not available in the categorical setting once ordered= is on). MLR is the robust ML that takes it, so MLM, MLMVS and MLMV are disabled while FIML is selected and a selection among them switches to MLR; a FIML request under a least-squares estimator falls back to pairwise deletion, which that family accepts, and under GLS to lavaan’s listwise default. Because the estimator can itself be substituted for ordinal indicators, a plain ML + FIML pick can land on pairwise without the user touching either control – so the summary reports the mode that ran with the request beside it, and the sample-size line follows the same resolution, since lavaan’s case count is every case carrying data under FIML and the complete cases under listwise.

The scaling radio rewrites every marker; a widget edit rewrites none. Factor scaling is a model-wide identification choice, so the radio replaces the whole F ~~ 1*F set – every latent, second-order factors included – and passes std.lv = TRUE beside it. The markers are otherwise user state rather than a projection of the radio: the buffer may pin a subset of the latents, or a second-order factor alone, and regenerating the list on every widget write would rewrite the model’s identification on an unrelated edit – so the list is only pruned to latents that still exist, and a deleted factor leaves no dangling marker.

Allow factors to correlate pins exogenous pairs only, and never an interaction latent. lavaan auto-correlates the exogenous latents alone – an endogenous factor’s disturbance is uncorrelated by default – so an orthogonal model needs 0* on the exogenous pairs and nothing else, and a 0* typed outside that set is kept as written. An interaction latent is exempt from the block: its covariances with its own components are what the product-indicator construction estimates rather than assumes, and pinning them to zero would bias the interaction path itself.

A structural model takes the group dimension without the cascade. A group variable on a measurement-only buffer runs the invariance cascade, because its levels – loadings, intercepts, residuals – are measurement constraints; once a ~ line is present the model is fitted once across the groups with whatever group.equal levels are ticked, so the Constrain across groups panel appears on that routing alone – the cascade picks its own levels rung by rung and would ignore it. The panel offers six of lavaan’s levels, reaching regressions, the one a structural comparison needs and the cascade never gets to; a group.partial release is the cascade’s alone. The groups are passed as group.label in their canonical level order, since lavaan otherwise orders them by first appearance and row order would decide the reference group. The group-size line counts complete cases under listwise deletion, since raw sizes would hold the adequacy rules against an n lavaan never sees; its free-parameter bar is the fit’s count divided by the number of groups – this group’s share, and an equated rung frees fewer, so the bar is the conservative side – and the refusal below two groups counts the column’s distinct non-blank values. With the visual editor disabled there is no spec to pick columns from, so the buffer goes to sem() text-driven and a selected group variable is dropped with a warning.

BCa is attempted only above the replication bar, and the cascade ends at percentile. A BCa interval needs more successful replications than observations, and lavaan raises rather than falling back when it has fewer – so boot.ci.type = "bca" is tried only when the replication count exceeds the analysed N, then bca.simple, the bias-corrected percentile interval, then perc, the universal floor; if that fails too the error reaches the card. The CI header names the type that was used, and each standardized solution retries independently, so a standardized column can sit on a percentile interval while the unstandardized one is BCa – the headers say which. The replication count is the shared bootstrap setting, used as it stands rather than overridden per module; the card instead reports what the chosen count costs – below 1000 replications, and below the N that BCa needs – and the fit runs under the bootstrap seed setting, so a rerun reproduces its bounds.

A robust estimator is bootstrapped as its base estimator with its own test, and the fit indices come from a companion fit. lavaan refuses every robust estimator under se = "bootstrap" – MLR, MLM, MLMV and MLMVS want ML, WLSMV DWLS, ULSMV ULS – so under Bootstrap each is spelled as its base estimator plus the robust test it names: MLR as ML with yuan.bentler.mplus, MLM, MLMV and MLMVS as ML with satorra.bentler, scaled.shifted and mean.var.adjusted and expected information, WLSMV and ULSMV as DWLS and ULS with scaled.shifted. That leaves every χ² test and scaled index as the named estimator reports it and trades only its sandwich standard errors for the bootstrap’s, which the card says beside the estimator. It does not leave the categorical robust indices alone: under a bootstrap lavaan forces observed information, which those indices invert, so they drift in the third or fourth decimal. Every fit index is therefore read from a second fit of the model without the bootstrap, under the estimator the card names, and the bootstrap fit supplies only the standard errors and intervals. The header names the estimator attempted even when the fit fails. The invariance cascade’s latent-mean level takes the same spelling, where a robust estimator used to fail and fall back to analytic standard errors without a word.

Validation. On the café’s six satisfaction items typed ordinal, the bootstrapped WLSMV fit’s indices equal a plain WLSMV fit’s to every digit, where the bootstrap fit alone would have printed TLI 1.0117 for 1.0112. On lavaan’s HolzingerSwineford1939 every ML-family spelling, MLR under FIML included, reproduced its estimator’s fit measures exactly.

Check data #

The battery checks what the fit will see. The diagnostics run on the columns the model names – indicators and the observed variables of the structural equations – as soon as two of them are numerical, so the sample-size, normality and outlier rows describe the fit’s own variables rather than the selection; every selected numerical variable is the fallback while there is no model yet. A non-numerical column the model names is excluded and named in a toast because R’s per-variable loop would die on it and the readback would return a truncated battery that reads as a clean run. Sample-size adequacy is graded per free parameter when the count is known – the 5–10 cases per free parameter guideline is what the literature states, and a model always estimates more parameters than it has variables, so cases per variable is only its upper bound – and the count is the same one the summary reports, unknown for a grouped model until lavaan’s lands. The factorability battery (KMO, Bartlett, determinant, rank) is skipped for a model-based run with no latents, since a path model has nothing to factor, and runs on mixed correlations when a checked column is typed Ordinal, the estimator the fit will use for it.

Validation rules #

Run is closed on what lavaan would die on, and warned on what it would fit badly. The rules mirror what the fit cannot survive: a factor with one indicator, a second-order factor with one first-order factor, a name lavaan resolves to nothing – the editor already underlines it, but without the gate the fit runs and dies inside lavaan – and a latent interaction whose product columns nothing can build, since the =~ line survives the measurement model outgrowing it and the run would otherwise reach lavaan naming columns that do not exist. A syntax error closes the run through its own flag rather than a rule, because the buffer is what runs and a failed parse leaves no model to validate. The advisory rules are Bollen’s (1989): a factor carries its own identification only with three indicators, or with two plus a free link to another latent – a regression, a second-order loading, a free covariance, or the automatic correlation among exogenous latents while Allow factors to correlate is on – so a two-indicator factor with no link warns rather than blocks, since a user constraint can still identify it, and a just-identified model warns that its fit indices may not be meaningful.

Free indicator intercepts and free latent means are jointly identified only across groups. In a single group the observed means are the only mean-level data, and observed intercepts together with latent means are one parameter set more than those moments can pin down, so the pair is refused whenever the model has latents – a pure path model has no latent means and the rule has nothing to bite on. The invariance cascade identifies its mean structure through the cross-group equality constraints instead, and pins meanstructure = TRUE, int.ov.free = TRUE, int.lv.free = FALSE for every rung rather than forwarding the widgets – freeing latent means in every group, or fixing the indicator intercepts, would change the df of every level – which is why the pair is a warning rather than a refusal on a measurement-only model with a group variable: the boxes have no effect on that path.

The df gate prefers lavaan’s count and abstains on a grouped model. Free parameters and df are read from lavaan’s own parameter table on every parse, the only source that agrees with what sem() estimates – == constraints, repeated labels and the mean structure included – and dataPoints follows from it, since lavaan’s df already nets out the moments fixed.x removes. A widget-side estimate mirrors lavaan’s defaults for the gap before the first parse lands: one variance per non-exogenous observed variable, one per latent under marker scaling, auto.fix.first on each factor’s first loading, auto-covariances among exogenous latents, one mean per modelled variable under a mean structure, a repeated label constraining to its first occurrence and := lines derived. The estimate is single-group only: nothing in it has a group axis, and answering anyway would put a grouped model’s single-group df in front of the gate, so it abstains and lavaan’s group-aware count fills in when the parse lands.

Reading results #

The χ² row’s df is the test’s; the summary’s is the model’s. Under MLMVS the mean- and variance-adjusted statistic is referred to a Satterthwaite-adjusted df – fractional, lavaan’s df.scaled – which is the df the scaled test was evaluated against and belongs beside it in the Details cell; the model’s df – how many claims the model makes, the count a report quotes and the one the N:parameter ratio rests on – is the summary’s line. Printing one number for both would misstate one of them, so the two are kept apart on the card.

CFI, TLI, RMSEA and SRMR read the fit-cutoff settings; the other graded indices read a fixed scale. GFI, AGFI, NFI, RFI and IFI were never calibrated on the simulations the CFI and TLI cutoffs come from (Hu & Bentler, 1999) – GFI in particular moves with sample and model size – so tying them to the CFI setting would let a deliberate change to one index quietly move five others; they are banded at 0.95 / 0.90 / 0.85 for excellent, acceptable and poor fit, the conventional values, and the setting leaves them alone. PNFI and PGFI have no accepted bands and the information criteria rank models rather than grade one, so none of those takes a verdict. A value is rounded to its displayed precision before it is banded, so a cell printed as 0.95 is never graded as falling short of 0.95.

A saturated model is not graded, and RMSEA is withheld below 50 df. With df ≤ 0 the model reproduces the observed moments by construction, so every index sits at its optimum whatever the model is worth: the Interpretation column is dropped and a note says there is no fit left to assess – the numbers lavaan computed still print. RMSEA and its close-fit test are biased upward and unstable in models with few degrees of freedom, the more so at modest N (Kenny, Kaneko & McCoach, 2015), so with df < 50 and N not above 200 their verdict is left blank with a note saying why, the values still shown.

Significance on the diagram follows the standardization on screen, and an estimate lavaan could not standardize is drawn as no estimate. A standardized coefficient is tested against its own standard error, so the dashes and stars of a standardized view read the standardized p rather than the unstandardized one – otherwise a path could be dashed in the picture and starred in the table beside it. Where standardizedSolution returns no value – a Heywood case, an indicator with no variance – the edge takes a neutral grey stroke of fixed width and no label, rather than borrowing the raw coefficient and colouring it on the standardized ramp, which would put an unstandardized number on a scale of correlations; its tooltip quotes the unstandardized value and says which one it is. A labelled free parameter shows its label in italics in place of the estimate, since a := line’s terms are read off the picture and the estimate is in the table.

Each standardization carries its own inference columns, and a grouped card shows the standardized estimates alone. Standardizing rescales a coefficient and its sampling distribution with it: lavaan’s standardizedSolution returns an SE, z, p and CI for each standardized estimate that are not the unstandardized ones rescaled, so the toolbar reveals a whole column group per standardization rather than a lone estimate column, and a standardized estimate is never read beside the unstandardized z – the two would mix two different tests. Under a grouped fit each group already spends six columns; the standardized SE, CI, z and p are dropped, since a per-group standardized Wald p tests no between-group difference – the score test in the invariance card does – and the two checkboxes become a one-at-a-time picker, so the column count stays readable.

A grouped card widens the per-parameter tables, repeats the matrices and keeps one fit table. A per-parameter family keeps one row per parameter and takes a column block per group, so the groups are read side by side; a := parameter is reported by lavaan once for the whole model, with no group of its own, so the defined-parameters table keeps its single-group shape rather than folding the value into the first group’s block. A matrix – discriminant validity, residual correlations – has no column dimension left to spend on groups, so it is repeated per group. The fit belongs to the model, not to a group, so the fit table stays single, with each group’s χ² contribution in the χ² row’s Details, since a per-group contribution says whether one group carries the misfit, and the group sizes stated above it, since that n is the one the adequacy rules are graded against. The diagram is drawn once per group on one shared colour ramp: per-panel ramps would give the same stroke width and colour a different meaning in each group, defeating the comparison the side-by-side layout exists for.

Defined parameters fall back to the delta method, and the decomposition’s interval is its own. Without bootstrap draws a := product’s SE is the delta method’s (Sobel), whose normal reference distribution a product of coefficients does not follow, so the card says so under the table and recommends bootstrapping – the Sobel test block has the argument. The effect decomposition computes its own interval from the route products – the percentile interval of the per-route bootstrap draws when bootstrapping is on, a delta-method interval otherwise – rather than reading lavaan’s, so its header must not inherit the estimate table’s BCa or BC claim, which this block never computed; its standardized columns are point estimates alone because standardizedSolution has no row for a quantity the model never named. The bootstrap SE divides by the replication count as lavaan’s does, so the two tables agree where they report the same quantity.

The decomposition enumerates every pair joined by an indirect route, and declines cycles and more than 500 routes. A pair whose only route is the direct path has nothing to decompose and gets no row; a pair with any indirect route reports the total over every route, the total indirect over the routes of two steps or more, one specific row per such route and the direct path when the model has it, in Mplus’s order. A non-recursive model has no finite set of routes – Kahn’s algorithm finds the cycle before any walk starts – so the card says so and points at :=. The route count grows exponentially with a dense model, and a card that takes minutes to draw is worse than one that says it gave up, so the walk stops at 500 routes and the card declines, the := route staying open for one specific effect out of a large model. A path fixed by the user contributes its pinned value to every route it lies on and no derivative.

The ~~ diagonal splits three ways. An endogenous latent’s variance parameter is its disturbance – what its predictors leave unexplained – not the variance of the factor; listed under Factor variances it overstates how much of that factor is unexplained, and in a second-order model every first-order factor would read as unexplained variability its higher-order factor accounts for. So the diagonal is bucketed by the variable’s role – exogenous latent, predicted latent, observed indicator – and each bucket is its own table. A negative estimate is a Heywood case in any bucket, flagged in bold red with a warning above the table; the Interpretation column that names it appears only when one exists, and travels inside each group’s block under a grouped fit, since a variance can go negative in one group and not another.

R² with ordinal indicators is a share of the underlying response’s variance. Under WLSMV’s delta parameterization the continuous response behind each ordered item has its variance fixed at 1, so the R² lavaan reports for an ordered indicator is the proportion of that latent variance the factor explains, not of the observed column’s; whenever the fit had ordered indicators (lavNames(fit, "ov.ord") is non-empty) the section heading and the diagram legend both say so, so the number is never read as the ordinary one. lavInspect(fit, "rsquare") returns an indicator’s communality and an endogenous variable’s structural R² in one list, read against different expectations, which is what the Type column is for.

Reliability is computed against each first-order factor as its own true score. Left to itself semTools::compRelSEM follows the model upstream, so an endogenous or second-order factor’s ω would be computed against the factor predicting it and collapse – 0.20 where the indicators load .70–.80 – so every first-order factor is named as its own true score, and ω is the reliability of that factor’s own indicators whatever predicts it. compRelSEM’s coefficient is congeneric ω, which under a congeneric model is also composite reliability, so one coefficient is reported once, under the ω header. A latent ~~ observed covariance row makes semTools return a collapsed ω and an NA AVE for the factor owning that indicator – the row’s presence does it, even fixed at zero – so both are withheld for that factor with a note; α is unaffected. A higher-order factor has no observed indicators of its own, so α, ω and AVE do not apply to it; ω hierarchical – the share of its indicators’ total observed variance the higher-order factor’s own variance carries, computed from the model matrices over the total loadings $\mathrm{\Lambda }{(I-B)}^{-1}$, which reproduces compRelSEM’s own number where both apply – is banded on Reise, Bonifay & Haviland (2013), a total score reflecting one general factor at 0.75 and above and too multidimensional to score that way below 0.50. AVE keeps its own two-band rule at 0.50 (Fornell & Larcker, 1981), half the indicator variance captured, where α and ω take the shared reliability bands.

A factor whose loadings disagree in sign has its minority reversed for reliability. A reverse-keyed item fitted as answered loads negatively, and summed as answered it cancels the other items: α falls below zero and ω collapses although the factor is measured well. So the reliability block reverses, per factor and per group, the indicators on the minority sign of the standardized loadings, and the negative ones on a tie, a two-indicator factor included. That is psych::omega()'s flip = TRUE convention, so the factor module and this card key the same items the same way. A factor whose loadings all share a sign is left as it is, since α and ω do not change when every indicator is reflected at once. The reversal gives exactly what a refit on the reverse-scored items would, without refitting. ω goes through compRelSEM with a −1 composite weight on each reversed indicator, which reflects both its common-variance numerator ${(\sum d_i{\lambda }_i)}^2\psi$ and its observed-variance denominator. α is computed from the same observed covariance matrix reflected as $DSD$, with $D$ the diagonal of ±1, because compRelSEM refuses α on unequal weights. ω hierarchical reflects its items by the signs of the higher-order factor’s standardized total loadings. The fit, the loadings table and AVE (squared loadings) are untouched, and HTMT is built from absolute correlations, so the reversal cannot move it. The note under the table names each reversed indicator with its loading, so a flip at a loading near zero, which points to a misspecified indicator rather than a reversed one, can be recognized as such.

Validation. On the café’s two-factor CFA with slow as given (standardized loading −0.59 on Service), Service reads α 0.670 and ω 0.673, equal to nine decimals to a fit on a copy with slow reverse-scored in the data; summed as answered, the same items would give α −0.52 and ω 0.21. Food, whose loadings agree, stays at α 0.711 and ω 0.714. Grouped by branch, both branches equal the reversed-data fit (Service α 0.596 and 0.693, ω 0.598 and 0.700). With taste and fresh reverse-scored in the data, a one-factor model of taste, fresh, value, friendly and quick loads the two against the other three and reverses them: α 0.678 and ω 0.664, as on the original items. A two-indicator factor value + fresh on the same copy is a tie, so the negatively loading fresh is reversed: α 0.550 and ω 0.652, as on the original items. A second-order factor over Food and Service reverses slow by its total loading (−0.45) and reads ω hierarchical 0.485 with slow as given or reverse-scored.

Modification indices are filtered at 3.84, bolded at 10, and their multiplicity is stated. 3.84 is χ²(1) at α = .05 – the drop in χ² one freed parameter must reach to be significant on its own test – and lavaan returns a candidate for every fixed parameter, so the R side keeps only those above it. 10 is Kline’s (2016) practical-significance rule, marked in bold rather than used as a second filter, and each bucket shows its 20 largest with the total stated when it was cut. Every fixed parameter was tested, so a single-test threshold is optimistic, and the note under the heading states how many were, once, rather than under each bucket. The buckets follow the operator and the model’s shape: a structural model keeps every ~~ pair in one bucket and adds ~ paths, which the measurement-only card has no use for, while the CFA-mode card splits residual from factor covariances so the buttonless “other” bucket never swallows the one modification an orthogonal model most wants; a =~ candidate whose right-hand side is a latent is a higher-order loading, not a cross-loading, and goes to the read-only bucket rather than being applied as one.

Residual correlations are flagged on |z| ≥ 1.96, and on |r| > 0.10 only when no z exists. A residual of 0.10 is a lot in a sample of 1000 and nothing in a sample of 60, so a fixed correlation cutoff flags noise in small samples and misses real misfit in large ones; the standardized residual scales the discrepancy by its own standard error, so the two-sided 5% flag means the same thing at any N. Some estimators return no standardized residuals, and there the matrix falls back to the conventional |r| > 0.10 rule, the legend saying so along with the caveat that it does not adjust for sample size.

Mediation walkthrough #

Bootstrap and FIML combine: each resample is refitted under FIML. lavaan’s ordinary bootstrap draws cases – incomplete ones included – and refits the model on every draw with the options of the original fit, so under missing = "ml" each replication evaluates the available-data likelihood on its own resample rather than on a completed dataset: the interval carries the uncertainty FIML accounts for, and nothing is imputed. The price is a FIML fit per replication, which is why the pair is allowed and the manual asks only for care with the replication count on a large model.

Model comparison #

Fit indices are ranked only on a common scale. A best-value highlight is a ranking, and two fits are on one scale only when they were fitted on the same observed columns, with the same estimator and missing-data handling as lavaan resolved them, the same mean-structure settings, the same ordinal set and the same number of cases. A model that drops an indicator is fitted on a different set of variables, and its AIC is smaller not because it fits better but because there is less to explain; different observed sets also give different n under listwise deletion, and a scaled statistic cannot be differenced against an unscaled one. The resolved pair is compared rather than the requested one, since two requests that both substitute to WLSMV with pairwise deletion produced fits on one scale. When any pair differs the highlight is withheld from the whole table, and when no pair is nested the closing note names what differs instead of recommending AIC and BIC, which are the entries least comparable across case sets. χ² is never ranked, since it falls mechanically as parameters are added and its smallest value always belongs to the least constrained model. The runner has the last word on a pair the nesting rule accepted: anova() is refused when the fits analysed different numbers of cases, carry equal degrees of freedom, or either did not converge.

The group variable is part of the common scale; the equated levels are read by the nesting rule. A grouped fit’s χ² sums each group’s discrepancy against that group’s own moments, so a fit grouped on one variable reproduces a different set of moments from a single-group fit or from one grouped on another variable – their χ² and information criteria are on different scales, and a Δχ² between them tests no restriction. Which Constrain across groups levels a fit equated is a restriction, not a scale, and two fits of one model that differ only in those levels are the textbook nested pair – a configural baseline against its constrained twin. A parameter table parsed from the syntax alone cannot see it, since the group.equal option never reaches the syntax, so a grouped model is parsed across its groups with its group.equal and group.partial applied: lavaan’s parameter table then carries each equated level as generated == rows between group-specific labels, which the module translates to group-qualified parameter keys, and the rule’s own subset, kept-restrictions and degrees-of-freedom tests decide the pair. Equality by shared label is read the same way in a single-group model, so a*x2 + a*x3 is seen as a restriction too. A single-group fit ignores the boxes, since the hidden panel keeps its checkbox state and those levels never reached lavaan.

Validation. On a café structural model grouped by branch, the configural fit contains the structural-paths fit, which contains the loadings-plus-paths one; a fit freeing one loading under equated loadings sits between the configural and the fully equated one; a loadings-plus-intercepts fit pairs with the configural but not with a structural-paths one, and a single-group fit pairs with none – each tested pair’s Δdf is the refits’ own.

Nesting is parameter identity, plus one re-parameterisation: merged factors. Two models are paired for a chi-square difference test when one’s parameter space is contained in the other’s, and the module reads that off the two parameter tables by name: the constrained model’s free parameters are a subset of the free model’s, it keeps every fixed value and every equality the free model imposes, and lavaan gives it strictly more degrees of freedom. Read literally, that rule refuses the most common CFA contrast of all – a one-factor model against a multi-factor model over the same items – because G =~ x1 + x2 and F1 =~ x1 are different names for what is, after the merge, the same loading. The models are nested: a k-factor CFA with its factor correlations fixed at 1 has exactly the implied covariance structure of the CFA that merges those factors, whichever scaling convention identifies them, so the likelihood-ratio test compares two genuine parameter spaces (Bollen, 1989, ch. 7). The module therefore reads one equivalence past the names: when both models are plain CFAs over the same indicator pool – every indicator on exactly one factor, no regressions, no fixed or equated loadings, no residual covariance in one model that the other lacks, no ~~ or intercept row naming a factor, no higher-order factor – and every factor of the free model lies inside exactly one factor of the constrained model, the constrained model is taken to be the free one with the correlations among the factors it merges fixed at 1, and the Note column says so. Merging three factors into one fixes two correlations and equates the cross-correlations the merged factors had with the rest; lavaan’s own degrees of freedom count both, so Δdf is right without the module counting anything. The other familiar twins are not recognised – a second-order factor against its correlated first-order factors imposes a proportionality constraint on the factor covariance matrix rather than a fixed value, and a bifactor model is not nested in its correlated-factors counterpart at all – so those pairs are refused unless written in the same parameter names, and a manual re-parameterisation stays the honest route: F1 ~~ 1*F2 fixes the covariance, which is the correlation only under Fixed variance scaling, so write the constrained model under that scaling or fix the variances too.

The Δχ² reads the main significance level; each model’s exact-fit p reads the assumption level. The difference test’s null is that the more constrained model fits as well as the free one, and which of the two to keep is the analyst’s hypothesis rather than a check on an assumption, so its p is judged where the results tables’ are. The p row of the fit table is each model’s own χ² test of exact fit, whose null is the reversed one, and it reads the assumption test significance level as the single-model fit table does. Under a robust estimator the difference is anova()'s, which dispatches to lavaan’s lavTestLRT and returns the Satorra–Bentler scaled difference test once the fits carry a scaled statistic, so it does not equal the difference of the two scaled χ² values printed above it.

ΔCFI and ΔRMSEA are highlighted together, on Chen’s band for the fits’ n. The Δχ² grows with the sample on its own, and the practical indices are the counterweight, both fits being already in hand. The rule is the invariance card’s, directional and corroborated: a pair is highlighted only when the constrained model is worse on both – ΔCFI at or below −0.010 and ΔRMSEA at or above 0.015 above 300 cases, −0.005 and 0.010 at or below (Chen, 2007) – since one index moving on its own is not evidence. ΔSRMR is not computed here, so the corroboration is RMSEA’s alone, and the highlight is a flag beside the test rather than a verdict. A tested pair has one n by construction, so the band is read off it.

The comparison refits without the bootstrap. The card shows fit indices and likelihood-ratio tests, which depend on the point estimates alone; resampling would only move standard errors and confidence intervals, neither of which the card displays, at the cost of a bootstrap per queued model.

Confirmatory factor analysis #

Data diagnostics #

KMO and Bartlett’s test are withheld on a singular matrix. Both are read off the inverse of the correlation matrix. When solve() fails, psych’s KMO() falls back to exactly 0.5 for every variable – a sentinel indistinguishable from a real value without its flag – and Bartlett’s statistic, built on −log|R|, diverges as the determinant goes to zero, so it is guaranteed significant for the very reason the model cannot be fitted. Printing either would put a number that means nothing beside the rank that explains it, so both read Unreliable (singular matrix): KMO whenever the inverse failed or a pair is perfectly correlated, Bartlett’s test also when the matrix is not positive definite or a variable is reproduced exactly by others. The recommendations follow the same rule – no low-KMO or non-significant-Bartlett line on a matrix the rank checks have already disqualified.

Outlier distances are robust by default, and the classical fallback takes its exact null. A classical Mahalanobis distance is measured from the mean and covariance the outliers themselves inflate, so clustered outliers mask one another; the minimum covariance determinant (Rousseeuw, 1984), via MASS::cov.rob, estimates the centre and scatter from the tightest half of the cases instead, under the reproducibility seed. MCD needs a margin of cases over variables and a non-degenerate spread in every column, so it runs only with more than twice as many complete cases as variables, and the classical estimator stands in – named on the card with its reason – when that fails, when the fit errors or when its scatter is singular. The threshold changes with the estimator: the χ² 0.999 quantile is the asymptotic null of a robust distance, whose centre and scatter do not depend on the row measured, while a classical distance is bounded by (n − 1)²/n – below that quantile for every indicator count a model reaches – and is graded on its exact scaled-beta null (Wilks, 1963). The case table is capped at the 50 largest, since a robust run can flag a large share of the sample.

The univariate bars are |skewness| > 2 and |excess kurtosis| > 7, on G1 and G2. They are the values past which West, Finch & Curran (1995) found ML’s χ² and standard errors going wrong, and they are stated for the bias-adjusted sample estimators, so the diagnostics use the same G1 and G2 as the rest of the app rather than the population forms, and a column’s skewness here matches the descriptive statistics module’s. A constant column has no shape and is left to the variance check.

Near-zero variance follows caret’s rule. A column is flagged when it is constant, or when its most frequent value outnumbers the second by 19 to 1 while fewer than 10% of its values are distinct – caret::nearZeroVar’s defaults (Kuhn & Johnson, 2013). Both conditions together pick out the column that is nearly constant rather than merely skewed or coarse, and the table prints both numbers so the flag can be checked.

Missing data warns on the grid and on any single variable. The listwise complete-case count the sample-size rows grade is driven by the worst variable, not by the average, so a grid that is 2% blank can still lose a fifth of its cases to one indicator. The recommendation therefore fires when the grid is more than 5% blank or any one variable more than 20%, and the per-variable table sets the second case in bold.

Single-group results #

The χ² and close-fit p cells read the assumption test significance level. Every SEM text reports the model χ² at .05 and reads it beside the fit indices rather than as a verdict, which was the case for leaving its cell at the main level. It is judged at the assumption test significance level all the same, because its null – exact fit – is the reversed one: a user who lowers the main level to 0.01 for the paths would otherwise see an unhighlighted p = 0.03 whisper “fits” in the misleading direction, and an exception for one table would need documenting as carefully as the rule it breaks. The RMSEA close-fit p (null: RMSEA ≤ 0.05) is the same kind of test and reads the same level; the p row of the model comparison card is the same χ² test per model and follows.

Measurement invariance testing #

The Δχ² of each invariance step reads the assumption test significance level. Its null is that the constraint holds. The literature’s preference for ΔCFI over the χ² difference test (Cheung & Rensvold, 2002; Chen, 2007) is about the criterion, not the α: a Δχ² kept in the verdict at all is a check on an assumption, and it is judged where the other checks are – the structural level, which the Δχ² decides alone, moves with it. Under a robust estimator the difference is anova()'s, which dispatches to lavaan’s lavTestLRT and returns the Satorra–Bentler scaled difference test once the fits carry a scaled statistic, so it does not equal the difference of the two scaled χ² values the table prints.

A level fails only when the Δχ² test and the practical criteria agree, and the practical criteria are directional and corroborated. The χ² difference test grows sensitive to trivial differences as the sample grows, and ΔCFI, ΔRMSEA and ΔSRMR do not, so neither is allowed to fail a level alone: both rejecting is Fail, one of them is Mixed. The practical criteria are Chen’s (2007): ΔCFI is the primary one and fails only when corroborated by ΔRMSEA or ΔSRMR, because ΔSRMR is noisy at the intercept and residual levels and flagging on any one of the three would let it veto a genuinely invariant model on its own. Only a worsening counts – a constrained model that happens to fit better is no evidence against invariance. Chen states one band for a total N above 300 – ΔCFI −0.010, ΔRMSEA 0.015, ΔSRMR 0.030 for loadings and 0.010 for intercepts, thresholds and residuals – and a stricter one at or below it – −0.005, 0.010, 0.025 and 0.005 – and asks for the stricter band also when group sizes are markedly unequal, which the module reads as the largest group more than twice the smallest, however large the total. The card’s note states the band in force and the N it was selected on. Cheung & Rensvold (2002) are cited beside the ΔCFI cutoff they proposed, not beside the χ² test.

The structural level is judged on its Δχ² alone, and nominates nothing to free. Chen’s cutoffs cover loadings and the intercept, threshold and residual family; no published ΔCFI, ΔRMSEA or ΔSRMR cutoff covers latent variances and covariances, and inheriting the residuals’ band would pass off an arbitrary number as a criterion. So the level takes no band, its Δ values are printed for reference without a verdict of their own, and it can pass or fail but never be Mixed. It offers no partial-invariance candidates either: a latent variance that differs across groups is a finding about the constructs to report, not a nuisance parameter to free on the way to a higher level.

A level that did not converge is a failed level. A fit that runs to its iteration cap is an ordinary lavaan object carrying a full set of fit indices, so every consumer downstream – the Δ comparison, the score tests, the latent means – would read it as though it had settled. The runner refuses it where it is fitted: the level is named in a message, excluded from the comparisons and the cascade stops there, as it does for a level lavaan could not fit at all.

Only the latent-mean level is bootstrapped. Nothing else the cascade reports depends on a standard error – χ², the Δχ² test, CFI, RMSEA, SRMR and the score tests are point estimates, and lavTestScore falls back to the ordinary score test on a non-analytic covariance matrix anyway. Bootstrapping every level would also hand lavaan a rank-deficient covariance matrix whenever the replication count is below the number of free parameters, which it reports as a model that is not identified – at ordinary replication counts every level came back with that warning. The latent-mean table is the one place a standard error is read, from one level, so that level alone is refitted with bootstrap standard errors, under the bootstrap seed.

Ordinal indicators take semTools’ cascade, under Wu & Estabrook’s identification. lavaan fixes the intercepts of ordinal indicators at 0, so a hand-built cascade that equates intercepts constrains nothing and only frees the latent means: its “scalar” model has fewer degrees of freedom than the metric model it is supposed to nest in, and the χ² difference test fails outright. The levels are therefore generated by semTools::measEq.syntax() under the Wu & Estabrook (2016) identification and the theta parameterization, which nests them properly, identified by fixed latent variances whatever scaling the model uses, since marker scaling leaves the latent means unidentified once the intercepts are fixed to scale the items. A binary item has one threshold, identified only together with its loading, so a threshold level is inserted only when some item has more than two categories, and an all-binary item set equates thresholds and loadings in one step, named for both and score-tested on both kinds; its scalar level is its own, and its strict and structural levels are dropped. The generator regenerates every parameter, so a := definition, an == constraint or a fixed or labelled modifier in the model reaches no cascade fit, and each such line is listed on the card rather than lost silently.

An all-binary cascade reads its latent means from a scalar level that equates the residual variances too. With binary items, equal thresholds and loadings alone are equivalent to the configural model (Wu & Estabrook, 2016, Proposition 6); the metric level measEq.syntax generates also holds the residual variances at 1 in every group under the theta parameterization, which makes it their Proposition 12 model – testable, and with the latent-response intercepts free in the groups after the first but every group’s latent means fixed at 0. Equating the intercepts alone, the obvious next step, is their Proposition 11 model, which frees the residual variances again: it has fewer degrees of freedom than the metric level, is not nested in it, and was weakly identified in a check. The scalar level therefore equates thresholds, loadings, intercepts and residual variances together – nested in metric with $(G-1)(J-m)$ more degrees of freedom for G groups, J items and m factors – and frees the latent means in every group but the reference. It is reported as equating the intercepts, since the residual variances were already equal at the metric level. The strict and structural levels stay dropped.

Validation. On simulated binary data – six items on two factors, two groups, a 0.5 SD shift on the first factor – the scalar level fitted at df 24 against metric’s 20, Δχ² 3.1 on 4 df (p .54), and read the shift as 0.621 (SE 0.127); over 60 replications at n 300 per group it averaged 0.477 on the shifted factor and 0.002 on the other, with 95% coverage and the scalar level rejected in 8% at .05.

Higher-order loadings stay free across groups. lavaan applies group.equal = "loadings" to latent-on-latent =~ rows as well, so a second-order model’s metric level would equate the higher-order loadings and test more than measurement invariance, which is defined over the indicator-to-factor loadings. The continuous cascade releases them through group.partial at every level and names them on the card; the marker loading stays fixed, or the higher-order factor would have no scale in the groups after the first.

Candidates are ranked by the score test on the constrained level, not by the raw difference. A per-group estimate difference is scale-bound – 0.1 is large for a standardized loading and trivial for an intercept on a 1–7 item – and read off the configural fit it absorbs any latent mean gap into the intercepts and cannot reach a marker loading at all, so it nominates invariant parameters. lavTestScore on the level’s own constrained fit asks how much the fit would improve if one equality were released, which is comparable across loadings, intercepts, thresholds and residual variances; a parameter equated across more than two groups carries several constraints and is tested on releasing them jointly. The per-group columns beside it are context, read from the level below, where the parameter was still free – read off the constrained level they would be identical, and on the ordinal cascade the identification constants, 0 in every intercept cell and 1 in every residual-variance cell.

The ordinal scalar and strict levels are refitted in a form the score test can read. lavTestScore tests only the == constraints of a fit, and measEq.syntax does not equate the latent-response intercepts at the scalar level or the residual variances at the strict level with ==: under Wu & Estabrook’s identification it fixes them at their reference values, 0 and 1, in every group. Left as generated, those levels offer the score test nothing to release. So in every group after the first each such fixed row is made free and pinned to its former value by an == constraint – the same model, with the same parameter count, χ², degrees of freedom and fit indices, but one whose constraints the score test can release, a parameter spanning several groups on as many degrees of freedom as the refit that frees it. modindices is no substitute, since it ignores equality constraints, and nor is lavTestScore(add = …), which merges the added row into the fixed one the model already holds and returns an empty test. Only the level’s own class is rewritten: the binary scalar level keeps its residual variances at the metric level’s 1, and a continuous item’s intercepts, which lavaan equates with == already, are left as they are.

Validation. On simulated ordinal data – 2 factors × 3 items, n 300 per group, one item’s intercept shifted in the groups after the first – the shifted item ranked first in 60 of 60 replications across binary and 4-category sets with 2 and 3 groups, and its score test came within 0.42 of the unscaled Δχ² of the refit that frees it; null sets gave 2 of 180 item tests below .05. On one 3-group, 4-category set the score test read 13.92 against the refit’s 13.88.

The expected change is computed for nominated candidates only. A score test says whether releasing a constraint would improve the fit, not by how much the parameter differs, which depends on sample size; the expected parameter change (Saris, Satorra & Sörbom, 1987) gives the size in the parameter’s own units, and its effect on the latent means (Oberski, 2014) says whether the difference matters for the comparison the cascade exists for. Each candidate’s expected change needs a lavTestScore call of its own, a few seconds on a three-group model, so only the candidates whose p is below the significance level pay for it – with more than two groups the joint test of a candidate is that call already. The latent-mean column is shown where the level’s fit estimates the latent means, and only for a parameter that measures one factor.

Validation. On the same 3-group set the expected change of the shifted intercept was 0.431 and 0.430 against the refit’s 0.381 and 0.379, and the factor’s expected means 0.195 and 0.207 against the refit’s 0.201 and 0.214.

Candidates are freed one at a time. Every score test in a table is computed on the fully constrained level, so each asks what releasing its constraint would gain while every other constraint holds – including the constraint on a parameter that really does differ. One non-invariant parameter therefore raises its neighbours’ statistics: at the scalar level the latent mean absorbs part of one intercept’s shift, and the other intercepts on that factor then look shifted the other way. Freeing the top-ranked parameter and re-ranking removes that spillover; freeing everything flagged in one pass releases invariant parameters and weakens the comparison the cascade exists for. The card says so above the tables, and the Apply buttons still allow several at once for a user with a substantive reason.

Validation. On simulated ordinal data – 3 groups, 4 categories, n 300 per group, one intercept shifted in groups 2 and 3 – the shifted item’s factor-mates reached p < .05 in 33–73% of replications; freeing the shifted item alone and re-ranking dropped the largest of them from χ² 15.1 to 2.4.

Every factor keeps two invariant indicators, unless the level equates residual variances alone. A factor’s scale is comparable across groups only while its referent and at least one other indicator are invariant at that level; below that the latent means and variances are no longer on one metric, so an Apply that would leave a factor with fewer is refused. A residual variance carries none of the factor’s scale, so at a level that equates residual variances and nothing else – the strict level – any number of them may be freed.

The re-run refits the model the parameters were nominated on, and freed parameters accumulate. A freed parameter is a name in that model’s parameter table, and lavaan’s group.partial ignores, without warning, an entry naming a parameter the model no longer has – so a re-run fitted to whatever the editor holds by then could quietly test a different model with fewer releases than it claims. The re-run replays the card’s own syntax, options and group variable instead. A parameter freed at one level stays free at every level after it, since each of those nests in the level where it was released.

Latent means are compared on the canonical reference group and, by default, on the pooled SD. lavaan orders the groups by their first appearance in the data unless told otherwise, so the cascade passes the group variable’s levels in their own order and the reference group is the first level whatever the row order – re-sorting the file cannot flip the sign of every difference. The means are read from the scalar level, the first that identifies them, and from nowhere else: every earlier level fixes each group’s latent means at 0 – the all-binary metric level included – so a cascade whose scalar level did not fit reports no means rather than a table of zeros that reads as estimates. Every level up to it has to hold: equal intercepts on top of unequal loadings leave the scales incommensurate, so the card warns when one did not pass. The standardized difference divides by the latent SD pooled across groups and weighted by each group’s n − 1 (Hancock, 2001), which puts every group on one metric; each group’s own SD, lavaan’s std.all, stops being a common metric as soon as the latent variances are free, which they are at every level below the structural one – hence the picker, and the structural level’s verdict as the test of the pooling assumption.