Package 'mfrmr'

Description

mfrmr provides estimation, diagnostics, and reporting utilities for many-facet ordered-response measurement models: the Rasch-family RSM / PCM route and a GPCM extension in which one selected facet supplies level-specific discriminations. MML permits a different facet to supply category steps; JML requires the same facet for both roles. confint.mfrm_fit() with parm = "slopes" supplies approximate relative-slope intervals for eligible MML fits; model ranking and matched PCM/GPCM tests use separate checks in compare_mfrm(). gpcm_capability_matrix() explains the supported uses. An experimental two-family GPCM route also provides estimation and conditional scoring, with explicit uncertainty limitations. Native RSM/PCM conditional maximum likelihood (CML) conditions on Person totals without assuming an ability distribution. Its separate result class supplies structural contrast intervals through mfrm_facet_intervals(), without Person abilities, a population SD or ordinary fit diagnostics. See fit_mfrm() for the supported conditional estimation scope.

For assessment planning from observed scores, mfrm_g_study() and mfrm_d_study() separate sources of variation and project dependability for specified future rating plans. For rater feedback and assessment decisions, mfrm_report(), interrater_agreement_table() and classification_table() provide their respective reports, agreement measures and fixed-plan classification summaries. These analyses have different assumptions and uncertainty targets; follow the linked help for the question you want to answer.

Details

Start with the complete script in the Examples section:

1. Load the package and example ratings with load_mfrmr_data().

2. Fit a model with fit_mfrm().

3. Draw the Wright map with plot(fit).

4. Save results <- summary(fit) and inspect results$person_overview and results$facet_overview.

Each data row is one rating event. The person, facets, and score arguments name columns. The example uses marginal maximum likelihood (MML) and a rating-scale model (RSM) with shared category thresholds. head(ratings) displays the first six rows. Only Person, Rater, Criterion, and Score are used by this model; the example's Study and Group columns are additional labels. ⁠<-⁠ saves an object, and $ selects a named part: results$person_overview displays one table from the saved summary. The overview tables summarize distributions; as.data.frame(fit) returns individual person, rater, and criterion estimates. The assignment ratings <- load_mfrmr_data("example_operational") reads an existing synthetic dataset supplied with mfrmr. You do not need to create ratings first; the name on the left stores what the function returns.

Before interpreting or reporting estimates, read results$decision and follow its NextAction. The default summary does not compute diagnostics. For your own data, first inspect the design and score categories with describe_mfrm_data().

Optional packages

The basic example below uses packages installed with mfrmr. Its base-R plots do not need ggplot2. Additional packages support particular features; you do not need to install all suggested packages to get started.

• as_ggplot() needs ggplot2 for ggplot2 rendering.

• as_kable() needs knitr and as_flextable() needs flextable for their respective table formats.

• fit_mfrm_random_rater() needs RTMB version 2.0 or later for shared random-rater estimation.

• Corrected JML in fit_mfrm(), requested with jml_correction_order, needs nleqslv to solve the adjusted equations.

• ordinal_logistic_table() needs ordinal for ordinal-logistic DIF.

These are common examples; other optional analyses name their dependencies in their own help. A missing-package message identifies what to install. For example, use install.packages("ggplot2") if you want as_ggplot().

A saved mfrm_results() bundle can be reopened with base::readRDS() and passed to mfrm_report() without refitting or installing optional estimation packages. Reporting reuses stored results and retains unavailable diagnostics as unavailable. Updating an analysis or computing new diagnostics is a separate operation and may need the original dependencies. Save the complete bundle with base::saveRDS() and retain its matching package version/source; see "Updating saved analyses for 0.2.4" in mfrmr_workflow_methods before reusing older saved analyses.

Measurement scale

Native response models use the logistic metric (L-metric), with scaling constant $D = 1$; no extra 1.7 normal-ogive approximation factor is used. This is distinct from choosing a normal ability distribution. See the "Logistic metric (D = 1)" section of fit_mfrm() for logit interpretation with fixed or varying slopes and comparison with other software.

Where to go next

• mfrmr_methods_reference for the question-to-calculation map, plain-language explanations and links to the methodological literature.

• mfrmr_workflow_methods and vignette("mfrmr-workflow", package = "mfrmr") for your own CSV, column mapping, input troubleshooting, and the full analysis workflow. If vignettes are not installed, start with describe_mfrm_data() and recode_missing_codes().

• mfrmr_visual_diagnostics for choosing follow-up figures.

• mfrmr_reporting_and_apa and mfrmr_reports_and_tables for reporting.

• mfrmr_linking_and_dff for linking and differential facet functioning.

• classification_table() for fixed-plan future classification accuracy, consistency and conditional precision under ordinary RSM/PCM calibration.

• gpcm_capability_matrix for the GPCM extension.

• mfrmr_output_guide() for the broader purpose-to-function map.

A printable reference card is available at system.file("cheatsheet", "mfrmr-cheatsheet.pdf", package = "mfrmr").

Portable fixed calibration

A saved, versioned calibration artifact can be created from an eligible one-scale RSM or PCM MML fit under the fixed standard-normal scoring basis. Use mfrm_calibration_capabilities() before extraction, then follow mfrm_calibration_workflow to review, validate, freeze, save, load, and score the artifact. Review the returned batch with mfrm_calibration_score_methods before using its estimates. One-family GPCM MML also supports an estimated intercept-only normal population, shared or separate slope/step owners, unit weights and no anchors/interactions, subject to conditional source and per-batch numerical scoring checks. Estimated- population RSM/PCM and latent-regression MML remain fitted-object-only. RSM/PCM JML supports portable post-hoc EAP with a standard-normal reference prior, unit weights, no anchors/interactions and passing finite identified source checks. Shared-owner GPCM JML also supports conditional portable EAP after its joint-likelihood/curvature checks, preserving incomplete global audits. Explicit corrected GPCM JML instead uses adjusted-equation checks and preserves the correction order in format 5; residual calibration bias may remain. That reference prior is not estimated by JML. Experimental two-family GPCM MML uses format 6 with both ordered slope owners and a fixed N(0,1) prior, separate source/batch checks and no prior overrides. Artifact score uncertainty is conditional on the frozen point calibration and its recorded prior; loading validates consistency but does not authenticate an untrusted file.

Advanced scope

After the basic route above:

• the package supports latent-regression MML for ordered-response RSM / PCM models with a one-dimensional conditional-normal population model and explicit one-row-per-person covariates expanded through stats::model.matrix()

• GPCM support is summarized by gpcm_capability_matrix()

• GPCM supports the core fit/summary/scoring/information path, direct Wright/pathway/CCC plots, residual-PCA follow-up, and the residual-based diagnostics tables/plots as exploratory tools

• posterior-predictive checks and MCMC estimation are not available for GPCM; use external Bayesian software when they are required

• direct GPCM data generation through build_mfrm_sim_spec(), extract_mfrm_sim_spec(), and simulate_mfrm_data() is available when the specification carries both thresholds and slopes with the same owner; general simulation/design workflows do not support separate owners

• slope-aware fair_average_table() and estimate_bias() are available for GPCM with explicit caveats; build_apa_outputs(), build_visual_summaries(), run_qc_pipeline(), build_mfrm_manifest(), build_mfrm_replay_script(), and export_mfrm_bundle() are available as caveated partial reporting/export surfaces; score-side FACETS compatibility remains limited, and fully arbitrary-facet planning is outside the documented planning scope

• predict_mfrm_population() remains a scenario-level forecast helper and should not be described as the latent-regression estimator itself

• the current simulation/planning layer remains role-based for two non-person facets rather than fully arbitrary-facet planning, with boundaries exposed through planner metadata such as planning_scope, planning_constraints, and planning_schema

• latent-class mixture models and response-time / careless-rating adjustment are not estimated by mfrmr; use residual, person-fit, local-dependence, and rater-drift diagnostics as screening layers rather than as mixture-model substitutes

Equal weighting versus GPCM

The package's operational reference route is the Rasch-family RSM / PCM branch. That route enforces fixed discrimination and therefore preserves an equal-weighting scoring interpretation across observed ratings.

GPCM is supported because some users want a slope-aware model- comparison or sensitivity layer inside the same many-facet workflow. However, the package does not treat GPCM as a universal replacement for the Rasch-family route. A better fit under GPCM should be read as evidence about discrimination-based reweighting, not as an automatic reason to discard the equal-weighting model.

Observation weights are a different concept again. Optional Weight columns change how observed rating events enter estimation and summaries, but they do not create a free-form facet-weighting scheme and do not alter the fixed-discrimination meaning of RSM / PCM.

Public entry map:

• First-screen results: mfrm_results(), summary(res)$next_actions, and mfrmr_output_guide("entry")

• Interactive exploration: mfrm_results_interactive() only when prompts are explicitly wanted at the console

Function families:

• Model fitting: fit_mfrm(), summary.mfrm_fit(), plot.mfrm_fit()

• Legacy-compatible workflow wrapper: run_mfrm_facets(), mfrmRFacets()

• Diagnostics: diagnose_mfrm(), summary(diag), analyze_residual_pca(), plot(x, type = "residual_pca")

• Bias and interaction: estimate_bias(), estimate_all_bias(), summary(bias), bias_interaction_report(), plot(x, type = "bias")

• Differential functioning: analyze_dff(), analyze_dif(), mantel_haenszel_table() for separately matched binary responses, ordinal_logistic_table() for matched cumulative-logit comparisons, dif_interaction_table(), plot(x, type = "heatmap"), dif_report()

• Design simulation: build_mfrm_sim_spec(), extract_mfrm_sim_spec(), simulate_mfrm_data(), simulate_mfrm_recovery(), assess_mfrm_recovery(), evaluate_mfrm_design(), evaluate_mfrm_signal_detection(), predict_mfrm_population(), predict_mfrm_units(), sample_mfrm_plausible_values() (including fit-derived empirical / resampled / skeleton-based simulation specifications; fitted-object posterior scoring supports MML fits directly, latent-regression MML fits through the fitted population model when scored units also provide one-row-per-person background data, and JML fits through a post hoc reference-prior EAP layer; fit-derived simulation specifications also support direct GPCM data generation, recovery checks, role-based design evaluation, population forecasting, diagnostic-screening, and signal-detection helpers with documented caveats; curve reports, graph-only exports, fair-average tables, and bias screening are also available for GPCM with documented caveats)

• Reporting: build_apa_outputs(), build_visual_summaries(), reporting_checklist(), apa_table() for the full RSM / PCM route; GPCM uses these as caveated partial surfaces and retains direct table, plot, checklist, and summary-appendix routes

• Weighting review: compare_mfrm(), build_weighting_review(), build_model_choice_review(), compute_information(), plot(x)

• Case review: build_misfit_casebook(), plot(x, type = "unexpected"), plot(x, type = "displacement"), plot(x, type = "marginal_fit"), plot(x, type = "marginal_pairwise")

• Linking and scale maintenance: review_mfrm_anchors(), detect_anchor_drift(), build_equating_chain(), build_linking_review(), plot(x)

• Dashboards: facet_quality_dashboard(), plot(x, type = "facet_dashboard")

• Export / reproducibility: build_mfrm_manifest(), build_mfrm_replay_script(), build_conquest_overlap_bundle(), normalize_conquest_overlap_exports(), normalize_conquest_overlap_files(), normalize_conquest_overlap_tables(), review_conquest_overlap(), export_mfrm_bundle() for the diagnostics-compatible Rasch-family route; GPCM supports caveated partial manifest, replay, and bundle output as well as export_summary_appendix()

• Equivalence: analyze_facet_equivalence(), plot(x, type = "facet_equivalence")

• Data and anchors: describe_mfrm_data(), review_mfrm_anchors(), make_anchor_table(), load_mfrmr_data()

Data interface:

• Input analysis data is long format (one row per observed rating).

• Required columns are one person column, one ordered score column, and one or more non-person facet columns named in facets = c(...).

• Score values should be ordered integer categories. Binary 0/1 or 1/2 input is supported as the two-category Rasch-family special case; by contrast, fractional score values should be recoded before fitting rather than relying on automatic coercion.

• If keep_original = FALSE, unused intermediate categories are collapsed to a contiguous internal scale and the mapping is stored in fit$prep$score_map.

• If the intended scale has unused boundary categories, such as a 1-5 scale with only 2-5 observed, set ⁠rating_min = 1, rating_max = 5⁠ so the zero-count boundary category remains explicit in the data-support review. A missing boundary remains review evidence for the separate element- boundary contract. If unused intermediate categories should remain in the original scale, set keep_original = TRUE; a retained zero-count internal category in a polytomous fitted ladder creates an unsupported adjacent- step contrast, and fit_mfrm() returns a structured pre-optimization error.

• summary(describe_mfrm_data(...)) reports retained zero-count categories in Notes, printed Caveats, and ⁠$caveats⁠; summary(fit) carries full structured rows into printed Caveats and ⁠$caveats⁠, with ⁠Key warnings⁠ as a short triage subset. Summary-table exports route those rows through score_category_caveats or analysis_caveats.

• A category can be observed globally but unused by one rater or other facet level. Review that local pattern with data_quality_report() rather than treating it as a globally missing step or an automatic reason to select GPCM.

• Optional columns such as Subset, Weight, and Group support linking, weighted analysis, and fairness-focused follow-up workflows.

• Packaged synthetic data is available via load_mfrmr_data() or data(). Use list_mfrmr_data(details = TRUE) to distinguish the applied teaching, idealized, planted-effect, and larger sparse-design datasets.

Interpreting output

Core object classes are:

• mfrm_fit: fitted model parameters and metadata.

• mfrm_diagnostics: fit, facet-level reliability, and flag diagnostics, plus inter-rater agreement when one facet is treated as a rater facet.

• mfrm_bias: interaction bias estimates.

• mfrm_dff / mfrm_dif: differential-functioning contrasts and screening summaries.

• mfrm_population_prediction: scenario-level forecast summaries for one future design.

• mfrm_unit_prediction: posterior summaries for future or partially observed persons under the fitted scoring basis.

• mfrm_plausible_values: posterior draws for future or partially observed persons under the fitted scoring basis.

• mfrm_bundle families: summary/report bundles and draw-free plot data.

Typical workflow

1. Prepare long-format data.

2. Fit with fit_mfrm().

3. For RSM / PCM, diagnose with diagnose_mfrm() and prefer diagnostic_mode = "both" for final MML runs.

4. For RSM / PCM, run analyze_dff() or estimate_bias() when fairness or interaction questions matter; GPCM also supports estimate_bias() as a conditional screening review.

5. For RSM / PCM, report with build_apa_outputs() and build_visual_summaries().

6. For design planning, move to build_mfrm_sim_spec(), evaluate_mfrm_design(), mfrm_generalizability(), mfrm_d_study(), and predict_mfrm_population(). GPCM also supports direct simulation via extract_mfrm_sim_spec() / simulate_mfrm_data() and caveated role-based evaluate_mfrm_design() / predict_mfrm_population() routes. Those helpers assume two non-person facet roles even though the estimation core supports arbitrary facet counts. Treat evaluate_mfrm_design() as Monte Carlo design evaluation, and use mfrm_d_study() for analytic G/Phi design projections. Always read the IdentificationStatus, GStatus, and PhiStatus columns before reporting those projections; boundary or singular mixed-model fits are design-identification warnings, not high-stakes-ready reliability evidence. predict_mfrm_population() remains the scenario-level forecast helper, not the latent-regression estimator.

7. For future-unit scoring, retain an MML calibration when you want the fitted marginal model directly, use an active latent-regression MML fit when scored units also provide one-row-per-person background data, or use a JML calibration when a post hoc fitted-object EAP layer is acceptable; then score with predict_mfrm_units() or sample_mfrm_plausible_values().

8. For GPCM, use summary.mfrm_fit(), diagnose_mfrm(), analyze_residual_pca(), predict_mfrm_units(), sample_mfrm_plausible_values(), compute_information(), plot(x, type = "qc"), plot.mfrm_fit(), category_structure_report(), category_curves_report(), graph-only facets_output_file_bundle(), direct simulation-spec generation/data generation, recovery checks, fair_average_table(), estimate_bias(), and the residual-based table helpers with their documented caveats. Caveated APA/QC/export bundles and exploratory linking review plus role-based design evaluation, population forecasting, diagnostic-screening, and signal-detection helpers are available, while full score-side FACETS review, posterior-predictive checks, and MCMC estimation are not available for GPCM. Use gpcm_capability_matrix() as the formal boundary statement.

Model formulation

The Rasch-family branch used by RSM and PCM extends the basic Rasch model by incorporating multiple measurement facets into a single additive linear predictor on the log-odds scale.

RSM/PCM adjacent-category equation

For an observation where person $n$ with ability $\theta_n$ is rated by rater $j$ with severity $\delta_j$ on criterion $i$ with difficulty $\beta_i$, the probability of observing category $k$ (out of $K$ ordered categories) is:

$P(X_{nij} = k \mid \theta_n, \delta_j, \beta_i, \tau) = \frac{\exp\bigl[\sum_{s=1}^{k}(\theta_n - \delta_j - \beta_i - \tau_s)\bigr]} {\sum_{c=0}^{K}\exp\bigl[\sum_{s=1}^{c}(\theta_n - \delta_j - \beta_i - \tau_s)\bigr]}$

