eval=FALSE to computationally
intensive code chunks in vignettes (GridSearch, IRM, LDLRA, LDLRA_PBIL,
LDB, BINET, ordinal/nominal Biclustering, ordinal/rated LRA). The
Japanese guide (guide-ja) now uses eval=FALSE
for all model estimation chunks to avoid duplicating computation from
English vignettes. Full rendered output is available on the pkgdown site.log(n + 1)
but should be log(n) + 1 per Bozdogan (1987, Psychometrika,
52(3), p.358, Proposition 2, Eq.44). The original Mathematica
implementation had this error (Log[nobs + 1]), and the R
port inherited it. Both the R version
(R/00_ModelFitModule.R) and the Mathematica version
(develop/mtmk15forVer13/mod/Module_ModelFit.nb) have been
corrected. The numerical difference is approximately 1 (constant), but
the corrected formula now matches the published definition:
CAIC(k) = -2 log L + k * (log(n) + 1). This affects all
models that compute fit indices: IRT, LCA, LRA (binary/ordinal/rated),
Biclustering (binary/ordinal/nominal), IRM, BNM, LDLRA, LDB, BINET, and
GRM. All 873 tests pass with the corrected formula.\donttest to
\dontrun: GRM’s multi-panel plot examples caused
“invalid graphics state” errors in non-interactive environments
(pkgdown, CI). Changed to \dontrun to prevent build
failures.CA
parameter: When CA (correct answer vector) was
provided but response.type was not explicitly specified,
dataFormat() incorrectly classified the data as
"ordinal" instead of "rated". This caused
$U (binary scoring matrix) to be NULL. The
auto-detection now checks for CA before ordinal/nominal
detection: binary → rated (if CA provided) → ordinal → nominal.id parameter not working for column
1: Changed id default from 1 to
NULL. Previously, id = 1 (both default and
explicit) triggered auto-detection heuristics, making it impossible to
explicitly specify the first column as the ID column. Now
id = NULL triggers auto-detection, and any numeric value
(including 1) forces that column to be used as ID.1:N) as ID regardless of
whether they also look like valid response values. Previously,
1:10 was misclassified as response data because
looks_like_response_data() returned TRUE.U
matrix: When the response matrix contained missing values
(-1), the binary scoring matrix U incorrectly
scored them as 0 (incorrect) instead of -1
(missing). Now U[i,j] = -1 when
Q[i,j] = -1.drop=FALSE missing in item
exclusion: When items were excluded due to invalid variance
(e.g., containing Inf), the column subsetting
response.matrix[, mask] dropped matrix dimensions if only
one item remained, causing a crash. Added
drop = FALSE.dataFormat() now reports via message() when it
detects problematic data:
g_list
/ adj_list Input Path Fixg_csv variable undefined error when using
g_list or adj_list input:
BINET() crashed with an undefined variable error at the
graph construction step when the DAG was specified via
g_list or adj_list parameters. The internal
variable g_csv (used to build the integrated graph object
all_g) was only defined in the adj_file code
path. Added logic to reconstruct g_csv from
adj_list for the g_list and
adj_list input paths, ensuring all_g is
correctly built with Field edge attributes regardless of input
method.g_list / adj_list length
validation: The length check for g_list and
adj_list incorrectly compared against ncls
(number of classes) instead of nfld (number of fields). In
BINET, each element of adj_list represents the DAG
structure at a specific field, so the list length should equal
nfld. Also fixed g_list type check from
g_list[[1]] (always checking the first element) to
g_list[[j]] (checking each element).... to
Biclustering_IRM.binary() signature: The S3
generic Biclustering_IRM(U, ...) requires all methods to
include ... in their formal arguments. The
.binary method was missing it, causing an R CMD check
WARNING. Added ... to match the generic and the
.nominal/.ordinal methods.msg Field Fixmsg field to correctly distinguish
Rank vs Class models: Binary IRM and Ordinal IRM are
Ranklustering models (ordered latent classes), so
msg = "Rank" is correct. Nominal IRM has no Ranklustering
concept (unordered latent classes), so msg = "Class" is
correct. Previously all three methods had inconsistent values; now
Biclustering_IRM.binary and
Biclustering_IRM.ordinal return msg = "Rank",
while Biclustering_IRM.nominal returns
msg = "Class". Also updated internal column names of
cls01 and FRP from "Class" to
"Rank" for binary/ordinal IRM. The shared Gibbs core
(irm_gibbs_core()) output column names were also
updated.LRD AliasLRD (Latent Rank Distribution) alias to
Biclustering_IRM.binary() return value: The binary
IRM was the only Biclustering model that did not include
LRD in its return value (it only had LCD).
Other Biclustering models (binary Biclustering, nominal IRM, ordinal
IRM) all return both LRD and LCD. Added
LRD = clsdist as an alias for LCD to ensure
consistency across all Biclustering-family models. This improves
interoperability with downstream packages (ggExametrika) that look up
LRD when msg == "Rank".Biclustering_IRM() seed default back
to 123: Ensures reproducibility by default.t(apply()) Dimension Drop Fixt(apply(log_S, 1, irm_log_to_prob))
dimension drop in Nominal/Ordinal IRM: When the Gibbs sampler
converged to ncls=1 or nfld=1,
apply() on a single-column matrix returned a vector instead
of a matrix, causing t() to produce a 1-row matrix instead
of an N-row matrix. This led to incorrect class/field membership
assignment. Replaced all 6 instances of
t(apply(mat, 1, fun)) with explicit for-loops in both
Biclustering_IRM.nominal() (3 locations: EM E-step, final
class membership, final field membership) and
Biclustering_IRM.ordinal() (3 locations: same). Consistent
with the project’s apply() caution policy.Biclustering_IRM() Gibbs sampler: Changed
exp(ptab - min(ptab)) to exp(ptab - max(ptab))
for numerical stability. The previous implementation subtracted the
minimum log-probability, which could cause overflow
(exp(large positive) = Inf) when the range of
log-probabilities was large, resulting in NaN after
normalization. Subtracting the maximum ensures the largest value becomes
exp(0) = 1 and all others underflow harmlessly to near-zero
values.Biclustering.nominal() using stale
log-likelihood in model fit indices: When the EM algorithm
exited due to log-likelihood decrease, BCRM was correctly
reverted to the previous iteration’s value, but
test_log_lik was not. The model fit section recalculated
the correct log-likelihood as testell, but
chi_A, model_log_like, log_lik,
and LogLik all referenced the stale
test_log_lik instead of testell. Also removed
a duplicated model fit code block (copy-paste error).J20S400 response type from
nominal to binary: The
J20S400.rda dataset was incorrectly stored with
response.type = "nominal" despite being a binary (0/1)
dataset with -99 as missing values. Missing values were not properly
handled (Z was all 1s, Q contained raw -99 values). Regenerated from the
original CSV source (develop/sampleData/J20S400.csv) with
dataFormat(..., na = -99), now correctly typed as
binary with 86 missing values properly masked in Z.Biclustering_IRM.nominal() for
nominal/polytomous data: Extends the Infinite Relational Model
(IRM) from binary to nominal scale data using a Dirichlet-Multinomial
collapsed Gibbs sampler. The Chinese Restaurant Process (CRP)
automatically determines the optimal number of classes and fields. After
the Gibbs sampling phase, small classes are consolidated and refined
with an EM algorithm. The Dirichlet prior concentration parameter
alpha controls smoothing of category probabilities.Biclustering_IRM.ordinal() for
ordinal/polytomous data: Extends IRM to ordinal scale data.
Shares the same Dirichlet-Multinomial collapsed Gibbs sampler as nominal
IRM (Phase 1), then applies ordinal-specific EM refinement (Phase 2)
with cumulative normalization to enforce monotonic category boundaries.
The mic parameter (default TRUE) enforces
monotone increasing class ordering by expected score sum. Reports BFRP
(Bicluster Field Reference Profile), FRPIndex, TRP, and Strongly/Weakly
Ordinal Alignment Conditions (SOACflg/WOACflg).Biclustering_IRM() is now an S3
generic: Dispatches to Biclustering_IRM.binary
(existing binary IRM), Biclustering_IRM.nominal (new), and
Biclustering_IRM.ordinal (new). Raw data is automatically
formatted and dispatched based on response.type.irm_gibbs_core() (in R/00_IRM_Gibbs_CORE.R),
used by both nominal and ordinal IRM. Helper functions
(irm_calc_Ufcq, irm_lmvbeta,
irm_log_to_prob, irm_bic_calc, etc.) are also
shared.U_fcq, computing only the contribution of the target
student/item rather than recalculating the full array at each step.LCA() and LRA() now support a
conf parameter for confirmatory analysis: The
conf argument accepts an IRP matrix (ncls/nrank x
testlength) where non-NA values are held fixed throughout EM estimation
and NA values are freely estimated. This enables test equating scenarios
where anchor items retain their known IRPs while new items are
calibrated against them.conf
has column names, items are matched by label rather than by position.
This allows conf to contain a subset of items (anchor items
only) or items in a different order than the data. Items in the data but
not in conf are automatically set to freely estimated.
Unknown labels in conf produce an informative error.somclus()
internal function: The Self-Organizing Maps estimation code
previously inlined in LRA.binary() (~100 lines) has been
extracted into a standalone internal function somclus() in
00_EMclus.R. This parallels emclus() (GTM/EM)
and returns the same structure. Also fixed a pre-existing typo
(h_cout → h_count) and another
(oldsBIC → oldBIC).tidyverse
and readxl for reading Excel-based Mathematica reference
data. Replaced with 24 modern test files using base R
read.csv() and test_path() to load CSV
fixtures from
tests/testthat/fixtures/mathematica_reference/. The new
test suite has zero external package dependencies beyond
testthat.readxl/tidyverse from test
dependencies: DESCRIPTION Suggests no longer
requires any packages beyond knitr, rmarkdown,
and testthat.test-irm-nominal.R): 18 test blocks for
Biclustering_IRM.nominal() using J20S600 data — basic
execution, dimensions, FRP validity, membership, fit indices, backward
compatibility, seed reproducibility, alpha validation.test-irm-ordinal.R): 25 test blocks for
Biclustering_IRM.ordinal() using J35S500 data — basic
execution, dimensions, FRP validity, expected scores, TRP, BFRP,
FRPIndex, SOAC/WOAC flags, mic parameter, fit indices, backward
compatibility, seed reproducibility, alpha validation.stem(nrs(dataFormat(J15S500))) to demonstrate
score distribution visualization.^tests$
and ^inst$ entries that were incorrectly excluding tests
and vignettes from the built package.index
Parameter FixesGridSearch() now
accepts common aliases for fit indices. "loglik",
"log_lik", "LogLik", and "LL" are
mapped to "model_log_like". "Chi_sq" and
"chi_sq" are mapped to "model_Chi_sq".
Previously, using these aliases caused silent NULL
extraction and eventual errors.GridSearch(), before the
computationally expensive grid search loop runs. The error message lists
all valid options including available aliases.model_log_like is now correctly treated as a maximization
target (larger log-likelihood = better fit). Previously, it was
incorrectly placed in the minimization group, which would have selected
the worst-fitting model.tolerance = 1e-2 in
expect_equal) to absolute difference comparison (threshold
0.005). Q3 residual correlations include values near zero where relative
tolerance comparisons are unreliable.layout and
plot.new from graphics to NAMESPACE. These
functions are used by the legend strip layout helpers for polytomous
Biclustering plots.apply(U$Q, 2, unique) returning matrix
instead of list: When all items have the same number of
response categories (e.g., all 5-point Likert), apply()
returns a matrix rather than a list. The subsequent
lapply() then iterates over individual elements instead of
per-column vectors, causing ncat to be a vector of 1s and
crashing the algorithm. Replaced with
lapply(seq_len(nitems), function(j) sort(unique(U$Q[, j])))
which always returns a proper list. Same fix applied to
catfreq999 computation. Both LRA.ordinal and
LRA.rated were affected.apply(tmp$Q, 2, table) returning matrix
instead of list: Same class of bug as the LRA.ordinal/LRA.rated
fix above. In GRM(), the null model log-likelihood
computation used apply(tmp$Q, 2, table) to compute per-item
category frequencies. When all items have the same number of response
categories (e.g., all 5-point Likert), apply() returns a
matrix instead of a list, causing response_list[[j]] to
extract a single number rather than the full frequency table. This
produced incorrect null_log_like values and consequently
wrong chi-square statistics, RMSEA, TLI, and CFI in
ItemFitIndices. Replaced with
lapply(seq_len(nitems), function(j) table(tmp$Q[, j]))
which always returns a proper list.LRA.ordinal() now raises an informative error when items
have different numbers of response categories (e.g., some items with 3
categories and others with 5). The internal matrix algebra uses
fixed-stride indexing that assumes uniform category counts. The error
message suggests alternatives (LRA.rated,
Biclustering.ordinal) that support mixed category counts
via list-based designs.plot(lca_result, type = "FRP") passed validation but failed
at runtime because the $FRP field does not exist in LCA/LRA
return values. Now properly rejected with an informative error message
at the validation stage.stat parameter supporting
"mean" (default), "median", and
"mode".style parameter supporting
"line" (default) and "bar" (stacked bar
chart).stat
parameter.FRPIndex (Field
Reference Profile indices) including location parameters (Alpha, Beta),
slope parameters (A, B), and monotonicity indices (Gamma, C).plot.exametrika()stat: Controls the summary statistic for FRP and RRV
plots ("mean", "median",
"mode").style: Controls the display style for FCRP plots
("line", "bar").#808080) to distinguish from white (incorrect) and black
(correct). Polytomous data uses black (#000000) to
distinguish from the category color palette.print() now
displays FRPIndex (Field Reference Profile Indices) with a note that the
values are based on normalized expected scores
(E[score]-1)/(maxQ-1).?Biclustering, including detailed definitions and
polytomous adaptation logic.Systematic unification of return value structures across all analysis functions for consistency and interoperability.
n_class (retaining
Nclass for backward compatibility)n_rank,
n_field (retaining Nrank,
Nfield)n_class,
n_field (retaining Nclass,
Nfield)n_class,
n_field, n_cycle (retaining
Nclass, Nfield, N_Cycle)n_class,
n_field, n_cycle (retaining
Nclass, Nfield, N_Cycle)n_cycle,
N_Cycle (retaining em_cycle,
EM_Cycle as IRM-specific aliases)log_lik Added to All Functionslog_lik
at the top level of the return object, matching
TestFitIndices$model_log_like.ModelFit class, matching the format used by IRT/LCA/LRA and
other functions. Previously these used bare
calcFitIndices() output (9 fields, no class). Added fields:
model_log_like, bench_log_like,
null_log_like, model_Chi_sq,
null_Chi_sq, model_df,
null_df.ModelFit class.TestFitIndices added as primary
name for multigroup fit indices (previously only
MG_FitIndices). MG_FitIndices retained as
backward-compatible alias. SM_FitIndices (saturated model)
remains unchanged.Estimate
column (most probable class assignment) to the Students matrix.Estimate
column (most probable class assignment) to the Students matrix.FRPIndex
(Field Reference Profile indices) for consistency with
Biclustering.binary and LDB.TRP from
matrix(1×ncls) to numeric vector, consistent
with all other functions.class = c("exametrika", "GridSearch") to return value for
method dispatch support.msg field assignment: Fixed
msg <- "Class" to msg = "Class" inside
structure() call (line 150 of 05_LCA.R). The
<- operator was being interpreted as a standalone
assignment rather than a named list element, causing the
msg field name to be empty.plot(r, type="RMP", students=1)). Added
drop = FALSE to prevent matrix-to-vector coercion when
extracting a single row from the Students matrix.TestFitIndices
log-likelihood fields: Fixed null_log_like which
incorrectly stored the saturated model log-likelihood
(log_lik_satu) instead of the null model log-likelihood
(sum(null_itemll)). Added missing
bench_log_like field containing the saturated model
log-likelihood (sum(satu_itemll2)). Note: the
chi-square-based fit indices (NFI, CFI, TLI, RMSEA, AIC, BIC, etc.) were
always computed correctly; only the stored log-likelihood labels were
affected.TestFitIndices and ItemFitIndices to
the standard 16-field structure with ModelFit class
(c("exametrika", "ModelFit")), matching all other analysis
functions. Added model_log_like,
bench_log_like, null_log_like to
ItemFitIndices. Removed ScoreRankCorr /
RankQuantCorr from TestFitIndices (already
available at the top level of the return object).test-12OLR.R and test-13NLR.R to handle the
new ModelFit class (add unclass() before
as.data.frame()) and adjusted column/index references to
match the unified 16-field structure.layout() with a
thin dedicated row (height ratio 0.2) to reduce visual clutter in data
panels. Added setup_legend_layout() and
draw_legend_strip() internal helper functions.P(Q>=1)=1.0 reference line at the top of FCBR plots for
visual completeness.dataFormat() now correctly identifies factor-type ID
columns before converting factors to numeric. Previously, the
factor-to-numeric conversion occurred before ID detection, causing
factor ID columns with many levels (>=20) to trigger a “Too many
categories” error instead of being recognized as IDs.is_response_data()) that was defined but never
called.obj$Q instead of
obj$U for test length calculationpmax(Ufcq_prior, 1e-10) to prevent division by zero and NaN
errorsThis release improves naming consistency across the package while maintaining full backward compatibility through a deprecation path.
All analysis functions now return results with snake_case field names for better consistency:
n_class - Number of latent classes (replaces
Nclass)n_field - Number of latent fields (replaces
Nfield)n_rank - Number of latent ranks (replaces
Nrank)n_cycle - Number of EM iterations (replaces
N_Cycle)log_lik - Log-likelihood value (replaces
LogLik)The old field names continue to work for backward compatibility:
Nclass, Nfield, Nrank,
N_Cycle, LogLik - Deprecated but
functional# Old code (deprecated but still works)
result <- LCA(data, ncls = 3)
n_classes <- result$Nclass
# New code (recommended)
result <- LCA(data, ncls = 3)
n_classes <- result$n_classProgress messages now display properly in R Markdown documents:
verbose parameter
(default: TRUE)
ncls = 2: nfld=2 nfld=3 nfld=4 ...verbose = FALSE to suppress all progress
messagesiter 1: match=0 nfld=15 ncls=30Adjusting classes: BIC=-99592.5 ncls=21 (min size < 20)verbose = TRUE to
verbose = FALSE (consistent with binary/rated
versions)Saturation Model - iter 1: log_lik=-1.234567Restricted Model - iter 1: log_lik=-0.987654verbose = TRUE to
verbose = FALSE (consistent with binary/ordinal
versions)All progress messages now use proper line breaks instead of carriage returns, ensuring clean output in R Markdown/knitr documents and web documentation.
sort(unique(...)) to ensure consistent ordering: 0
(white/incorrect) → 1 (black/correct)testEll →
test_log_lik)log_lik
terminologyconverge variable to all EM-based functions to
indicate algorithm convergence status
Added implementation of latent rank model for ordinal scale data
New function LRA() now supports ordinal response data
Bug fixes and improvements
Standardized terminology: unified the usage of “class” and “rank” throughout the package
Various minor bug fixes