where $\tau_s$ are the Rasch-Andrich threshold (step) parameters in the RSM reference case and $\sum_{s=1}^{0}(\cdot) \equiv 0$ by convention. Additional facets enter as additive terms in the linear predictor $\eta = \theta_n - \delta_j - \beta_i - \ldots$.

This additive predictor generalises to any number of facets; the facets argument to fit_mfrm() accepts an arbitrary-length character vector.

Rating Scale Model (RSM)

Under the RSM (Andrich, 1978), all levels of the step facet share a single set of threshold parameters $\tau_1, \ldots, \tau_K$.

Partial Credit Model (PCM)

Under the PCM (Masters, 1982), each level of the designated step_facet has its own threshold vector on the package's common observed score scale. In the current implementation, threshold locations may vary by step-facet level, but the fitted score range is defined by one global category set taken from the observed data.

Generalized Partial Credit Model (GPCM)

With one slope family under GPCM (Muraki, 1992), the adjacent-category partial-credit kernel is multiplied by a positive slope $\alpha_g$ for the designated slope-facet level $g$:

$\ln\frac{P(X_{nij} = k)}{P(X_{nij} = k-1)} = \alpha_g(\theta_n - \delta_j - \beta_i - \tau_{hk}).$

MML allows a slope owner $g$ distinct from the step owner $h$; JML requires slope_facet == step_facet. Slopes are identified on the log scale with geometric mean 1. This makes GPCM a slope-aware sensitivity/extension route, not a replacement for the equal-weighting RSM/PCM interpretation. A separate provisional two-family MML route fits the multiplicative Uto–Ueno equation with explicit fixed-N(0,1) identification, exactly two facets and second-owner steps. The first family's slopes have geometric mean one and the second family's slopes are free. Fixed-grid EM and adaptive direct MML supply summaries and conditional fitted curves without intervals. Both support separately checked experimental component log-Wald and location/contrast intervals and descriptive posterior residuals; component profiles preserve the fitted fixed-grid or adaptive integration and require separate numerical checks. predict_mfrm_units() and the portable calibration workflow supply experimental conditional new-Person EAP and posterior intervals, excluding calibration uncertainty. Coverage and population transport remain unqualified. The one-family diagnostic and inference routes below do not apply to it. See the Two slope families section of fit_mfrm(). Unit slopes reduce to the equal-discrimination PCM kernel. Under default one-family MML, an intercept-only person distribution $N(\beta_0,\sigma^2)$ is estimated. The geometric-mean-one slopes are relative discriminations and $\sigma\alpha_g$ are their equivalent fixed-latent-standard-deviation values, so the conventional common discrimination degree of freedom is retained. The legacy gpcm_mml_identification = "fixed_standard_normal" mode fixes both $\sigma=1$ and the slope geometric mean to one and is therefore a narrower relative-discrimination model. Under JML, geometric-mean-one is required to resolve the scale of jointly estimated person coordinates. The former "bounded" label described limited model/workflow scope, not finite optimizer box constraints. Use GPCM with explicit scope restrictions. The GPCM JML objective is unpenalized. Certified extreme-person or facet recession is reported through typed primary boundary results; a finite optimizer iterate is retained only as a numerical trace and is not a finite JML maximum.

Ordered-response scope

The implemented response-model scope is ordered categorical only. Binary responses are the $K = 1$ special case of the same formulation, so they are handled through the ordinary ordered-score interface. This means mfrmr supports ordered binary and ordered polytomous data under RSM and PCM, plus GPCM with one designated slope_facet and, under MML, a possibly different step_facet. Unordered nominal/multinomial response models are outside the documented model scope, as are Poisson, negative-binomial, and grouped binomial-trial count-response families. A positive observation weight weights the conditional ordered-rating contribution; it is not a general collapsed-person frequency table and does not change the response family or model dependence among replicated ratings. Under MML, powering conditional responses within one Person is not the same as replicating a complete Person pattern after marginalization. Consequently, integer count values supplied as scores are modeled as ordered category codes, not as counts from a Poisson or related distribution.

Estimation methods

Marginal Maximum Likelihood (MML)

MML integrates over the person ability distribution using Gauss-Hermite quadrature, in the broader marginal-likelihood framework introduced by Bock & Aitkin (1981) for IRT:

$L = \prod_{n} \int P(\mathbf{X}_n \mid \theta, \boldsymbol{\delta}) \, \phi(\theta) \, d\theta \approx \prod_{n} \sum_{q=1}^{Q} w_q \, P(\mathbf{X}_n \mid \theta_q, \boldsymbol{\delta})$

where $\phi(\theta)$ is the assumed normal prior and $(\theta_q, w_q)$ are quadrature nodes and weights. Person estimates are obtained post-hoc via Expected A Posteriori (EAP):

$\hat{\theta}_n^{\mathrm{EAP}} = \frac{\sum_q \theta_q \, w_q \, L(\mathbf{X}_n \mid \theta_q)} {\sum_q w_q \, L(\mathbf{X}_n \mid \theta_q)}$

MML does not estimate one fixed person parameter per respondent and therefore avoids that JML incidental-parameter mechanism. Its structural inference is nevertheless conditional on the specified population distribution, response model, quadrature approximation, and regularity conditions; it is not an automatic small-sample guarantee.

Note: Bock & Aitkin (1981) is the canonical citation for the Gauss-Hermite-quadrature MML framework. The default mfrmr engine (mml_engine = "direct") optimises this marginal log-likelihood by direct gradient methods (BFGS / L-BFGS-B), not by Bock & Aitkin's signature EM algorithm. The "em" and "hybrid" engines do follow the EM template but use the selected gradient-based M-step rather than B&A's probit IRLS, because the target is the polytomous Rasch family rather than B&A's 2PL probit model.

Conditional Maximum Likelihood (CML)

For supported RSM/PCM models, CML conditions on each Person's total and observed assignments, removing ability parameters without specifying their population distribution. A generalized elementary-symmetric-function recursion evaluates the conditional likelihood and information. Identification and finite-MLE checks precede structural contrast intervals; their Wald interpretation is asymptotic, not a general finite-sample coverage guarantee. The separate mfrm_cml_fit and comparison-only results supply no Person abilities or population SD. Use mfrm_facet_intervals() for their native-scale intervals. CML does not support GPCM or ordinary diagnostics.

Joint Maximum Likelihood (JML)

JML estimates all person and facet parameters simultaneously as fixed effects by maximising the joint log-likelihood $\ell(\boldsymbol{\theta}, \boldsymbol{\delta} \mid \mathbf{X})$ directly. It does not assume a parametric person distribution, which can be advantageous when the population shape is strongly non-normal. Its incidental-parameter concern is not simply a small-person-sample issue: structural-parameter bias can persist as the number of persons increases when the number of ratings or items per person remains fixed (Neyman & Scott, 1948). The package still accepts "JMLE" as a backward-compatible alias, but user-facing summaries and documentation use "JML" as the public label.

See fit_mfrm() for practical guidance on choosing an estimator. General GPCM inference uses MML; JML remains useful for descriptive point estimates and FACETS comparisons.

Strict marginal diagnostics

For RSM / PCM, diagnose_mfrm(..., diagnostic_mode = "both") returns two complementary targets: the legacy residual / EAP diagnostics and a marginal_fit layer whose expected counts and pairwise summaries are integrated over the posterior quadrature bundle rather than plugged in at the EAP point. The screen is structured as limited-information evidence (Orlando & Thissen, 2000; Haberman & Sinharay, 2013; Sinharay & Monroe, 2025), not as an omnibus accept / reject test, and it complements rather than replaces separation / reliability and inter-rater agreement summaries. The full derivation, with notation and pairwise local-dependence events, lives in vignette("mfrmr-mml-and-marginal-fit", package = "mfrmr").

Statistical background

Key statistics reported throughout the package:

Infit (Information-Weighted Mean Square)

Weighted average of squared standardized residuals, where weights are the model-based variance of each observation:

$\mathrm{Infit}_j = \frac{\sum_i Z_{ij}^2 \, \mathrm{Var}_i \, w_i} {\sum_i \mathrm{Var}_i \, w_i}$

Expected value is 1.0 under model fit. Values below 0.5 suggest overfit (Mead-style responses); values above 1.5 suggest underfit (noise or misfit). Infit is most sensitive to unexpected patterns among on-target observations (Wright & Masters, 1982).

Note: The 0.5–1.5 range is the general "productive for measurement" band given by Linacre (2002, RMT 16(2), 878). Context-specific bands come from Wright & Linacre (1994, RMT 8(3), 370): 0.8–1.2 for high-stakes MCQ, 0.7–1.3 for run-of-the-mill MCQ, 0.6–1.4 for rating-scale surveys, 0.5–1.7 for clinical observation, and 0.4–1.2 for judged performance. See also Bond & Fox (2015) for textbook summaries of these conventions.

Outfit (Unweighted Mean Square)

Simple average of squared standardized residuals:

$\mathrm{Outfit}_j = \frac{\sum_i Z_{ij}^2 \, w_i}{\sum_i w_i}$

Same expected value and flagging thresholds as Infit, but more sensitive to extreme off-target outliers (e.g., a high-ability person scoring the lowest category).

ZSTD (Standardized Fit Statistic)

Wilson-Hilferty (1931) cube-root transformation that converts the mean-square chi-square ratio to an approximate standard normal deviate:

$\mathrm{ZSTD} = \frac{\mathrm{MnSq}^{1/3} - (1 - 2/(9\,\mathit{df}))} {\sqrt{2/(9\,\mathit{df})}}$

Values near 0 indicate expected fit. The conventional $|\mathrm{ZSTD}| \ge 2$ and $|\mathrm{ZSTD}| \ge 3$ cutoffs are heuristic two- and three-standard-deviation reference bands (Wright & Linacre, 1994; see also Wilson & Hilferty, 1931), not calibrated 5\ tests. Parameter estimation, sparse cells, the selected df convention, and repeated screening across elements all affect that interpretation. Diagnostics default to FACETS fourth-moment df, admitting df > 0 and capping ZSTD at +/-9. The explicit engine convention withholds df < 1. whexact = FALSE uses the expression above (FACETS YES); TRUE removes its additive q/3 term (FACETS NO). Both use a cube root. The former linear TRUE transformation is available only as whexact = "legacy_linear".

Residual basis under MML vs JMLE engines

For method = "MML" fits, the standardized residuals behind Infit, Outfit, and ZSTD are evaluated at EAP person measures, which are shrunken toward the population mean. JMLE engines such as FACETS evaluate the same formulas at unshrunken JMLE estimates, so MnSq and ZSTD values are not numerically interchangeable across the two residual bases, most visibly for extreme-scoring persons. Use method = "JML" when an external FACETS fit comparison requires a JMLE-style residual basis, and see facets_fit_df_guide() for the separate standardization-side df/ZSTD conventions.

PTMEA (Point-Measure Correlation)

Pearson correlation between observed scores and estimated person measures within each facet level. A positive value is directionally consistent with the fitted orientation, but it is not a confirmatory test of alignment. A negative value prompts review of orientation, response patterns, and model specification; it does not by itself establish a scoring error.

Separation

Package-reported separation is the ratio of adjusted true standard deviation to root-mean-square measurement error:

$G = \frac{\mathrm{SD}_{\mathrm{adj}}}{\mathrm{RMSE}}$

where $\mathrm{SD}_{\mathrm{adj}} = \sqrt{\mathrm{ObservedVariance} - \mathrm{ErrorVariance}}$. Higher values indicate the facet discriminates more statistically distinct levels along the measured variable. In mfrmr, Separation is the model-based value and RealSeparation provides a more conservative companion based on RealSE.

Reliability

$R = \frac{G^2}{1 + G^2}$

Analogous to Cronbach's alpha or KR-20 for the reproducibility of element ordering. In mfrmr, Reliability is the model-based value and RealReliability gives the conservative companion based on RealSE. For MML, these are anchored to observed-information ModelSE estimates for non-person facets; JML keeps them as exploratory summaries.

The facet's role determines the interpretation. High Person separation reliability indicates distinguishable abilities. High rater separation reliability indicates distinguishable severities, not agreement between raters. When raters are intended to be interchangeable, reliable severity differences merit review rather than a high-reliability quality label. Low rater separation reliability can also result from imprecise estimates; it does not establish agreement or equivalence. The index neither tests every rater pair nor quantifies the practical size of severity differences. Its Person interpretation is analogous to test reliability, but it is not numerically identical to Cronbach's alpha.

For the person facet under MML, the same $G$ and $R$ formulas are applied to EAP person measures with posterior SDs in the error slot. EAP measures are shrunken, so their observed variance is already deflated (approximately the true variance times the reliability), and subtracting the mean posterior variance deflates it again. The reported MML person separation/reliability is therefore a conservative summary: it is systematically lower than the IRT empirical-reliability convention $\mathrm{Var}(\mathrm{EAP}) / (\mathrm{Var}(\mathrm{EAP}) + \overline{\mathrm{PSD}^2})$ and is not numerically comparable to JMLE-based person separation reliability from FACETS. The gap is small when measurement is precise and grows as precision drops. Person rows can still carry the model-based precision tier because posterior SDs are model-based quantities; that tier describes the SE source, not FACETS comparability. Use method = "JML" when a FACETS-style person separation table is required, and treat MML person rows as conservative summaries.

This is a Rasch/FACETS-style separation reliability on the fitted logit scale, not an intra-class correlation. Use compute_facet_icc() only when you want the complementary random-effects variance-share view on the observed-score scale; for non-person facets, large ICC values indicate systematic facet variance rather than desirable measurement reliability.

Strata

Number of statistically distinguishable groups of elements:

$H = \frac{4G + 1}{3}$

Three or more strata are commonly used as a practical target (Wright & Masters, 1982), but in this package the estimate inherits the same approximation limits as the separation index.

Key references

• Andrich, D. (1978). A rating formulation for ordered response categories. Psychometrika, 43, 561–573.

• Bond, T. G., & Fox, C. M. (2015). Applying the Rasch model (3rd ed.). Routledge.

• Bock, R. D., & Aitkin, M. (1981). Marginal maximum likelihood estimation of item parameters: Application of an EM algorithm. Psychometrika, 46, 443–459.

• Burnham, K. P., & Anderson, D. R. (2002). Model selection and multimodel inference: A practical information-theoretic approach (2nd ed.). Springer. (AIC / BIC weights and Delta-IC bands used by compare_mfrm().)

• Drasgow, F., Levine, M. V., & Williams, E. A. (1985). Appropriateness measurement with polychotomous item response models and standardized indices. British Journal of Mathematical and Statistical Psychology, 38(1), 67–86. (Source for the lz person-fit statistic implemented in compute_person_fit_indices().)

• Haberman, S. J., & Sinharay, S. (2013). Generalized residuals for general models for contingency tables with application to item response theory. Journal of the American Statistical Association, 108, 1435–1444.

• Eckes, T. (2005). Examining rater effects in TestDaF writing and speaking performance assessments: A many-facet Rasch analysis. Language Assessment Quarterly, 2, 197–221.

• Muraki, E. (1992). A generalized partial credit model: Application of an EM algorithm. Applied Psychological Measurement, 16(2), 159–176. (Source for the GPCM response-probability kernel used by the package's bounded fit_mfrm(model = "GPCM") route. The slope-aware fair_average_table() and estimate_bias() extensions are package-specific, caveated adaptations built from that kernel, not procedures proposed in Muraki 1992.)

• Muraki, E. (1993). Information functions of the generalized partial credit model. Applied Psychological Measurement, 17(4), 351–363. (Companion paper to Muraki 1992 that derives the GPCM item information identity $I_j(\theta) = D^2 a_j^2 \mathrm{Var}(T \mid \theta)$ via Samejima's (1974) polytomous information formula. This is the canonical reference for compute_information() under GPCM.)

• Samejima, F. (1974). Normal ogive model on the continuous response level in the multidimensional latent space. Psychometrika, 39, 111–121. (General polytomous information formula that Muraki 1993 specializes to the GPCM.)

• Snijders, T. A. B. (2001). Asymptotic null distribution of person fit statistics with estimated person parameter. Psychometrika, 66(3), 331–342. (Source for the lz_star correction in compute_person_fit_indices() when person estimates come from the JML/fixed-effect route. MML/EAP person scores are left uncorrected because EAP does not satisfy the Snijders estimating-equation setup.)

• Linacre, J. M. (1989). Many-facet Rasch measurement. MESA Press.

• Linacre, J. M. (2002). What do Infit and Outfit, mean-square and standardized mean? Rasch Measurement Transactions, 16(2), 878.

• Masters, G. N. (1982). A Rasch model for partial credit scoring. Psychometrika, 47, 149–174.

• Orlando, M., & Thissen, D. (2000). Likelihood-based item-fit indices for dichotomous item response theory models. Applied Psychological Measurement, 24, 50–64.

• Orlando, M., & Thissen, D. (2003). Further investigation of the performance of S-X2: An item fit index for use with dichotomous item response theory models. Applied Psychological Measurement, 27, 289–298.

• Sinharay, S., Johnson, M. S., & Stern, H. S. (2006). Posterior predictive assessment of item response theory models. Applied Psychological Measurement, 30, 298–321.

• Sinharay, S., & Monroe, S. (2025). Assessment of fit of item response theory models: A critical review of the status quo and some future directions. British Journal of Mathematical and Statistical Psychology, 78, 711–733.

• Wright, B. D., & Masters, G. N. (1982). Rating scale analysis. MESA Press.

• Wright, B. D., & Linacre, J. M. (1994). Reasonable mean-square fit values. Rasch Measurement Transactions, 8(3), 370.

• Wilson, E. B., & Hilferty, M. M. (1931). The distribution of chi-square. Proceedings of the National Academy of Sciences of the United States of America, 17(12), 684-688.

Model selection

RSM vs PCM

The Rating Scale Model (RSM; Andrich, 1978) assumes all levels of the step facet share identical threshold parameters. The Partial Credit Model (PCM; Masters, 1982) allows each level of the step_facet to have its own set of thresholds on the package's shared observed score scale. Use RSM when the rating rubric is identical across all items/criteria; use PCM when category boundaries are expected to vary by item or criterion. In the current implementation, PCM assumes one common observed score support across the fitted data, so it should not be described as a fully mixed-category model with arbitrary item-specific category counts.

MML, CML and JML

Marginal Maximum Likelihood (MML) integrates over the person ability distribution using Gauss-Hermite quadrature and does not directly estimate person parameters; person estimates are computed post-hoc via Expected A Posteriori (EAP). Joint Maximum Likelihood (JML) estimates all person and facet parameters simultaneously as fixed effects; "JMLE" remains a backward-compatible alias.

Estimator choice should follow the inferential target and assumptions. MML avoids estimating one fixed effect per person, but its structural inference depends on the specified population distribution and response model. JML does not assume a parametric person distribution and can be lighter computationally, but structural-parameter bias can persist as the number of persons increases when the number of ratings or items per person remains fixed. Sample size alone therefore does not determine the preferred route. For supported RSM/PCM structural contrasts, CML instead conditions out Person ability without a population model. It estimates neither abilities nor a population SD; same-fit population-SD standardization requires MML.

See fit_mfrm() for usage.

Fixed-calibration scoring after fitting

predict_mfrm_units() and sample_mfrm_plausible_values() score future or partially observed persons on a quadrature grid under the fitted scoring basis. For ordinary MML fits, these summaries inherit the fitted marginal calibration directly. For latent-regression MML fits, they use the fitted one-dimensional conditional normal population model and therefore require one-row-per-person background data for the scored units when the fitted population model includes covariates. Intercept-only latent-regression fits (population_formula = ~ 1) can reconstruct that minimal person table from the scored person IDs. For JML fits, mfrmr uses the fitted facet and step parameters together with a standard normal reference prior introduced only for the post hoc scoring layer. This is useful for practical fixed-scale scoring, but it should still be described as a limited approximation rather than as full ConQuest-style population modeling.

ConQuest overlap

The package supports a scoped latent-regression MML model, but the overlap with ConQuest should still be described conservatively. The documented comparison is binary RSM / PCM, one latent dimension, exactly one item facet, complete rectangular responses, and exactly one numeric person covariate beyond the intercept in a conditional-normal person population model. This is a scoped file-based comparison, not a claim of broad ConQuest numerical equivalence for polytomous, multidimensional, incomplete, or arbitrary imported designs.

Author(s)

Maintainer: Ryuya Komuro [email protected] (ORCID) [copyright holder]

Authors:

• Ryuya Komuro [email protected] (ORCID) [copyright holder]

See Also

Useful links:

• https://ryuya-dot-com.github.io/mfrmr/

• https://github.com/Ryuya-dot-com/mfrmr

• Report bugs at https://github.com/Ryuya-dot-com/mfrmr/issues

Examples

# Load the package
library(mfrmr)

# Load example ratings and look at the first six rows
ratings <- load_mfrmr_data("example_operational")
head(ratings)

# Fit the model
fit <- fit_mfrm(
  data = ratings,
  person = "Person",
  facets = c("Rater", "Criterion"),
  score = "Score",
  rating_min = 1, rating_max = 4, category_policy = "preserve",
  method = "MML",
  model = "RSM"
)

# Plot the results (Wright map)
plot(fit)

# Save the summary, then display its tables
results <- summary(fit)
results$person_overview # One row summarizing person ability estimates
results$facet_overview  # One row per facet: number of levels, mean, SD, range

# Check the interpretation status and recommended next step
results$decision

Description

Compares observed-minus-expected scores between groups, or describes linked subgroup facet estimates after refitting. Residual differences do not isolate differential functioning and are returned without tests or classifications.

analyze_dif() is retained for compatibility with earlier package versions. In many-facet workflows, prefer analyze_dff() as the primary entry point.

Usage

analyze_dff(
  fit,
  diagnostics,
  facet,
  group,
  data = NULL,
  focal = NULL,
  method = c("residual", "refit"),
  min_obs = 10,
  p_adjust = "holm",
  comparison_review = NULL
)

analyze_dif(...)

Arguments

Details

Differential facet functioning (DFF) occurs when the difficulty or severity of a facet element differs across subgroups of the population, after controlling for overall ability. In an MFRM context this generalises classical DIF (which applies to items) to any facet: raters, criteria, tasks, etc.

Differential functioning is a threat to measurement fairness: if Criterion 1 is harder for Group A than Group B at the same ability level, the measurement scale is no longer group-invariant.

Differences between group ability distributions are not themselves DFF. Residual screens inherit the fitted model's ability and population assumptions. With fixed-standard-normal RSM/PCM MML, subgroup refits retain that population assumption; linking anchors do not estimate subgroup ability distributions. A residual difference can therefore reflect an inadequately represented group difference as well as differential facet functioning. For MML, observation expectations are evaluated at each Person's EAP ability, rather than integrated over the conditional ability distribution. This nonlinear substitution can produce different mean residuals across groups even under a correctly specified no-DFF model. A zero mean-residual contrast and absence of differential functioning are therefore distinct null hypotheses; correcting the SE alone does not make them equivalent.

Two methods are available:

Residual method (method = "residual"): Uses the existing fitted model's observation-level residuals. For each facet-level $\times$ group cell, the observed and expected score sums are aggregated and a standardized residual is computed as:

$z = \frac{\sum (X_{obs} - E_{exp})}{\sqrt{\sum \mathrm{Var}}}$

Pairwise contrasts compare the mean observed-minus-expected score in Group1 minus that in Group2. Positive values mean higher residual scores in Group1, not greater rater leniency for that group. Differences use score units. StdResidual is a descriptive scaling by model response variance, not a t-statistic. No residual-contrast SE, confidence interval, p-value or positive/negative classification is provided. Compatibility columns SE, t, df, p_value and p_adjusted contain NA.

Refit method (method = "refit"): Subsets the data by group, refits the MFRM model within each subset, anchors all non-target facets back to the baseline calibration when possible, and compares the resulting facet-level estimates. When linking, convergence, model-based MML precision, finite positive SEs, disjoint Person rosters and cell-count screens all pass, a conditional plug-in Welch statistic is also shown:

$t = \frac{\hat{\delta}_1 - \hat{\delta}_2} {\sqrt{SE_1^2 + SE_2^2}}$

This provides group-specific point estimates on a common scale when linking anchors are available. However, the plug-in standard error conditions on the baseline anchors and omits baseline-anchor uncertainty and cross-refit covariance. Refit statistics therefore remain screening evidence and set FormalInferenceEligible, SupportsFormalInference, PrimaryReportingEligible, and ETS_Eligible to FALSE. The Welch degrees of freedom use a response-count approximation, not independently estimated degrees of freedom for model-based SEs. Disjoint Person IDs do not establish independence or valid matching.

N_Group1 / N_Group2 count retained responses for the particular target level in each group; N_Subgroup1 / N_Subgroup2 count all retained subgroup responses. Failed or skipped fits have unknown retained counts (NA), not zero; ⁠$group_fits[[group]]$N_Input⁠ retains the input count. An absent target in a successfully prepared subgroup has a known count of zero. All baseline target levels remain in the result. Failed/skipped fits and absent levels have missing estimates and contrasts, with ComparisonStatus and ComparisonReason; they never inherit the baseline estimate. Finite estimates from sparse cells or incomplete diagnostics are retained, while ConditionalScreenReason explains withheld screening uncertainty. NPersons_Group1 / NPersons_Group2 count Persons in each target cell, and NPersons_Overlap counts shared Persons across the entire two subgroup fits. Any overlap withholds the independent-subgroup screen.

ContrastDirection describes the location parameter only: larger_location_in_group1, smaller_location_in_group1, or equal_locations. Translating this into higher/lower expected scores requires the facet orientation and the remaining response parameters. Subgroup step ladders or slopes can differ; a location difference alone does not order their entire response curves. The subgroup fits replay the baseline response family, resolved score range, step/slope facets, weighting, repeated-rating independence declaration, optimizer, MML engine, and numerical controls; the non-target anchors are the intentional linking change. Direct anchors on the target facet and JML Person anchors from the baseline are preserved. A target fixed by constraints in either subgroup keeps its point estimates, with Fixed1 / Fixed2 identifying the fixed values; its contrast, SE and screening tests are unavailable. Equal supplied anchors do not demonstrate absence of differential functioning. Use residual DFF to describe response departures, or explicitly fit a baseline with different target constraints before requesting free-parameter comparisons. Saved refit results created before this preservation rule require rerunning analyze_dff() from the baseline fit; loading an old result does not repair its subgroup estimates. Refit analysis fails closed when the baseline uses a user-specified active latent-regression population model, facet interactions, or group-anchor constraints, because those structures do not yet have a complete subgroup-replay and linking contract. The automatically generated intercept-only population model used for default GPCM-MML scale identification is replayed and is not treated as a substantive covariate model.

When facet refers to an item-like facet (for example Criterion), this recovers the familiar DIF case. When facet refers to raters or prompts/tasks, the same machinery supports DRF/DPF-style analyses.

Refit contrasts are not assigned ETS A/B/C labels. Even when subgroup point estimates share a linked logit scale, a validated joint, bootstrap, or replicate covariance contract is required before formal refit inference.

Refit screening tail areas are adjusted using Holm's step-down procedure by default, jointly across the returned facet-level/group-pair rows in this call. Holm's family-wise error control requires valid unadjusted p-values; the approximate screening tail areas here have not been shown to meet that requirement. Adjustment therefore does not establish an error-rate guarantee or formal inference eligibility. Alternative adjustments can be specified via p_adjust; see stats::p.adjust() for their assumptions.

Value

An object of class mfrm_dff (with compatibility class mfrm_dif) with:

• dif_table: data.frame of differential-functioning contrasts. Refit rows include Fixed1 and Fixed2; fixed-target rows have ReportingUse = "not_applicable_fixed_target" and no contrast or test.

• cell_table: (residual method) per-cell detail table.

• summary: counts by method-appropriate screening classification.

• group_fits: (refit method) per-group facet estimates.

• gpcm_boundary: for GPCM fits, a capability-boundary table.

• config: list with facet, group, method, min_obs, p_adjust settings.

Response matching for residual analyses

data selects fitted response rows using Person, all fitted facets, the observed score and fitted weight. Source column names, identifier trimming, the saved original/internal category map and the saved missing-code policy are respected. A changed score does not inherit a residual from another response with the same Person and facets. Rows with missing groups, invalid responses or no fitted match are excluded with reasons in ⁠$cases⁠; unmatched complete rows also produce a message. Requesting more copies of a response than were fitted is an error. No requested match is an error.

Repeated identical response rows are matched by their multiplicities, with a canonical group order. They must have identical fitted expectations and variances. FitRow then identifies an exchangeable fitted copy for this calculation, not a recovered occasion/event identity. This permits descriptive repeated-rating summaries without first-row group assignment or replication of unselected rows. It does not supply cluster-adjusted inference.

⁠$cases⁠ contains InputRow (row of supplied data, or fit$prep$data when data = NULL), FitRow (row of fit$prep$data), GroupValue, Included and Reason. Preserve the data and fit alongside this audit. N and N_Group1/N_Group2 count retained response rows; NPersons and NPersons_Group1/NPersons_Group2 count distinct Persons. The same Person may occur in multiple groups. Neither count establishes independence. Aggregates remain unweighted row summaries even when fitting used weights. Subsetting data retains the original model's fitted expectations; it does not recalibrate the model on the subset.

Choosing a method

Use method = "residual" to describe group differences in model residuals without fitting subgroup models. This does not test differential functioning. Use method = "refit" when you specifically want group-specific parameter estimates and can tolerate extra computation. Agreement between methods is not guaranteed by a universal per-group sample-size threshold: stability and detection depend jointly on effect size, category support, response-pattern overlap, linking strength, slope heterogeneity, and estimator behavior. min_obs is only a cell-computability/sparsity guard; it is not evidence of power, parameter stability, or sample-size adequacy.

Comparing with MH or ordinal logistic DIF

mantel_haenszel_table() supplies a separate binary conditional comparison with explicit matching strata and one response per Person-target cell. ordinal_logistic_table() supplies a separate cumulative-logit comparison with uniform/nonuniform terms and a proportional-odds diagnostic. Neither is a method of this function. A comparison needs the same persons, target response, group orientation and defensible matching variable. With a common binary Rasch item set and invariant other items, the total including the studied item removes the ability distribution from the conditional item response; a rest score generally does not. This result does not justify pooling different rater/task plans at the same raw total.

Ordinary MH targets a common conditional odds ratio and can miss opposing effects across matching strata. Ordinal logistic models distinguish group and matching-by-group terms, but their cumulative proportional-odds assumption is not generally implied by RSM/PCM. Repeated ratings from one Person cannot be counted as independent persons for either ordinary method. Agreement or disagreement with these procedures does not by itself validate a DFF conclusion. See the linking-and-DFF vignette for the comparison scope and primary references.

Interpreting output

• ⁠$dif_table⁠: one row per facet-level x group-pair with contrast, a residual mean difference or linked subgroup difference. Residual rows have no SE, test or binary classification. Includes Method, N_Group1, N_Group2, EffectMetric, ClassificationSystem, ContrastBasis, SEBasis, StatisticLabel, ProbabilityMetric, DFBasis, ReportingUse, PrimaryReportingEligible, and sparse columns.

• ⁠$cell_table⁠: (residual method only) per-cell detail with N, NPersons, ObsScore, ExpScore, ObsExpAvg, StdResidual and an interpretation note.

• ⁠$cases⁠: (residual method only) requested-row matching and exclusion audit, also retained by summary().

• ⁠$summary⁠: available/unavailable residual comparisons or linked- screening, insufficient-linking, fixed-target and unavailable rows (method = "refit").

• ⁠$group_fits⁠: (refit method only) list of per-group facet estimates and subgroup linking diagnostics, with EstimateStatus / EstimateReason and cell, subgroup and input counts.

• ⁠$gpcm_boundary⁠: for GPCM fits, a capability-boundary table marking the DFF/DIF output as caveated screening evidence.

• ⁠$comparison_review⁠: when supplied, the separate prospective interaction question, weights, counts, assumptions and identification states. Preserved by summary(), dif_report() and both DFF plot views. Its model can differ from the fitted screening model; read its fit_arguments and notes.

GPCM boundary

For GPCM, DFF/DIF rows are available as slope-aware screening evidence over the fitted expected-score and residual scale. Keep residual-method contrasts and interaction cells in screening language. Refit contrasts require explicit subgroup linking and precision support for conditional screening, but remain in screening language.

Saved results

Residual results, summaries and reports saved by earlier versions must be recomputed with this function using the fitted model and original data. This does not require refitting the model. Previously exported tables and figures should also be regenerated. Refit results created before the subgroup availability/cell-count checks must instead be recomputed with method = "refit". This reruns the subgroup fits from the saved baseline fit and original data. Reading an old result cannot repair borrowed estimates, dropped target levels or incorrect counts; summary, report and plotting methods reject these earlier refit objects.

Typical workflow

1. Fit a model with fit_mfrm(). For RSM / PCM fairness review, prefer method = "MML".

2. Run diagnose_mfrm() and, for RSM / PCM, prefer diagnostic_mode = "both" so legacy and strict marginal screens remain visible together.

3. Run analyze_dff(fit, diagnostics, facet = "Criterion", group = "Gender", data = my_data).

4. Inspect ⁠$dif_table⁠ for contrasts and their reporting status, and ⁠$summary⁠ for method-specific counts. Residual rows have no flags.

5. Use dif_interaction_table() when you need cell-level diagnostics.

6. Use plot(x, type = "heatmap") or dif_report() for communication.

See Also

fit_mfrm(), estimate_bias(), compare_mfrm(), dif_interaction_table(), plot(x, type = "heatmap"), dif_report(), subset_connectivity_report(), mfrmr_linking_and_dff

Examples

# Long-running workflow; run explicitly when needed.
# example("analyze_dff", run.dontrun = TRUE)
## Not run: 
toy <- load_mfrmr_data("example_bias")

fit <- fit_mfrm(toy, "Person", c("Rater", "Criterion"), "Score",
                 method = "MML", model = "RSM", quad_points = 7, maxit = 30)
diag <- diagnose_mfrm(fit, residual_pca = "none", diagnostic_mode = "both")
dff <- analyze_dff(fit, diag, facet = "Rater", group = "Group", data = toy)
dff$summary
# Read the residual mean difference and the number of observations per group.
head(dff$dif_table[, c("Level", "Group1", "Group2", "Contrast",
                       "N_Group1", "N_Group2")])
# The residual contrast is an observed-minus-expected average contrast
# between groups. It is useful for screening, but it is not an ETS
# A/B/C logit-delta classification.
dff_refit <- analyze_dff(fit, diag, facet = "Rater", group = "Group",
                         data = toy, method = "refit")
unique(dff_refit$dif_table$ClassificationSystem)
# Look for: "descriptive". Linked refit contrasts remain screening-only
#   because their plug-in uncertainty omits anchor uncertainty and
#   cross-refit covariance.
sc <- subset_connectivity_report(fit, diagnostics = diag)
plot(sc, type = "design_matrix", draw = FALSE)
if ("ScaleLinkStatus" %in% names(dff_refit$dif_table)) {
  unique(dff_refit$dif_table$ScaleLinkStatus)
}
# Look for: "linked" in `ScaleLinkStatus` confirms the focal and
#   reference groups share enough common elements for a comparable
#   contrast; "demoted_*" rows lose linking under the refit branch
#   and should be read as exploratory.

## End(Not run)

Description

Analyze practical equivalence within a facet

Usage

analyze_facet_equivalence(
  fit,
  diagnostics = NULL,
  facet = NULL,
  equivalence_bound = 0.5,
  ci_level = 0.95,
  conf_level = NULL
)

Arguments

Details

Pair differences use the full constrained covariance from the MML observed information: $\mathrm{Var}(A-B) = \mathrm{Var}(A) + \mathrm{Var}(B) - 2\mathrm{Cov}(A,B)$. No model is refitted. JML, inference-ineligible fits, missing or regularized covariance, and singular contrast covariance stop with an error. Levels with fixed or unavailable contrasts are not silently dropped. Known anchors are treated as fixed; their uncertainty is excluded. Non-unit observation weights are inference-ineligible, including weights normalized to mean one. Older bundles must retain the current readiness contract as well as the covariance basis before they can be displayed. A numerically verified inverse of ill-conditioned information is retained with a warning, not silently treated as strong evidence. When applicable, cautions, information_review and table column InferenceCaution preserve that warning through summaries and plots. Numerical verification does not establish the accuracy of the normal or chi-square approximations.

The heterogeneity table uses a joint Wald chi-square test of equality of the facet levels. Non-significant heterogeneity is neither necessary nor sufficient for practical equivalence. FixedChiSq, FixedDF, and FixedProb retain their column names but use this joint contrast test. Separation and reliability remain descriptive summaries.

GrandMean is the equally weighted mean of the facet estimates. For each deviation from that mean, uncertainty includes the covariance with the estimated mean. ROPEPct is the mass of its normal confidence distribution inside the practical bound; it is descriptive, not a Bayesian posterior probability or a separate equivalence decision.

The former BIC/Bayes-factor heuristic is unavailable: a Wald statistic is not a fitted likelihood comparison. For compatibility, BF01 is NA and BF01Label states why no value is supplied.

Value

A named list with class mfrm_facet_equivalence.

What this analysis means

The analysis asks whether differences between facet levels fall within a prespecified practical bound under the fitted model. These are asymptotic normal-approximation tests; numerical eligibility does not establish finite-sample coverage or adequacy of the rating design.

What this analysis does not justify

A non-significant difference is not evidence of equivalence. Pairwise conclusions are unadjusted for multiplicity: selecting some positive pairs does not provide family-wise error control for that selected set. Estimated linking or anchor uncertainty, population transport, and model misspecification require separate evaluation.

Decision rule

Each pair uses two one-sided normal tests at alpha 0.05. Equivalent is true when both tests reject non-equivalence, equivalently when its 90% interval lies strictly inside the bound. Decision summarizes all pairs as "all_pairs_equivalent", "partial_pairwise_equivalence", or "no_pairwise_equivalence_established". No pair may be omitted from the all-pairs summary. Heterogeneity and ROPE summaries do not change this rule. The all-pairs conclusion is an intersection-union decision: if any pair is truly outside the bound, declaring all pairs equivalent requires falsely passing that pair's test. Its asymptotic type-I error is therefore at most alpha when each pair's approximation is valid. This does not provide simultaneous confidence intervals or error control for selected pairs.

Interpreting output

Start with summary$Decision and examine the corresponding pairwise differences, SEs and intervals. A negative result can reflect imprecision or a material difference. chi_square addresses exact equality, while rope and forest describe proximity to the facet mean.

How to read the main outputs

• summary: pairwise decision, covariance basis and multiplicity convention.

• pairwise: differences, covariance-aware SEs, 90% intervals and TOST tests.

• chi_square: joint Wald heterogeneity test and descriptive separation.

• rope / forest: Measure, marginal SE and CI_Lower/CI_Upper, plus Deviation, DeviationSE and DeviationCI_Lower/DeviationCI_Upper for the difference from the estimated facet mean. Fixed marks an anchored or constraint-fixed location: its conditional marginal SE is zero and its marginal interval is absent. Its deviation from an estimated mean, or contrast with an estimated level, can still have nonzero uncertainty. The forest plot displays these deviations, not the anchored locations. No uncertainty in supplied anchor values is propagated.

• CIReason, DeviationCIReason, and pairwise CI90Reason distinguish "available", "fixed", "numerical_resolution" (the requested width rounds away), and "outside_numeric_range". Unavailable endpoints are NA; estimates and SEs remain. ROPEStatus = "unavailable" means the deviation interval could not be displayed, not that it is inside the bound. The fixed-alpha TOST p-values do not depend on the display ci_level.

Recommended next step

Review numerical integration and the model's uncertainty assumptions before interpreting a borderline result. Sensitivity to another practical bound should be reported transparently, without selecting a bound to obtain a desired decision.

Typical workflow

1. Fit and review an MML model with fit_mfrm().

2. Prespecify the practical bound and run analyze_facet_equivalence().

3. Read summary and pairwise.

4. Use plot(x, type = "facet_equivalence") for descriptive grand-mean proximity.

Output

A bundle with summary, chi_square, pairwise, rope, forest, and settings. Older bundles without the current inference/covariance basis must be recomputed before using summary(), print(), or plotting.

References

Schuirmann, D. J. (1987). A comparison of the two one-sided tests procedure and the power approach for assessing the equivalence of average bioavailability. Journal of Pharmacokinetics and Biopharmaceutics, 15(6), 657-680.

See Also

facets_chisq_table(), fair_average_table(), plot(x, type = "facet_equivalence")

Examples

toy <- load_mfrmr_data("example_core")
fit <- fit_mfrm(toy, "Person", c("Rater", "Criterion"), "Score",
                method = "MML", quad_points = 31, maxit = 150)
eq <- analyze_facet_equivalence(fit, facet = "Rater")
eq$summary[, c("Facet", "Elements", "Decision", "MeanROPE")]
head(eq$pairwise[, c("ElementA", "ElementB", "Equivalent")])

Description

One-stop review that combines the nesting, cross-tabulation, ICC, and design-effect reports into a single object. Its summary feeds into reporting_checklist() after checking the reviewed data against the fit.

Usage

analyze_hierarchical_structure(
  data,
  facets = NULL,
  person = "Person",
  score = "Score",
  compute_icc = TRUE,
  ci_method = c("none", "boot"),
  ci_level = 0.95,
  ci_boot_reps = 1000L,
  ci_boot_seed = NULL,
  igraph_layout = TRUE,
  icc_ci_method = NULL,
  icc_ci_level = NULL,
  icc_ci_boot_reps = NULL,
  icc_ci_boot_seed = NULL,
  missing = c("error", "omit")
)

Arguments

Value

A list of class mfrm_hierarchical_structure with:

• nesting: output of detect_facet_nesting().

• crosstabs: list of pairwise observation-count data.frames (long format, suitable for heatmap plotting).

• icc: output of compute_facet_icc() when requested.

• design_effect: output of compute_facet_design_effect() when requested.

• connectivity: named list with bipartite-graph component summary when igraph is available.

• summary: one-row summary used by downstream reporting helpers.

• facets: character vector of facet names that were reviewed (echoed for downstream reporting helpers that need to label rows by review scope).

• data_identity: reviewed column roles, input row count and a compact fingerprint for reuse checks. The optional person column also contributes to nesting pair counts; summary$NFacets counts the requested facets.

plot(x, type = "icc") draws variance shares and any available intervals without an embedded title or footnote. Save its returned object as p and use plot_data(p, "notes") to read the interval method or the reason bounds are unavailable. These notes survive saveRDS()/readRDS(); they do not change the underlying ICC estimates.

Interpreting output

• nesting: a detect_facet_nesting() object with every facet pair classified as Crossed / Partially / Near-perfectly / Fully nested.

• crosstabs: list of ⁠(LevelA, LevelB, N)⁠ long-format tables, one per facet pair. Plot via plot(x, type = "crosstab", pair = "FacetA__FacetB").

• icc: per-facet variance shares. See compute_facet_icc() for the two-scale interpretation.

• design_effect: per-facet design-effect approximations and descriptive equivalent row counts; not precision estimates for the full design.

• connectivity: number of bipartite components linking Person x facet levels. A single component is required for a common measurement scale; multiple components indicate a disconnected design.

Typical workflow

1. Optional: fit the MFRM with fit_mfrm().

2. Call analyze_hierarchical_structure(fit) (or on the raw data).

3. Read summary(x) for the condensed view.

4. Pass it as hierarchical_structure to reporting_checklist() and save the review with saveRDS() when its full tables are needed later. build_mfrm_manifest() records a separate fit-level rating-count review; it does not accept or archive this structural-review object.

Saved-result reuse

A reporting checklist checks the reviewed column roles and retained rating rows, including scores when an ICC table was computed. It accepts a subset of fitted facets and reports that scope explicitly. Reviews from raw data also work when the selected columns and rows match, using either canonical Person/Score names or the fit's original source-column names. Row order, unused factor levels, unused columns and MFRM estimates do not affect this match. Nesting and row counts are unweighted; when ICC is absent, the match does not depend on scores. These rules do not establish that a partial review covers all facets or that the design supports valid inference. A review of additional rows, other strata or columns absent from the fit remains usable on its own, but cannot mark this fit's checklist item as available. Recreate older reviews without data_identity from the saved fit or matching retained rows and original options; no MFRM refit is needed. See mfrmr_compatibility_layer.

References

McEwen, M. R. (2018). The effects of incomplete rating designs on results from many-facets-Rasch model analyses (Doctoral thesis, Brigham Young University). https://scholarsarchive.byu.edu/etd/6689/

Linacre, J. M. (2026). A User's Guide to FACETS, Version 4.5.0. Winsteps.com.

Kish, L. (1965). Survey Sampling. New York: Wiley.

Koo, T. K., & Li, M. Y. (2016). A guideline of selecting and reporting intraclass correlation coefficients for reliability research. Journal of Chiropractic Medicine, 15(2), 155-163.

See Also

detect_facet_nesting(), facet_small_sample_review(), compute_facet_icc(), compute_facet_design_effect(), reporting_checklist(), build_mfrm_manifest(), fit_mfrm().

Examples

toy <- load_mfrmr_data("example_core")
hs <- analyze_hierarchical_structure(toy,
                                     facets = c("Rater", "Criterion"),
                                     compute_icc = FALSE,
                                     igraph_layout = FALSE)
summary(hs)

# Full review when lme4 and igraph are available.
if (requireNamespace("lme4", quietly = TRUE) &&
    requireNamespace("igraph", quietly = TRUE)) {
  hs_full <- analyze_hierarchical_structure(toy,
                                            facets = c("Rater", "Criterion"))
  summary(hs_full)
  plot(hs_full, type = "icc")
}

Description

Legacy-compatible residual diagnostics can be inspected in two ways:

1. overall residual PCA on the person x combined-facet matrix

2. facet-specific residual PCA on person x facet-level matrices

Usage

analyze_residual_pca(
  diagnostics,
  mode = c("overall", "facet", "both"),
  facets = NULL,
  pca_max_factors = 10L,
  parallel = FALSE,
  parallel_reps = 200L,
  parallel_quantile = 0.95,
  parallel_method = c("residual_permutation", "model_bootstrap"),
  seed = NULL
)

Arguments

Details

The function works on standardized residual structures derived from diagnose_mfrm(). When a fitted object from fit_mfrm() is supplied, diagnostics are computed internally.

Conceptually, this follows the Rasch residual-PCA tradition of examining structure in model residuals after the primary Rasch dimension has been extracted. In mfrmr, however, the implementation is an exploratory many-facet adaptation: it works on standardized residual matrices built as person x combined-facet or person x facet-level layouts, rather than reproducing FACETS/Winsteps residual-contrast tables one-to-one.

Residual PCA should therefore be reported as residual-structure evidence, not as a formal proof of unidimensionality. It also should not be described as DIMTEST or UNIDIM: those essential-unidimensionality tests require a separate item-response-layer definition that is not uniquely determined by a many-facet long data set. In applied MFRM reporting, residual PCA is best triangulated with global residual fit, element fit, and Q3-style local-dependence screens.

Correlations use persons with observed residuals for each pair. Every column and pair must have defined correlations, and the resulting matrix must be positive semidefinite (allowing numerical roundoff). Missing correlations are not set to zero, and invalid matrices are not smoothed. Unavailable analyses retain an explanation in errors. Older stored PCA calculations are recomputed from the supplied observations when this function is called.

Output tables use:

• Component: principal-component index (1, 2, ...)

• Eigenvalue: eigenvalue for each component

• Proportion: component variance proportion

• Cumulative: cumulative variance proportion

When parallel = TRUE, the variance tables additionally include reference summaries (the method is retained in ParallelMethod):

• ParallelMean: mean reference eigenvalue

• ParallelCutoff: parallel_quantile cutoff of reference eigenvalues

• ExcessOverParallelCutoff: observed eigenvalue minus the cutoff

• ExceedsParallelCutoff: whether the observed eigenvalue exceeds the selected reference cutoff

With parallel_method = "residual_permutation", the comparison conditions on the fitted residuals and missingness pattern; it does not simulate responses or refit the model. It does not account for fitted-parameter uncertainty and is not a calibrated dimensionality test. If any requested permutation has unavailable or invalid correlations, the comparison is withheld rather than conditioning on successful permutations. parallel_status retains the successful count and the reason. More permutations improve numerical stability of the reference quantile, but do not establish error-rate control for the fitted model.

For mode = "facet" or "both", by_facet_table additionally includes a Facet column.

summary(pca) is supported through summary(). plot(pca) dispatches on class mfrm_residual_pca. Use type = "scree", "parallel_scree", "parallel_excess" or "loadings", with mode = "overall" or ⁠mode = "facet", facet = "Criterion"⁠. Existing prefixed spellings such as "facet_scree" remain supported. Plotting a saved PCA does not generate or refit a parallel reference. Save the complete result with base::saveRDS(), including errors, warnings, parallel_status and any bootstrap_trials; export the desired tables separately with utils::write.csv(). as_ggplot() converts a saved plot payload without redoing PCA. A blank unavailable plot is not evidence of an absence of residual structure.

Value

A named list with:

• mode: resolved mode used for computation

• facet_names: facets analyzed

• overall: overall PCA bundle (or NULL)

• by_facet: named list of facet PCA bundles

• overall_table: variance table for overall PCA

• by_facet_table: stacked variance table across facets

• parallel_settings, parallel_overall_table, parallel_by_facet_table, and parallel_status: returned for every call; the parallel tables are populated when parallel = TRUE

• bootstrap_trials, bootstrap_settings: attempted replicate/scope records and generation settings for the model bootstrap; otherwise NULL.

• errors: explanations for unavailable overall or facet PCA analyses; unavailable analyses have empty variance tables.

• warnings: named list of non-fatal PCA warnings captured from the underlying PCA engine. These indicate exploratory boundary conditions, not confirmatory evidence.

• InferenceTier, SupportsFormalInference, PrimaryReportingEligible, ReportingUse, and DecisionUse: machine-readable guards that keep this route exploratory and prohibit an automatic dimensionality or subscore decision.

Interpreting output

Use overall_table first:

• early components with noticeably larger eigenvalues or proportions suggest residual structure that may deserve follow-up. Small early components do not establish unidimensionality; inspect the residual aggregation, available person overlaps, fit and local-dependence screens.

Then inspect by_facet_table:

• describes patterns within each facet-level residual matrix. Different numbers of columns, aggregation and overlaps prevent treating eigenvalues across facets as a ranking of their contributions to a single decomposition.

Finally, inspect loadings via plot(x, type = "residual_pca") to identify which variables/elements drive each component.

Model-generated reference

parallel_method = "model_bootstrap" is an exploratory parametric-bootstrap reference for native RSM/PCM MML fits with a fixed standard-normal population, fixed quadrature, additive facets and unit weights. Anchors, dummy/positive facets, shrinkage, estimated population models, GPCM, JML and extended testlet/shared-rater models are not supported by this reference yet.

Each replicate draws one independent standard-normal ability per Person and shares it across that Person's retained rating rows. Fitted facet and step parameters generate conditionally independent ratings. The model is refitted with the original identification, category map, quadrature and optimization settings. Standardized residuals use the same scoring and aggregation as the observed analysis. One refit supplies every requested PCA scope.

The reference conditions on the analyzed assignment and observation pattern: it does not fill unassigned or missing ratings or simulate informative missingness. Sparse designs remain usable only when the required residual correlations are defined and their matrix is positive semidefinite. Source and refitted models must pass numerical convergence checks. The observed PCA itself can be unavailable despite successful fitting, for example because the assignment lacks shared-Person pairs. This is reported separately from refit failure. Every planned replicate is retained in bootstrap_trials, including warnings and failures. An affected scope has no reference cutoff if any replicate fails; successful replicates are not substituted or selected to form its reference. Raw eigenvalue draws, with missing rows for failures, are stored under each PCA bundle's parallel$draws.

Refitting includes variation from parameter estimation under the fitted null; generating parameters themselves are plug-in estimates, not posterior draws. Cutoffs are componentwise reference quantiles, not multiplicity-adjusted tests. Scanning every component does not preserve the nominal false-flag probability of one component. Absence of a flag is not evidence that all model assumptions hold. Few replicates give unstable tail quantiles. Neither a fixed number of replicates nor this algorithm establishes calibrated false-flag rates or proves unidimensionality. Inspect convergence, quadrature sensitivity, overlap and competing dependence explanations before substantive use.

Dimensionality and network follow-up

Checking a one-dimensional model does not require fitting a multidimensional alternative first. Residual structure can also reflect testlet, rater or assignment effects; neither the component count nor the count plus one is an estimated number of substantive abilities. Combine the screen with q3_statistic(), marginal-fit review and substantive task/criterion evidence. The current permutation reference does not replace a fitted-model parametric bootstrap. No native residual-network/EGA or calibrated dimensionality decision is supplied here. Residual communities would describe remaining associations, not automatically latent dimensions. See vignette("mfrmr-visual-diagnostics") for alternatives and their scope.

Typical workflow

1. Fit model and run diagnose_mfrm() with residual_pca = "none" or "both".

2. Call analyze_residual_pca(..., mode = "both").

3. Review summary(pca), then plot scree/loadings.

4. Cross-check with fit/misfit diagnostics before conclusions.

References

The residual-PCA idea follows the Rasch residual-structure literature, especially Linacre's discussions of principal components of Rasch residuals. The current mfrmr implementation should be interpreted as an exploratory extension for many-facet workflows rather than as a direct reproduction of a single FACETS/Winsteps output table.

The optional parallel analysis follows Horn's data-driven eigenvalue comparison logic and later recommendations to compare observed eigenvalues with high quantiles of a reference distribution. Here that reference is generated by within-column permutation of standardized residuals by default. The optional model bootstrap instead regenerates ratings and refits the supported model. The cited factor-retention literature does not calibrate either many-facet residual comparison as a fitted-model test.

• Horn, J. L. (1965). A rationale and test for the number of factors in factor analysis. Psychometrika, 30, 179-185. doi:10.1007/BF02289447.

• Glorfeld, L. W. (1995). An improvement on Horn's parallel analysis methodology for selecting the correct number of factors to retain. Educational and Psychological Measurement, 55, 377-393. doi:10.1177/0013164495055003002.

• Hayton, J. C., Allen, D. G., & Scarpello, V. (2004). Factor retention decisions in exploratory factor analysis: A tutorial on parallel analysis. Organizational Research Methods, 7, 191-205.

• Timmerman, M. E., & Lorenzo-Seva, U. (2011). Dimensionality assessment of ordered polytomous items with parallel analysis. Psychological Methods, 16, 209-220.

• Linacre, J. M. (1998). Structure in Rasch residuals: Why principal components analysis (PCA)? Rasch Measurement Transactions, 12(2), 636. https://www.rasch.org/rmt/rmt122m.htm

• Linacre, J. M. (1998). Detecting multidimensionality: Which residual data-type works best? Journal of Outcome Measurement, 2(3), 266-283.

Applications and broader measurement discussions below provide context; they do not validate this implementation's residual aggregation, reference cutoffs or automatic dimensionality decisions:

• Eckes, T. (2005). Examining rater effects in TestDaF writing and speaking performance assessments: A many-facet Rasch analysis. Language Assessment Quarterly, 2(3), 197-221.

• Yamashita, T. (2024). An application of many-facet Rasch measurement to evaluate automated essay scoring: A case of ChatGPT-4.0. Research Methods in Applied Linguistics, 3(3), 100133.

• Uto, M. (2021). A multidimensional generalized many-facet Rasch model for rubric-based performance assessment. Behaviormetrika, 48(2), 425-457.

• Aryadoust, V., Ng, L. Y., & Sayama, H. (2021). A comprehensive review of Rasch measurement in language assessment: Recommendations and guidelines for research. Language Testing, 38(1), 6-40.

• Tseng, W.-T. (2016). Measuring English vocabulary size via computerized adaptive testing. Computers & Education, 97, 69-85.

See Also

diagnose_mfrm(), plot(x, type = "residual_pca"), mfrmr_visual_diagnostics

Examples

# Saved synthetic JML/RSM analysis of example_core.
# For your own data, fit a model first; see ?fit_mfrm.
example_fit <- readRDS(system.file("examples", "core-jml.rds", package = "mfrmr"))
fit <- example_fit$fit
diag <- diagnose_mfrm(fit, residual_pca = "both")
pca <- analyze_residual_pca(diag, mode = "both")
pca2 <- analyze_residual_pca(fit, mode = "both")
summary(pca)
p <- plot(pca, type = "residual_pca", mode = "overall", view = "scree", draw = FALSE)
p$data$plot
head(p$data)
pca_pa <- analyze_residual_pca(diag, mode = "overall", parallel = TRUE, parallel_reps = 10)
head(pca_pa$overall_table)
head(pca$overall_table)

Description

Re-estimates a fitted many-facet model on new data while holding selected facet parameters fixed at the values from a previous (baseline) calibration. This transfers baseline coordinates into a new fit through direct equality constraints. Common-scale interpretation is conditional on a defensible baseline and invariant cross-run element identity. For GPCM, treat this as direct exploratory anchor/drift support rather than as the package's formal linking-synthesis route.

Usage

anchor_to_baseline(
  new_data,
  baseline_fit,
  person,
  facets,
  score,
  anchor_facets = NULL,
  include_person = FALSE,
  weight = NULL,
  model = NULL,
  method = NULL,
  anchor_policy = "warn",
  mml_integration = NULL,
  centering_scope = NULL,
  ...
)

## S3 method for class 'mfrm_anchored_fit'
print(x, ...)

## S3 method for class 'mfrm_anchored_fit'
summary(object, ...)

## S3 method for class 'summary.mfrm_anchored_fit'
print(x, ...)

Arguments

Details

This function automates the baseline-anchored calibration workflow:

1. Extracts anchor values from the baseline fit using make_anchor_table().

2. Re-estimates the model on new_data with those anchors fixed via fit_mfrm(..., anchors = anchor_table).

3. Runs diagnose_mfrm() on the anchored fit.

4. Compares the fixed non-Person anchor values in the new fit with their supplied values. This checks constraint consistency, not invariance.

The helper refuses a baseline unless its current readiness record has InferenceReady = TRUE. Also confirm compatible model, score coding, orientation, and population conventions, and document that anchor labels denote the same elements with invariant meaning. The function checks object shape and receiving-data anchor compatibility, but does not establish those substantive conditions. Fixing parameters can identify coordinates; it does not create observed common-rating links.

Only model, method, the centering convention and MML integration mode inherit baseline settings automatically. Supply other settings explicitly through ..., including rubric bounds, category policy, facet orientation, PCM/GPCM facet roles and any population model. RSM/PCM MML without a population model uses a fixed standard-normal ability distribution; anchors do not automatically estimate a new cohort's mean or variance.

The table fixes locations only: thresholds and GPCM slopes are re-estimated. Use extract_mfrm_calibration() and score_mfrm_calibration() when the question is scoring against a complete fixed calibration. Partial anchors retain the centering restrictions described in fit_mfrm(); this is not an unrestricted fixed-common-parameter linking implementation for arbitrary new levels. Elements present in the anchor table but absent from the new data are handled according to anchor_policy: "warn" (default) emits a message, "error" stops execution, and "silent" ignores silently.

The returned drift table contains only non-Person anchor rows that match the new data. Unanchored common levels and Person levels are absent. The compared values are fixed by construction, so zero differences do not show that raters or tasks remained stable. For separate cross-wave screening, fit the waves independently and use detect_anchor_drift(), retaining its offset and uncertainty limitations.

Differences compare each retained fixed value to the supplied anchor:

$\Delta_e = \hat{\delta}_{e,\text{new}} - \hat{\delta}_{e,\text{base}}$

An element is flagged when $|\Delta_e| > 0.5$ logits or $|\Delta_e / SE_{\Delta_e}| > 2.0$, where $SE_{\Delta_e} = \sqrt{SE_{\mathrm{base}}^2 + SE_{\mathrm{new}}^2}$. SE_Diff is a plug-in independence calculation. It does not propagate baseline-anchor uncertainty or cross-fit covariance, and anchored estimates in the new fit are constrained rather than independently re-estimated. Consequently, Drift_SE_Ratio and Flag are descriptive consistency screens, not formal z tests or calibrated error-rate decisions.

Value

Object of class mfrm_anchored_fit with components:

fit

The anchored mfrm_fit object.

diagnostics

Output of diagnose_mfrm() on the anchored fit.

baseline_anchors

Anchor table extracted from the baseline.

drift

Tibble of element-level drift statistics.

Which function should I use?

• Use anchor_to_baseline() when you have one new dataset and want to place it on baseline coordinates within the documented constraint system.

• Use detect_anchor_drift() when you already have multiple fitted waves and want to compare their stability.

• Use build_equating_chain() when you need cumulative offsets across an ordered series of waves.

Interpreting output

• ⁠$drift⁠: one row per retained non-Person anchor with columns Facet, Level, Baseline, New, Drift, SE_Baseline, SE_New, SE_Diff, Drift_SE_Ratio, and Flag. Read this as a constraint-consistency table. Zero differences are expected by construction and provide no independent evidence of stability.

• ⁠$fit⁠: the full anchored mfrm_fit object, usable with diagnose_mfrm(), measurable_summary_table(), etc.

• ⁠$diagnostics⁠: pre-computed diagnostics for the anchored calibration.

• ⁠$baseline_anchors⁠: the candidate table fed to fit_mfrm(). Some rows may be excluded when anchor_policy permits unmatched IDs. Inspect ⁠$fit$config$anchor_review$anchors⁠ for the retained direct constraints.

Typical workflow

1. Fit the baseline model: fit1 <- fit_mfrm(...).

2. Collect new data (e.g., a later administration).

3. Call res <- anchor_to_baseline(new_data, fit1, ...).

4. Inspect summary(res) for numerical readiness and fixed-value consistency; neither establishes anchor invariance or a successful substantive link.

5. For multi-wave drift monitoring, fit waves separately and pass the fits to detect_anchor_drift() or build_equating_chain().

See Also

fit_mfrm(), make_anchor_table(), detect_anchor_drift(), diagnose_mfrm(), build_equating_chain(), mfrmr_linking_and_dff

Examples

# Long-running workflow; run explicitly when needed.
# example("anchor_to_baseline", run.dontrun = TRUE)
## Not run: 
# Packaged synthetic teaching data; no data-generation code is needed.
ratings <- load_mfrmr_data("example_core")
people <- unique(ratings$Person)
# Two balanced 24-Person waves retain every Rater and Criterion.
d1 <- ratings[ratings$Person %in% people[1:24], , drop = FALSE]
d2 <- ratings[ratings$Person %in% people[25:48], , drop = FALSE]
fit1 <- fit_mfrm(d1, "Person", c("Rater", "Criterion"), "Score",
                 method = "MML", rating_min = 1, rating_max = 4,
                 category_policy = "preserve")
res <- anchor_to_baseline(d2, fit1, "Person",
                          c("Rater", "Criterion"), "Score",
                          anchor_facets = "Criterion", anchor_policy = "error",
                          rating_min = 1, rating_max = 4,
                          category_policy = "preserve")
summary(res)
# These fixed-value differences are not a test of Criterion stability.
head(res$drift[, c("Facet", "Level", "Drift", "Flag")])
res$baseline_anchors[1:3, ]

## End(Not run)

Description

Build an APA-oriented table handoff using base R structures

Usage

apa_table(
  x,
  which = NULL,
  diagnostics = NULL,
  digits = 2,
  caption = NULL,
  note = NULL,
  bias_results = NULL,
  context = list(),
  whexact = FALSE,
  branch = c("apa", "facets")
)

Arguments

Details

This helper avoids styling dependencies and returns a reproducible base data.frame plus manuscript-oriented metadata. It does not claim complete APA 7 or JARS compliance: outside the saved decisions table, digits applies the same rounding rule to every numeric column, so statistic-specific formatting (for example, exact p-value, confidence-interval, and effect-size conventions) and the target journal's final typography still require human review.

Supported which values:

• For mfrm_fit: "summary", "person", "facets", "steps"

• For summary() outputs or mfrm_summary_table_bundle: names listed in build_summary_table_bundle(x)$table_index

• For diagnostics list: "overall_fit", "measures", "fit", "reliability", "facets_chisq", "bias", "interactions", "interrater_summary", "interrater_pairs", "obs"

• For bias-result list: "table", "summary", "chi_sq"

• For ordinary score_mfrm_persons() output: "decisions". predict_mfrm_units() and score_mfrm_calibration() results and their summaries also accept this selector. Cuts must have been requested.

• For RSM/PCM fixed-facet intervals: "intervals" (default), "settings", "contrasts", and "clusters" when present. Method, confidence level and unavailable reasons remain with the selected estimates and bounds.

• For GPCM inference: "intervals" or "curves"; bootstrap results also retain "trials", "checks" and "source_checks" when recorded, "sampling", and "availability" for slope intervals or "test" for a null-model LRT. Extended results also expose "settings", and "clusters"/"contrasts" when present. Saved profile intervals also expose "profile", "profile_endpoints", "profile_checks" and "wald". Target and method columns are preserved.

Value

A list of class apa_table with fields:

• table (data.frame)

• which

• caption

• note

• digits

• branch, style

Interpreting output

Saved ordinary Person-score, unit-prediction and portable-score results with cut probabilities support which = "decisions". Numeric probabilities, cut boundaries and priors are retained without decimal rounding in the table. They condition on calibration and the scoring prior; the cut labels do not establish CEFR alignment or a decision rule.

• table: plain data.frame ready for export or further formatting.

• which: source component that produced the table.

• caption/note: manuscript-oriented metadata stored with the table.

Typical workflow

1. Build table object with apa_table(...).

2. Inspect quickly with summary(tbl).

3. Render base preview via plot(tbl, ...) or export tbl$table.

See Also

fit_mfrm(), diagnose_mfrm(), build_apa_outputs(), reporting_checklist(), mfrmr_reporting_and_apa, predict_mfrm_units()

Examples

library(mfrmr)

# Saved synthetic MML/RSM analysis of example_operational.
# For your own data, fit a model first; see ?fit_mfrm.
example_fit <- readRDS(system.file("examples", "operational-mml.rds", package = "mfrmr"))
fit <- example_fit$fit

# Turn one summary table into a table with a caption and note
results <- summary(fit)
tbl <- apa_table(results, which = "facet_overview",
                 caption = "Distribution of estimates within each facet")
tbl # Prints the table, caption, and note

# Extract the ordinary data frame for further formatting or export
tbl$table
# This table summarizes facets; use as.data.frame(fit) for individual estimates

Description

Post-hoc normal-model shrinkage helper that augments an mfrm_fit with descriptive empirical-Bayes adjustments for each non-person facet. The shrinkage variance $\hat{\tau}^2$ is estimated by method of moments from the facet-level point estimates and their standard errors:

$\hat{\tau}^2 = \max\!\left(0, \frac{1}{K}\sum_{j=1}^{K}\hat{\delta}_j^{2} - \overline{\mathrm{SE}^2}\right),$

where zero is the chosen shrinkage target and the moment is calculated on complete finite estimate/positive-SE pairs. A sum-to-zero identification constraint does not establish a known population mean or independent errors. The shrinkage factor is $B_j = \mathrm{SE}_j^2 / (\hat{\tau}^2 + \mathrm{SE}_j^2)$, with $\hat{\delta}_j^{EB} = (1-B_j)\hat{\delta}_j$ and $\mathrm{ShrunkSE}_j = \sqrt{(1-B_j)\mathrm{SE}_j^2}$.

Usage

apply_empirical_bayes_shrinkage(
  fit,
  facet_prior_sd = NULL,
  shrink_person = FALSE
)

Arguments

Details

ShrunkSE is a plug-in normal-model quantity conditional on the selected prior variance. It omits uncertainty in that variance, cross-level covariance and the fitted identification/anchor constraints. It is neither a calibrated SE nor a guaranteed lower bound. Unequal shrinkage can break the original sum-to-zero constraint. Original estimates, anchors and predictions are unchanged; this helper does not refit a hierarchical MFRM. The normal-normal formula is exact for an independent Gaussian observation with a known Gaussian prior variance. Applying it to jointly calibrated MFRM estimates is a descriptive approximation; a James-Stein improvement in risk is not established for this adjustment.

A zero estimated prior variance fully pools eligible estimates at zero and gives zero plug-in SE; this does not establish perfect precision. Optional plot whiskers are descriptive normal bands. Their legacy CI column names do not confer confidence-interval coverage; SupportsFormalInference is false. Invalid estimate/SE pairs retain their original values with unavailable shrinkage factors. With fewer than three valid pairs no shrinkage is applied. Unrepresentable moment variances also leave the original values unchanged, with unavailable factors and an explanatory Note in the report.

The report records complete-pair counts and mean shrinkage over those pairs. EffectiveDF sums retained weights over eligible levels, conditional on the chosen prior variance; it is not model degrees of freedom for testing or IC. Reapply this helper to existing fits to refresh old reports; no MFRM refit is needed. Regenerate previously saved plots and exports as well. Updated reports carry CalculationVersion = 2: positive prior variances no longer produce zero plug-in SE merely because a shrinkage factor rounded to one. Reapplication replaces the previous adjustment; shrink_person = FALSE removes any previous Person shrinkage columns. The original fitting inputs remain recorded separately from the latest post-fit shrinkage settings. build_mfrm_replay_script() replays the adjustment after fitting, preserving whether diagnostic SEs were attached by fit_mfrm(). Manually edited SEs or other table edits require their own reproducible editing steps.

Value

The same mfrm_fit, with augmented columns and a new shrinkage_report list entry, and with fit$config$facet_shrinkage set to "empirical_bayes".

Typical workflow

1. Fit the model as usual with fit_mfrm().

2. Call apply_empirical_bayes_shrinkage(fit) when small-N facets are present and zero-centered pooling is substantively defensible (see facet_small_sample_review()); small counts alone do not require it.

3. Report both the original and shrunk estimates in the manuscript, citing Efron & Morris (1973). build_apa_outputs() will add the sentence automatically when fit$config$facet_shrinkage is set.

References

Efron, B., & Morris, C. (1973). Combining possibly related estimation problems. Journal of the Royal Statistical Society: Series B, 35(3), 379-402.

Efron, B. (2021). Empirical Bayes: Concepts and methods (Technical report). Department of Statistics, Stanford University. https://efron.ckirby.su.domains/papers/2021EB-concepts-methods.pdf

Morris, C. N. (1983). Parametric empirical Bayes inference: Theory and applications. Journal of the American Statistical Association, 78(381), 47-55.

See Also

fit_mfrm() (which accepts facet_shrinkage directly), facet_small_sample_review(), compute_facet_icc().

Examples

# Saved synthetic JML/RSM analysis of example_core.
# For your own data, fit a model first; see ?fit_mfrm.
example_fit <- readRDS(system.file("examples", "core-jml.rds", package = "mfrmr"))
fit <- example_fit$fit
fit_eb <- apply_empirical_bayes_shrinkage(fit)
fit_eb$shrinkage_report
# Look for:
# - `Tau2` is the estimated between-level prior variance per facet.
#   `Tau2 = 0` fully pools eligible estimates at zero
#   (`MeanShrinkage = 1`). Zero plug-in SE is not perfect precision.
# - `MeanShrinkage` near 0 = little movement, near 1 = heavy pooling
#   toward the chosen zero target; inspect the SEs and prior variance.
# - `EffectiveDF` sums retained weights over eligible levels, conditional
#   on the chosen prior variance; it is not degrees of freedom for testing.
head(fit_eb$facets$others[, c("Facet", "Level", "Estimate",
                               "ShrunkEstimate", "ShrinkageFactor")])
# Look for: rows where `ShrinkageFactor` is large (close to 1) had
#   their estimates pulled most strongly toward the facet mean (0).

Description

Generic for converting objects to a flextable

Usage

as_flextable(x, ...)

Arguments

Value

A flextable object (concrete return type from the underlying method, e.g. ⁠[as_flextable.apa_table()]⁠ returns a flextable ready for flextable::save_as_docx()).

See Also

as_flextable.apa_table() for the apa_table method; as_kable() for a knitr::kable-targeted alternative; apa_table() for constructing an apa_table in the first place.

Examples

# Optional integration example (loads additional modelling or rendering packages).
# example("as_flextable", run.dontrun = TRUE)
## Not run: 
tbl <- structure(
  list(
    table = data.frame(Term = c("Rater A", "Rater B"), Estimate = c(-0.12, 0.18)),
    caption = "Facet estimates",
    note = "Toy values for formatting only."
  ),
  class = "apa_table"
)
if (requireNamespace("flextable", quietly = TRUE)) {
  invisible(as_flextable(tbl))
}

## End(Not run)

Description

Produces a Word / PowerPoint-friendly flextable with the caption and note wired in. Requires flextable (in Suggests).

Usage

## S3 method for class 'apa_table'
as_flextable(x, ...)

Arguments

Value

A flextable object, or a message when flextable is unavailable.

See Also

as_kable.apa_table(), apa_table().

Examples

tbl <- structure(
  list(
    table = data.frame(Term = c("Rater A", "Rater B"), Estimate = c(-0.12, 0.18)),
    caption = "Facet estimates",
    note = "Toy values for formatting only."
  ),
  class = "apa_table"
)
if (requireNamespace("flextable", quietly = TRUE)) {
  invisible(as_flextable(tbl))
}

Description

as_ggplot() is an optional renderer for an mfrm_plot_data object or an object whose plot() method supports draw = FALSE. Base graphics remain the default; the returned ggplot can be restyled or composed downstream.

Usage

as_ggplot(x, type = NULL, component = NULL, ...)

## Default S3 method:
as_ggplot(x, type = NULL, component = NULL, ...)

## S3 method for class 'mfrm_design_evaluation'
as_ggplot(x, type = NULL, component = NULL, ...)

## S3 method for class 'mfrm_design_evaluation_plot_data'
as_ggplot(x, type = NULL, component = NULL, ...)

## S3 method for class 'mfrm_signal_detection'
as_ggplot(x, type = NULL, component = NULL, ...)

## S3 method for class 'mfrm_signal_detection_plot_data'
as_ggplot(x, type = NULL, component = NULL, ...)

## S3 method for class 'mfrm_plot_data'
as_ggplot(x, type = NULL, component = NULL, ...)

## S3 method for class 'mfrm_slope_intervals'
as_ggplot(x, type = NULL, component = NULL, ...)

## S3 method for class 'mfrm_curve_intervals'
as_ggplot(x, type = NULL, component = NULL, ...)

## S3 method for class 'mfrm_gpcm_bootstrap'
as_ggplot(x, type = NULL, component = NULL, ...)

## S3 method for class 'mfrm_results'
as_ggplot(x, type = NULL, component = NULL, ...)

Arguments

Details

Bubble conversion uses the saved radii, facet colours, facet order and reference lines. It preserves relative circle radii, rather than replacing them with equal-sized points. Base circles use plot units; ggplot point sizes use physical units, so absolute sizes need not match across devices. Payloads lacking matching saved radii or facet colours must be recreated with plot(x, type = "bubble").

Dedicated conversions are provided for Wright maps, theta-to-expected-score pathways, fit-statistic-to-measure pathways, category characteristic curves, bubble charts, DIF/DFF summaries and heatmaps, portable-calibration score review plots, testlet calibration, shared-rater severity and bootstrap intervals, conditional Person scores, and multivariate D-study comparisons. Matched ordinary/extended facet comparisons preserve centered effects, paired/difference views, excluded rows, display controls and alternative text; see plot.mfrm_extended_comparison(). Portable-calibration score views preserve their selected axes, stored endpoints, review symbols and unplotted dispositions with either the default or component = "data". Use plot_data() for other components; conversion does not rescore Persons or change their interval level. Agreement conversions retain the selected exact/correlation/difference metric, saved pair ordering and review thresholds, model-expected overlay, and missing comparisons. Agreement heatmaps use the saved facet labels and palette with a fixed metric scale. Use plot_data() for other tables in an agreement payload; automatic conversion does not guess a different statistic from the first numeric column. Linking conversions preserve the selected drift, wave estimates and stored interval endpoints, cumulative offsets, or common-element counts. Missing offsets do not become zero or connect a line across an unavailable link. Drift heatmaps use their saved residual matrix. Use plot_data() for other tables; conversion does not recompute offsets or their uncertainty. Threshold-sensitivity tile and curve payloads also preserve their selected view, display controls and text alternatives; see plot.mfrm_screening_sensitivity(). Response residual displays preserve their posterior or conditional plug-in probability definition and lack of reference cutoffs. Paired views retain available Infit when Outfit is undefined; see plot.mfrm_response_diagnostics(). Testlet conversions retain unavailable rows, prior-only symbols and the conditional-interval note; they do not estimate diagnostics or add calibration uncertainty. Extended-model interval, precision and empirical-distribution views preserve display settings, data exclusions and text alternatives; see plot.mfrm_testlet_scores(). Bootstrap conversion retains infinite endpoints as arrows and ordinary intervals as a dashed comparison. D-study conversions preserve G/Phi or SEM panels, fixed-count groups, score units, and unavailable estimates. They retain supplied D-study coefficient ribbons and interval gaps. Converting a saved plot object does not refit models or calculate intervals. G-study component conversions retain raw negative points, component intervals and square outlines for unavailable intervals. For a G-study or D-study input, intervals = TRUE is forwarded to its plot method and computes intervals before conversion; also supply assumption = "normal". For multivariate D-study plots, an explicit component = "series" keeps this dedicated conversion; use plot_data() for other tables. Difference-interval plots from mfrm_multivariate_d_compare() use their base plot() method or plot_data(); automatic conversion is not supported. Severity profiles preserve their saved normal bounds, ordered levels, guide bands, fixed points, unavailable states and precision notes with the default or component = "data". The word "profile" does not imply a likelihood-profile interval. Conversion does not recalculate bounds or qualify exploratory JML or fallback/regularized uncertainty for inference. Residual-PCA scree, reference, excess and loadings plots have dedicated conversions with the default or component = "data". They preserve saved eigenvalues, reference cutoffs, signed loadings and display encodings; method identifiers are not inferred as plotted variables. No PCA is recomputed. External-feature PCA scree, scores and loadings views have dedicated conversions. They preserve selected axes, retained-component symbols, saved group colours and shapes, labels and transformation metadata. Scores retain equal axis units; loading coefficients are not correlations. The default and component = "table" use the same complete PCA view. Use ggplot2::labs(title = NULL, subtitle = NULL) to remove headings from the returned ggplot. No PCA or clustering is refitted. External-feature dendrogram conversion retains the stored merges, heights, leaf order and group boxes, including tied-height cuts. The default and component = "tree" use the same complete view. Labels follow the saved setting; no tree is refitted or partition selected. Source metadata remain available with plot_data(). Heights are not significance or branch support. Imputation co-membership heatmaps have dedicated conversion with the default or component = "matrix". Values retain their all-imputation denominator, selected ID order and fixed zero-to-one scale. Unavailable cells use grey fill and crosses; zero is a valid fraction. The saved legend, labels and imputation count are reused, without recomputing or pooling partitions. Silhouette conversion preserves widths on a fixed minus-one-to-one scale, their saved order and overall-mean line. Numeric profiles preserve original-unit means, medians and counts, using separate symbols with slight vertical offsets. Categorical profiles preserve the category order, unused levels, group counts and within-group proportions on a zero-to-one scale. The default and component = "table" keep the full selected view; categorical profiles also accept component = "matrix". These summaries do not estimate uncertainty, rater quality or sampling stability. For main-effects mfrm_d_study results, use the base plot() method or plot_data() for custom graphics; automatic conversion is refused because generic column selection does not preserve those design comparisons. Fit-measure measure_ci plots retain each facet/level row and its saved normal-band endpoints, including exploratory JML bands. Open diamonds mark fixed values and crosses mark estimates without intervals. Conversion does not recompute bounds or change their inferential qualification. The default and component = "table" preserve the complete view and its returned notes. Pooled fixed-facet MI intervals have dedicated conversion with the default or component = "table". Saved t-interval endpoints, degrees of freedom, contrasts and complete-data information cautions are retained, with no repooling or interval recalculation. Open diamonds mark fixed targets without intervals, crosses mark missing intervals, and arrows mark infinite bounds. Arrow tips are plotting limits, not replacement finite confidence limits. Selecting component does not enable unsupported conversions: it cannot preserve D-study differences through generic column selection. Use plot_data() to extract the table and specify the axes, intervals and grouping explicitly in custom graphics. Residual Q-Q plots require the base renderer or explicitly specified custom quantile axes; automatic conversion is unavailable, even with component. Residual-correlation heatmaps retain the selected facet on both axes, a fixed -1 to 1 colour scale, and grey unavailable cells; their pair table and notes retain overlap counts and interpretation qualifications. Saved rating-scale counts and thresholds, and category-structure counts, retain the selected values, missing outcomes, expected-count overlay and saved palette. Use the default or the matching category_table / threshold_table component. Threshold figures show point estimates; stored SEs and CI-use flags remain in the payload, not as curve or gap intervals. Category-curve bundles (including their overview), category-structure boundary views and threshold ladders require their base renderer or custom graphics from plot_data(); automatic conversion is unavailable even with component. Generic column selection would mix distinct curve targets or lose the selected view. Fitted CCC plots from plot(fit, type = "ccc") retain their dedicated conversion.

Bias scatter/ranked views preserve their axes, saved review rules and conditional interval bounds, including dashed search-limited bounds. Bias heatmaps, histograms and facet profiles require base/custom graphics; the generic converter cannot combine their distinct layouts safely. Halo distributions, selected per-rater summaries and signed retained-edge matrices preserve the selected metric, missing comparisons and caveats. Halo network conversion is unavailable. Save the complete analysis result and recreate its plot to obtain current interpretation notes.

Shrinkage plots from plot(fit, type = "shrinkage") preserve both estimates for each facet level, their adjustment direction and saved descriptive bands. Use the default or component = "data". Conversion does not turn plug-in bands into calibrated confidence intervals. Recreate older plots from a fit with a refreshed shrinkage report first.

Other draw-free payloads use a conservative tabular fallback; inspect plot_data_components() when automatic inference is not appropriate. Titles, subtitles and captions are omitted by default and retained as metadata. Use plot_data(p, "notes") to retrieve descriptions and cautions. Explicitly enabled text is wrapped at 72 columns; use ggplot2::labs() for custom labels or narrower exports. Fit plots created with show_title = FALSE or show_notes = FALSE retain those settings on conversion. To set them when converting a fit directly, pass the flags through .... Interpretation and display notes remain available in attr(plot, "mfrmr_notes"), when present in the source payload. Wright, pathway, and CCC conversions share the fit-family series palette, honour monochrome presets and supplied palette overrides, and preserve line-type or point-shape distinctions. CCC categories retain their source order. When slope_aes = "colour", slopes instead use a continuous viridis scale (a grey gradient in monochrome); categories still have line types.

Fixed-facet RSM/PCM intervals retain the selected covariance method, level, status symbols and optional ordinary comparison. Titles, subtitles, captions, reference lines and legends can be changed through the source plot arguments. Methods use both line types and vertical offsets, including in monochrome.

Value

A ggplot2 plot object, with a mfrmr_notes attribute when the source plot payload contains a notes table. Paired Wright/CCC payloads from plot(reference, comparison = comparison) retain their selected comparison or difference view, group selection and monochrome panel policy.

Interactive views

The visual gallery shows how to pass a supported ggplot to plotly::ggplotly() for browser-based hovering, zooming and legend filtering. Plotly is an optional downstream renderer. Keep the original mfrmr plot and its plot_data() notes alongside the widget; conversion does not recalculate estimates or intervals. Check labels, custom geoms and interval layers after conversion. The gallery includes static figures, numeric tables and an HTML export recipe.

See Also

mfrmr_output_guide() with scope = "plots" for selected purpose-based routes, conversion status and alternatives.

Examples

# Long-running workflow; run explicitly when needed.
# example("as_ggplot", run.dontrun = TRUE)
## Not run: 
toy <- load_mfrmr_data("example_core")
# A balanced slice retains every Rater and Criterion while running quickly.
toy <- toy[toy$Person %in% unique(toy$Person)[1:12], , drop = FALSE]
fit <- fit_mfrm(toy, "Person",
                c("Rater", "Criterion"), "Score", maxit = 30)
as_ggplot(fit, type = "wright")
as_ggplot(fit, type = "fit_pathway", include_person = TRUE)
as_ggplot(plot(fit, type = "ccc", draw = FALSE))

## End(Not run)

Description

Generic for converting objects to a knitr::kable

Usage

as_kable(x, ...)

Arguments

Value

A knitr::kable object (concrete return type from the underlying method, e.g. ⁠[as_kable.apa_table()]⁠ returns a kableExtra object when the package is installed).

See Also

as_kable.apa_table() for the apa_table method; as_flextable() for a flextable-targeted alternative; apa_table() for constructing an apa_table in the first place.

Examples

tbl <- structure(
  list(
    table = data.frame(Term = c("Rater A", "Rater B"), Estimate = c(-0.12, 0.18)),
    caption = "Facet estimates",
    note = "Toy values for formatting only."
  ),
  class = "apa_table"
)
if (requireNamespace("knitr", quietly = TRUE)) {
  invisible(as_kable(tbl))
}

Description

Renders the table data for direct inclusion in RMarkdown, Quarto, or HTML reports, wiring the caption and note slots into the standard APA placement (caption above, note below). When kableExtra is installed the note is attached as a footer for HTML or LaTeX. Otherwise the complete rendered table is collapsed to one knitr_kable string before a single note is appended.

Usage

## S3 method for class 'apa_table'
as_kable(x, format = c("pipe", "html", "latex"), digits = 3L, ...)

Arguments

Value

A knitr_kable object ready to be printed inline in a report, or a message when knitr is unavailable.

See Also

as_flextable.apa_table(), apa_table().

Examples

tbl <- structure(
  list(
    table = data.frame(Term = c("Rater A", "Rater B"), Estimate = c(-0.12, 0.18)),
    caption = "Facet estimates",
    note = "Toy values for formatting only."
  ),
  class = "apa_table"
)
if (requireNamespace("knitr", quietly = TRUE)) {
  invisible(as_kable(tbl))
}

Description

Returns all facet-level estimates (person and others) in a single tidy data.frame. Person rows retain their original identifiers; review or transform them before writing the result outside a controlled analysis environment.

Usage

## S3 method for class 'mfrm_fit'
as.data.frame(x, row.names = NULL, optional = FALSE, ...)

Arguments

Details

This method returns four columns (Facet, Level, Estimate, Extreme) so that the result is easy to inspect, join, or write to disk.

Value

A data.frame with columns Facet, Level, Estimate, and Extreme. The Extreme column is populated for person rows from the extreme-score flag ("Min" / "Max" / NA); non-person facet rows carry NA in that column by design. Corrected-JML fits additionally retain estimator, correction order, point/covariance status, uncertainty target and structural RootSE columns. Person RootSE is missing; extreme Person flags are "low" / "high" / "none". The estimation_note attribute explains the local-root uncertainty and remaining bias. Ordinary diagnostic attachments do not apply to this estimator; see the Corrected JML section of fit_mfrm.

Interpreting output

Person estimates are returned with Facet = "Person". All non-person facets are stacked underneath in the same schema.

Typical workflow

1. Fit a model with fit_mfrm().

2. Convert with as.data.frame(fit) for a compact long-format export.

3. Join additional diagnostics later if you need SE or fit statistics.

See Also

fit_mfrm, export_mfrm_tables

Examples

library(mfrmr)

# Saved synthetic MML/RSM analysis of example_operational.
# For your own data, fit a model first; see ?fit_mfrm.
example_fit <- readRDS(system.file("examples", "operational-mml.rds", package = "mfrmr"))
fit <- example_fit$fit

# Extract estimates and select the rows to display
estimates <- as.data.frame(fit)
head(subset(estimates, Facet == "Person")) # First six persons
subset(estimates, Facet == "Rater")       # All raters
subset(estimates, Facet == "Criterion")   # All criteria

Description

Converts the numerical output from simulate_mfrm_recovery() into a reviewer-facing adequacy checklist. The goal is not to impose one universal pass/fail rule; it is to make the main user questions explicit: Did the runs finish? Did the fitted models converge? Are uncertainty summaries available? Are coverage and Monte Carlo precision plausible? If practical RMSE or bias limits are supplied, which parameter groups need follow-up? For GPCM, which slope-regime generator condition frames the recovery evidence? This only reads an existing simulation result. It does not generate data, refit models, add replications or repair unavailable uncertainty estimates.

Usage

assess_mfrm_recovery(
  x,
  min_reps = 30,
  min_success_rate = 0.95,
  min_convergence_rate = 0.95,
  min_se_available = 0.8,
  coverage_target = 0.95,
  coverage_tolerance = 0.05,
  max_mcse_rmse_ratio = 0.25,
  max_rmse = NULL,
  max_abs_bias = NULL,
  top_n = 6,
  digits = 3,
  ...
)

## S3 method for class 'mfrm_recovery_assessment'
plot(
  x,
  y = NULL,
  type = c("status", "metrics"),
  metric = c("rmse", "bias", "coverage", "se_available", "mcse_rmse"),
  draw = TRUE,
  ...
)

Arguments

Details

RMSE and bias adequacy depends on the substantive scale and the use case, so the function does not mark them as failed unless the user supplies max_rmse or max_abs_bias. Without those limits, the corresponding rows are marked not_assessed and the next action asks the user to set practical thresholds when a decision depends on the metric. BiasBasis explains groups whose mean bias is fixed at zero by identification; even a supplied bias limit cannot turn that identity into recovery evidence. Threshold names are matched from the most specific group to its type or default; see max_rmse for examples.

Metrics are recomputed from row-level evidence. Coverage is conditional on the available centered facet intervals, not all planned targets. uncertainty_review retains counts, CoverageAvailabilityRate, McseCoverage95 and CIReasons alongside its descriptive checklist. See simulate_mfrm_recovery() for the estimands, independent-trial MCSE, unsupported interval groups and handling of legacy saved results. Checklist tolerances do not replace the Monte Carlo uncertainty or qualify an estimator.

The condition_review table is generator metadata for interpreting the recovery run. For GPCM, GPCMSlopeRegime, StressLevel, and generated score-category support describe the data-generating condition; they are not model-fit tests and they are not literature-derived adequacy cut points. condition_reporting_notes turns those generator conditions into reporter-facing caveats, such as high-dispersion slope stress or sparse generated score support.

The optional diagnostic_review table is available when simulate_mfrm_recovery() was called with include_diagnostics = TRUE. It summarizes fit and separation operating characteristics as diagnostic context only. Its availability fields do not mean that fit or separation values are adequate, and those rows do not enter the recovery adequacy status. diagnostic_reporting_notes should be read first when drafting fit/separation language because it separates zero separation/reliability, absolute fit-ZSTD flags, and df-sensitive ZSTD flags from recovery criteria.

plot.mfrm_recovery_assessment() is a user-facing review aid. Use type = "status" first to see where checklist attention is needed, then type = "metrics" to inspect the parameter groups behind RMSE, bias, coverage, standard-error availability, or Monte Carlo precision statuses. The intended reading order is summary(recovery_review), then condition_reporting_notes and condition_review, then diagnostic_reporting_notes and diagnostic_review, then the status plot, then the metric plot, then the row-level recovery table for the parameter groups that need follow-up. When draw = FALSE, the plot data also include reading_order, guidance, condition/diagnostic handoff tables, and user-facing plot tables such as section_status for status plots and metric_review for metric plots.

Value

An object of class mfrm_recovery_assessment with:

• overview: compact run-level status.

• checklist: reviewer-facing adequacy checks.

• condition_review: generator-condition metadata, including GPCM slope-regime interpretation and generated score-category support when available.

• condition_reporting_notes: reporter-facing generator-condition caveats separated from parameter-recovery conclusions.

• diagnostic_reporting_notes: reporter-facing fit/separation caveats retained as diagnostic context rather than recovery criteria.

• diagnostic_review: optional fit/separation operating-characteristic context when retained by simulate_mfrm_recovery().

• metric_review: parameter-group metric checks.

• uncertainty_review: coverage / SE availability interpretation with expected/available counts, coverage MCSE and unavailable reasons.

• reading_order: recommended first-read order for the summary, condition, plot, and row-level recovery outputs.

• next_actions: short action list sorted by severity.

• thresholds: thresholds used for the assessment.

See Also

simulate_mfrm_recovery(), plot.mfrm_recovery_simulation()

Examples

# Long-running workflow; run explicitly when needed.
# example("assess_mfrm_recovery", run.dontrun = TRUE)
## Not run: 
rec <- simulate_mfrm_recovery(
  n_person = 12,
  n_rater = 2,
  n_criterion = 2,
  reps = 1,
  maxit = 30,
  seed = 123
)
assess_mfrm_recovery(rec, min_reps = 1, max_rmse = 1)

# Read the GPCM generator condition separately from recovery adequacy.
gpcm_spec <- build_mfrm_sim_spec(
  n_person = 14,
  n_rater = 2,
  n_criterion = 2,
  raters_per_person = 2,
  model = "GPCM",
  step_facet = "Criterion",
  slope_facet = "Criterion",
  slopes = c(0.85, 1.15),
  assignment = "crossed"
)
gpcm_rec <- suppressWarnings(simulate_mfrm_recovery(
  sim_spec = gpcm_spec,
  reps = 1,
  fit_method = "MML",
  quad_points = 5,
  maxit = 12,
  include_diagnostics = TRUE,
  include_person = FALSE,
  seed = 456
))
gpcm_review <- assess_mfrm_recovery(
  gpcm_rec,
  min_reps = 1,
  max_rmse = c(slope = 2),
  max_abs_bias = c(slope = 1),
  min_se_available = NULL,
  max_mcse_rmse_ratio = NULL
)
gpcm_review$condition_reporting_notes[, c(
  "ConditionArea", "ReportingAttention", "ConditionFinding"
)]
gpcm_review$condition_review[, c(
  "Model", "GPCMSlopeRegime", "StressLevel", "ScoreSupportStatus"
)]
gpcm_review$diagnostic_reporting_notes[, c(
  "Facet", "ReportingAttention", "DiagnosticFinding"
)]
summary(gpcm_review)$reading_order

## End(Not run)

Description

Prefer score_comparison_table() for new analyses. This earlier interface retains its argument names, automated/human output columns and default Williamson (2012) numerical screen. The generic interface uses score/reference names and applies no thresholds unless requested.

Usage

automated_scoring_table(
  automated,
  human,
  human_reference = NULL,
  categories,
  automated_categories = automated,
  smd = c("pooled", "human"),
  qwk_min = 0.7,
  correlation_min = 0.7,
  max_degradation = 0.1,
  max_abs_smd = 0.15,
  fit = NULL,
  diagnostics = NULL,
  rater_facet = NULL,
  automated_rater = NULL
)

Arguments

Value

The same list structure as score_comparison_table(), retaining the earlier automated/human column and category-table names. Saved results created through the earlier interface can still be read directly.

See Also

score_comparison_table()

Description

Build a bias-cell count report

Usage

bias_count_table(
  bias_results,
  min_count_warn = 10,
  branch = c("original", "facets"),
  fit = NULL
)

Arguments

Details

This helper summarizes how many observations contribute to each bias-cell estimate and flags sparse cells.

Branch behavior:

• "facets": keeps legacy manual-aligned column labels (Sq, ⁠Observd Count⁠, ⁠Obs-Exp Average⁠, ⁠Model S.E.⁠) for side-by-side comparison with external workflows.

• "original": keeps compact field names (Count, BiasSize, SE) for custom QC workflows and scripting.

Value

A named list with:

• table: cell-level counts with low-count flags

• by_facet: named list of counts aggregated by each interaction facet

• by_facet_a, by_facet_b: first two facet summaries (legacy compatibility)

• summary: one-row summary

• thresholds: applied thresholds

• branch, style: output branch metadata

• fit_overview: optional one-row fit metadata when fit is supplied

Interpreting output

• table: cell-level contribution counts and low-count flags.

• by_facet: sparse-cell structure by each interaction facet.

• summary: overall low-count prevalence.

• fit_overview: optional run context (when fit is supplied).

Low-count cells should be interpreted cautiously because bias-size estimates can become unstable with sparse support.

Typical workflow

1. Estimate bias with estimate_bias().

2. Build bias_count_table(...) in desired branch.

3. Review low-count flags before interpreting bias magnitudes.

Further guidance

For a plot-selection guide and a longer walkthrough, see mfrmr_visual_diagnostics and vignette("mfrmr-visual-diagnostics", package = "mfrmr").

Output columns

The table data.frame contains, in the legacy-compatible branch:

FacetA, FacetB

Interaction facet level identifiers; generic names for the two interaction facets.

Sq

Sequential row number.

Observd Count

Number of observations for this cell.

Obs-Exp Average

Observed minus expected average for this cell.

Model S.E.

Standard error of the bias estimate.

Infit, Outfit

Fit statistics for this cell.

LowCountFlag

Logical; TRUE when count < min_count_warn.

The summary data.frame contains:

InteractionFacets

Names of the interaction facets.

Cells, TotalCount

Number of cells and total observations.

LowCountCells, LowCountPercent

Number and share of low-count cells.

See Also

estimate_bias(), unexpected_after_bias_table(), build_fixed_reports(), mfrmr_visual_diagnostics

Examples

# Saved synthetic JML/RSM analysis of example_bias.
# For your own data, fit a model first; see ?fit_mfrm.
example_fit <- readRDS(system.file("examples", "bias-jml.rds", package = "mfrmr"))
fit <- example_fit$fit
diag <- example_fit$diagnostics
bias <- estimate_bias(fit, diag, facet_a = "Rater", facet_b = "Criterion", max_iter = 2)
t11 <- bias_count_table(bias)
t11_facets <- bias_count_table(bias, branch = "facets", fit = fit)
summary(t11)
p <- plot(t11, draw = FALSE)
p2 <- plot(t11, type = "lowcount_by_facet", draw = FALSE)
if (interactive()) {
  plot(
    t11,
    type = "cell_counts",
    draw = TRUE,
    main = "Bias Cell Counts (Customized)",
    palette = c(count = "#2b8cbe", low = "#cb181d"),
    label_angle = 45
  )
}

Description

Bundles the ranked flagged-cells view of a bias-interaction run for downstream printing and plotting. The three sibling reports in this family are intentionally distinct:

• bias_interaction_report() (this one) = FACETS Table 13: a ranked list of interaction cells with t, ⁠bias size⁠, and screening tail area – use when reviewing which ⁠(facet_a, facet_b)⁠ cells deserve follow-up.

• bias_iteration_report() = iteration history / convergence trace for the bias recalibration (FACETS Table 9 territory) – use when diagnosing whether the bias run itself stabilised.

• bias_pairwise_report() = pairwise contrast table for a target facet (FACETS Table 14 territory) – use when comparing levels within a facet while controlling for the other.

Usage

bias_interaction_report(
  x,
  diagnostics = NULL,
  facet_a = NULL,
  facet_b = NULL,
  interaction_facets = NULL,
  max_abs = 10,
  omit_extreme = TRUE,
  max_iter = 4,
  tol = 0.001,
  top_n = 50,
  abs_t_warn = 2,
  abs_bias_warn = 0.5,
  p_max = 0.05,
  sort_by = c("abs_t", "abs_bias", "prob")
)

Arguments

Details

Preferred bundle API for interaction-bias diagnostics. The function can:

• use a precomputed bias object from estimate_bias(), or

• estimate internally from mfrm_fit + facet specification.

Save the complete result with base::saveRDS(); plot(out, type = "scatter") and type = "ranked" reuse its flags and thresholds. Rebuild the report from the saved bias estimates to change screening rules; conflicting plotting-time thresholds are rejected. Both views support as_ggplot(), including conditional bounds requested with show_ci = TRUE. Heatmaps, histograms and facet profiles use base graphics or custom plot_data() graphics. Export the bundle's data frames with utils::write.csv().

Value

A named list with bias-interaction plotting/report components. Class: mfrm_bias_interaction.

Interpreting output

Focus on ranked rows where multiple screening criteria converge:

• large absolute t statistic

• large absolute bias size

• small screening tail area

The bundle is optimized for downstream summary() and plot(x, type = "bias") views.

Typical workflow

1. Run estimate_bias() (or provide mfrm_fit here).

2. Build bias_interaction_report(...).

3. Review summary(out) and visualize with plot(x, type = "bias").

See Also

estimate_bias(), build_fixed_reports(), plot(x, type = "bias")

Examples

# Saved synthetic JML/RSM analysis of example_bias.
# For your own data, fit a model first; see ?fit_mfrm.
example_fit <- readRDS(system.file("examples", "bias-jml.rds", package = "mfrmr"))
fit <- example_fit$fit
diag <- example_fit$diagnostics
bias <- estimate_bias(fit, diag, facet_a = "Rater", facet_b = "Criterion", max_iter = 2)
out <- bias_interaction_report(bias, top_n = 10)
summary(out)
p_bi <- plot(out, draw = FALSE)
p_bi$data$plot

Description

This report is NOT an alias of bias_interaction_report() despite the similar name. It focuses on the recalibration path of a bias run: iteration table, convergence summary, and orientation review. Use this to confirm that the bias recalibration itself converged; use bias_interaction_report() to review the ranked flagged cells from the converged run.

Usage

bias_iteration_report(
  x,
  diagnostics = NULL,
  facet_a = NULL,
  facet_b = NULL,
  interaction_facets = NULL,
  max_abs = 10,
  omit_extreme = TRUE,
  max_iter = 4,
  tol = 0.001,
  top_n = 10
)

Arguments

Details

This report focuses on the recalibration path used by estimate_bias(). It provides a package-native counterpart to legacy iteration printouts by exposing the iteration table, convergence summary, and orientation review in one bundle.

Value

A named list with:

• table: iteration history

• summary: one-row convergence summary

• orientation_review: interaction-facet sign review

• settings: resolved reporting options

• direction_note: one-line interpretive note describing which direction the iteration moved (carried from the bias estimator; empty string when the underlying estimator does not emit one)

• recommended_action: one-line recommended action label (e.g. "converged", "increase max_iter"); empty string when the underlying estimator does not emit one

See Also

estimate_bias(), bias_interaction_report(), build_fixed_reports()

Examples

# Saved synthetic JML/RSM analysis of example_bias.
# For your own data, fit a model first; see ?fit_mfrm.
example_fit <- readRDS(system.file("examples", "bias-jml.rds", package = "mfrmr"))
fit <- example_fit$fit
diag <- example_fit$diagnostics
out <- bias_iteration_report(fit, diagnostics = diag, facet_a = "Rater", facet_b = "Criterion")
summary(out)

Description

Compare each target element's local measure across pairs of contexts. With target Rater and context Criterion, a row compares the same rater's local severity on two criteria. With target Criterion and context Rater, it compares the same criterion's local difficulty across two raters. The calculation uses conditional bias cells and does not establish numerical equivalence to FACETS Table 14.

Usage

bias_pairwise_report(
  x,
  diagnostics = NULL,
  facet_a = NULL,
  facet_b = NULL,
  interaction_facets = NULL,
  max_abs = 10,
  omit_extreme = TRUE,
  max_iter = 4,
  tol = 0.001,
  target_facet = NULL,
  context_facet = NULL,
  top_n = 50,
  p_max = 0.05,
  sort_by = c("abs_t", "abs_contrast", "prob"),
  df_method = c("facets", "classical")
)

Arguments

Details

This helper exposes the pairwise contrast table that was previously only reachable through fixed-width output generation. It is available only for 2-way interactions. Rows without available orientation or precision are retained with missing screening results. top_n limits displayed rows; the summary always describes the complete pairwise family.

Value

A named list with:

• table: pairwise contrast rows

• summary: one-row contrast summary

• orientation_review: interaction-facet sign review

• settings: resolved reporting options

• direction_note: one-line interpretive note describing the dominant pairwise-contrast direction (carried from the underlying bias estimator; empty string when not applicable)

• recommended_action: one-line recommended-action label (e.g. routing the user to follow-up review of the largest flagged pairs); empty string when the underlying estimator does not emit one

Interpreting output

• table: one row per target element and context pair, with Target, Context1, Context2, ⁠Local Measure1⁠, ⁠Local Measure2⁠, Contrast, SE, t, d.f., Prob., Flag, and calculation/status metadata. The default sort places the largest absolute t first and unavailable values last.

• summary: full-family Contrasts, ScreenAvailable, ScreenUnavailable, Flagged, MeanAbsContrast, and MeanAbsT. Flagged counts known screen-positive rows and is missing when none can be screened. A zero count does not clear unavailable comparisons.

• orientation_review carries the same facet-orientation sign review as the parent estimate_bias() run.

• Flags are unadjusted exploratory screens. This report does not apply multiplicity correction or establish the absence of differential bias.

Typical workflow

1. Fit and diagnose the model.

2. Run estimate_bias() to get the underlying interaction effects.

3. Pass that result to bias_pairwise_report() with the intended target and context facets.

4. Read out$summary for availability and out$table for direction, magnitude and uncertainty of the local comparisons.

5. For the ranked flagged-cells view (which (rater, item) pairs have the largest local bias), use bias_interaction_report() on the same estimate_bias() output.

Standard-error caveat

The contrast standard error is computed as SE(b_i - b_j) = sqrt(SE_i^2 + SE_j^2) – the independence approximation conditional on the calibration. The general variance of a difference also subtracts 2 * Cov(b_i, b_j): a negative covariance increases that variance, while a positive covariance decreases it. Neither its sign nor its magnitude is established here. There is no general claim that these screens are conservative. The cell shifts are estimated separately with fixed baseline parameters, not by a jointly estimated sum-to-zero bias model.

References

• Linacre, J. M. (1989). Many-Facet Rasch Measurement. MESA Press.

• Eckes, T. (2005). Examining rater effects in TestDaF writing and speaking performance assessments: A many-facet Rasch analysis. Language Assessment Quarterly, 2(3), 197-221.

• Myford, C. M., & Wolfe, E. W. (2003). Detecting and measuring rater effects using many-facet Rasch measurement: Part I. Journal of Applied Measurement, 4(4), 386-422.

• Myford, C. M., & Wolfe, E. W. (2004). Detecting and measuring rater effects using many-facet Rasch measurement: Part II. Journal of Applied Measurement, 5(2), 189-227.

Pairwise calculation

For each target element, compare its local measure in two contexts. If the target has predictor sign s, its local measure is target + bias/s; the contrast is Context1 minus Context2. A negatively oriented rater therefore has local severity severity - bias. TargetSign comes from the saved bias result's orientation_review. Missing orientation withholds local measures, contrasts and tail areas; it is not inferred from a facet name.

SE1 and SE2 are conditional cell-information SEs with baseline calibration held fixed. The contrast SE is sqrt(SE1^2 + SE2^2) under conditional independence of the two cells. Calibration uncertainty and cross-cell covariance are not propagated. ⁠Target S.E.⁠ describes the original target estimate separately and is not added to the conditional local SE; a fixed target can still have uncertain local bias.

The default df_method = "facets" uses the modified Rasch-Welch formula described in the FACETS t-statistics manual (von Davier, 2025). Writing v1 = SE1^2, v2 = SE2^2 and n1, n2 for weighted counts, its df is (v1+v2)^2 * (n1+n2-2) / ((n1+n2+2) * (v1^2/(n1+1)+v2^2/(n2+1))). The tail area uses its integer part, matching executed FACETS 4.5.0 unequal-count examples. Values are stabilized to ten decimal places before truncation to avoid numerical noise at integer boundaries. DFUnrounded retains the original continuous value. df_method = "classical" retains the former Welch-Satterthwaite reference with n1-1, n2-1 and no integer truncation. DFMethod and DFBasis identify the choice; equal counts alone cannot distinguish the formulas. ObsN1 and ObsN2 are raw row counts. Both weighted counts must exceed one, both SEs must be positive and finite, and the final df must be positive. Otherwise the screen remains unavailable. Fractional weights follow the same rule but are outside the executed comparison. These are conditional screens, not calibrated fairness tests or general FACETS equivalence.

Saved pairwise results

Current tables carry CalculationVersion = 3. Recreate older fixed or pairwise reports using build_fixed_reports() or bias_pairwise_report(). If the bias bundle predates baseline-score corrections or lacks orientation, recreate estimate_bias() from the saved fit and current diagnostics first. No main calibration refit is needed. Keep historical values separately; summaries and fixed-report plots reject older pairwise calculations.

See Also

estimate_bias(), bias_interaction_report(), build_fixed_reports()

Examples

# Saved synthetic JML/RSM analysis of example_bias.
# For your own data, fit a model first; see ?fit_mfrm.
example_fit <- readRDS(system.file("examples", "bias-jml.rds", package = "mfrmr"))
fit <- example_fit$fit
diag <- example_fit$diagnostics
out <- bias_pairwise_report(fit, diagnostics = diag, facet_a = "Rater", facet_b = "Criterion")
s <- summary(out)
s$summary
# Check ScreenUnavailable before reading the known Flagged count.
# Zero flags do not establish equivalence or absence of bias.
head(out$table)
# For target Rater, compare that rater's local severity across criteria.
# Read Context1 minus Context2 and the conditional SE before follow-up.

Description

Simulate Persons and scores on the analyzed rating assignment and refit the model. With no null_fit, saved draws support slope intervals. With a matched PCM null_fit, simulate under PCM and refit both models for a bootstrap LRT.

Usage

bootstrap_mfrm_gpcm(fit, nsim = 499L, seed, null_fit = NULL)

## S3 method for class 'mfrm_gpcm_bootstrap'
confint(
  object,
  parm = "slopes",
  level = 0.95,
  scale = c("relative", "standardized"),
  contrasts = NULL,
  contrast_scale = c("ratio", "difference"),
  simultaneous = c("none", "bonferroni"),
  ...
)

## S3 method for class 'mfrm_gpcm_bootstrap'
summary(object, ...)

## S3 method for class 'mfrm_gpcm_bootstrap'
print(x, ...)

Arguments

Details

One ability is drawn per Person from the fitted conditional-normal population, shared by all that Person's rows. Facets and covariates are fixed. The observed assignment is preserved; omitted ratings remain omitted. No assignment or missingness mechanism is simulated. Population coefficients and variance are reestimated if estimated in the source model. Retained formula, factor coding and person data must reproduce the population design.

confint() uses basic bootstrap errors on the log scale for slopes/ratios, and the identity scale for differences. The point estimate minus reversed empirical error quantiles (type 1) supplies the limits. Failed refits and rejected returned estimates are never removed or replaced: their unknown errors are placed at both extremes to enclose the empirical limits for every completion. Limits can be zero, infinite or unbounded. availability and expected_tail_draws describe this uncertainty. Bonferroni covers only the requested finite family, assuming adequate marginal bootstrap approximations; it is not a coverage guarantee.

For an LRT, test reports (1 + exceedances)/(nsim + 1). If any replicate test is unresolved, PValue is missing and PValueLower/PValueUpper bound all completions. Monte Carlo binomial bounds and resolution are also retained. This is fitted-model calibration, not an exact finite-sample test or a remedy for misspecification, dependence between Persons, or informative assignment. Increasing nsim improves simulation precision, not the underlying model. A profile-likelihood method is not used by this function.

Value

A serializable result containing source, optional null_fit, all trials with errors/warnings and seeds, draws on the free parameter scale, and settings. refit_draws additionally retains returned alternative-model parameter vectors before eligibility checks; rejected rows are diagnostic values, never inputs to confint(). trials identifies the last stage and whether each model refit returned. checks and source_checks retain category counts/states, numerical status and the information diagnostics already computed; missing fields mean not recorded, not a successful check. No extra information calculation is performed for recording. Older saved results may lack these fields. LRT results also have comparison and test. Printing or changing the interval level never refits. Failed replicates remain present. summary() returns planned/available/unresolved counts, settings and either the saved likelihood-ratio test or default slope intervals, plus retained numerical cautions. It does not refit or discard unresolved trials. If a saved analysis records selected refit updates in settings$recheck, printing and slope intervals retain a caution that it is not a complete rerun under one procedure. Sampling tables preserve that record and any settings$diagnostic_checks_scope. Missing history means not recorded; it does not certify that all draws used the currently installed estimator. After changing estimation or acceptance rules, run a separate complete bootstrap to evaluate the changed procedure; retain the original result.

Eligibility and unresolved replicates

A returned fit is not the same as an accepted replicate. For slope bootstrapping, a singleton score category (observed once) is a warning, not by itself a reason to exclude the point estimate. This applies to the source fit and refits only when a fresh category audit confirms all score categories are observed in every step scope, with no unsupported step coordinates or category contrasts. Model identity, numerical convergence, reevaluated likelihood/gradient and unregularized joint-information checks must still pass. The saved fit's readiness is not changed.

Basic bootstrap quantiles do not themselves require a variance estimate for every replicate. The retained information check is a conservative solution-quality restriction, not a requirement of the basic formula. Singular information or unstable numerical refinement remains unresolved, with returned estimates saved for diagnosis. Positive ill-conditioned information can be used with a caution after refinement, unregularized inversion and scaled-gradient checks pass. This is distinct from a random-effect variance estimated at zero; neither proves optimization failed by itself. Wald intervals, information-criterion comparisons and LRT eligibility share the numerical information review but retain their other checks.

checks separates BootstrapEligible from WaldEligible and retains BootstrapCaution. Accepted singleton or weak-information cases produce an aggregate warning; cautions also follow confint(), its printed output, default plot subtitle, APA tables and reports. A custom subtitle (including NULL) overrides the plot text, not the saved diagnostics. source_checks records the source decision. Use apa_table(result, which = "checks") to inspect refits. InformationRefinementVerified, InformationRelativeChange, InformationInverseResidual and InformationScaledGradient record a numerical refinement when needed; missing values mean it was not run. Do not drop unresolved replicates or substitute diagnostic refit_draws for accepted draws. A larger nsim does not remove an unresolved estimation problem or guarantee coverage.

PopulationSD, MinimumStandardizedSlope and MaximumStandardizedSlope describe the returned optimizer solution. Standardized slopes multiply relative slopes by the fitted population SD (one for a fixed standard-normal population). With covariates, this is the residual population SD. These diagnostics help distinguish scale changes from a slope approaching zero; they do not certify a boundary solution or label rater quality. A missing field means the required estimate was not retained. A near-zero slope, an unused score category and an ill-conditioned information matrix can have different consequences for different parameters: stable slopes do not establish finite thresholds or valid Wald intervals. Use checks to identify cases needing further numerical review before interpreting an unresolved case. APA check tables display scale, slope, gradient and information diagnostics with significant digits so that small positive values are not rounded to zero. The original checks fields remain unrounded numeric values.

References

Chalmers, R. P. (2012). mirt: A Multidimensional Item Response Theory Package for the R Environment. Journal of Statistical Software, 48(6), 1–29. doi:10.18637/jss.v048.i06.

Davison, A. C., and Hinkley, D. V. (1997). Bootstrap Methods and Their Application. Cambridge University Press, Chapter 5.

See Also

confint.mfrm_fit(), compare_mfrm(), mml_quadrature_sensitivity()

Examples

# After fitting compatible MML models:
# boot <- bootstrap_mfrm_gpcm(gpcm_fit, nsim = 499, seed = 92401)
# confint(boot, scale = "standardized")
# test <- bootstrap_mfrm_gpcm(gpcm_fit, nsim = 499, seed = 92402,
#                             null_fit = pcm_fit)
# test$test

Description

Build APA text outputs from model results

Usage

build_apa_outputs(
  fit,
  diagnostics,
  bias_results = NULL,
  context = list(),
  whexact = FALSE
)

Arguments

Details

context is an optional named list for narrative customization. Frequently used fields include:

• assessment, setting, scale_desc

• rater_training, raters_per_response

• rater_facet (used for targeted reliability note text)

• line_width (optional text wrapping width for report_text; default = 92)

Output text includes residual-PCA screening commentary if PCA diagnostics are available in diagnostics.

For GPCM, this helper returns a caveated partial reporting bundle over supported diagnostics, direct tables, and plots. It also includes a gpcm_boundary table. Treat the output as slope-aware sensitivity-reporting text, not FACETS score-side equivalence, automatic operational scoring, or design-forecasting evidence.

By default, report_text includes:

• model/data design summary (N, facet counts, scale range)

• optimization/convergence metrics (Converged, Iterations, LogLik) plus the eligible canonical MML AIC/Person-BIC/SABIC panel or its explicit ineligibility and legacy-descriptive status; q<31 MML criteria are explicitly labelled screening/review-only rather than automatic ranking

• anchor/constraint summary (noncenter_facet, centering convention, anchored levels, group anchors, dummy facets)

• latent-regression population-model wording when fit has an active population_formula

• category/threshold diagnostics (including disordered-step details when present)

• overall fit, misfit count, and top misfit levels

• facet reliability/separation, residual PCA summary, and bias-screen counts

• stored inter-rater agreement, its matched-opportunity denominator and unavailable outcomes, separately from separation reliability

Agreement wording uses the rater facet in the stored diagnostics summary; context$rater_facet does not recalculate or relabel that summary. The model-expected baseline uses fitted category probabilities with calibration fixed, not marginal-frequency chance agreement. Adjacent agreement includes exact matches and differences of at most 1 on the compared score scale. Missing expected values remain explicit, including the stored number of opportunities without valid probabilities. An absent summary does not establish agreement. See interrater_agreement_table() for the matching, repeated-rating and denominator conventions.

Value

An object of class mfrm_apa_outputs with:

• report_text: APA-oriented Method/Results draft prose

• decision: plain-language source-fit interpretation plus the separate precision-contract decision; FormalInference is "Yes" only when both fit readiness and contract$precision$supports_formal_inference support it

• fit_readiness, fit_readiness_components, and fit_readiness_parameters: exact source-fit readiness provenance

• table_figure_notes: consolidated draft notes for tables/visuals

• table_figure_captions: draft caption candidates without figure numbering

• section_map: package-native section table for manuscript assembly

• contract: structured APA reporting contract used for downstream checks