Factor rotation is a crucial step in interpreting the results of exploratory factor analysis. Several rotation methods have been developed for simple structure solutions, but their extensions to bi-factor analysis are often not well established. In this article, we propose a mathematical framework that incorporates customized factor structure as a regularization to produce the optimal orthogonal or oblique rotation. We demonstrate the utility of the framework using examples of simple structure rotation and bi-factor rotation. Through detailed simulations, we show that the new method is accurate and robust in recovering the factor structures and latent correlations when bi-factor analysis is applied. The new method is applied to a test data and a Quality of Life survey data. Results show that our method can reveal bi-factor structures that are consistent with the theories.

Nonparametric item response theory (IRT) models consist of assumptions that restrict the joint item-score distribution. These assumptions imply stochastic ordering properties that allow ordering of respondents and items using the simple sum score and item mean score, respectively, and imply observable data properties that are useful for investigating model fit. In this paper, we investigate these properties for two-level nonparametric IRT. We introduce four two-level nonparametric IRT models. Two models pertain to respondents nested in groups: The MHM-1, useful for ordering respondents and groups, and the DMM-1, useful for ordering respondents, groups, and items. Two models pertain to groups rated by multiple respondents: The MHM-2, useful for ordering groups, and the DMM-2, useful for ordering groups and items. We define the model assumptions, derive implied stochastic ordering properties, and derive observable data properties that are useful for model fit investigation. Relations between models and properties are also presented.

This paper presents a model specification for group comparisons regarding a functional trend over time within a trial and learning across a series of trials in intensive binary longitudinal eye-tracking data. The functional trend and learning effects are modeled using by-variable smooth functions. This model specification is formulated as a generalized additive mixed model, which allowed for the use of the freely available mgcv package (Wood in Package ‘mgcv.’ https://cran.r-project.org/web/packages/mgcv/mgcv.pdf, 2023) in R. The model specification was applied to intensive binary longitudinal eye-tracking data, where the questions of interest concern differences between individuals with and without brain injury in their real-time language comprehension and how this affects their learning over time. The results of the simulation study show that the model parameters are recovered well and the by-variable smooth functions are adequately predicted in the same condition as those found in the application.

Assessing fit in common factor models solely through the lens of mean and covariance structures, as is commonly done with conventional goodness-of-fit (GOF) assessments, may overlook critical aspects of misfit, potentially leading to misleading conclusions. To achieve more flexible fit assessment, we extend the theory of generalized residuals (Haberman & Sinharay, 2013), originally developed for models with categorical data, to encompass more general measurement models. Within this extended framework, we propose several fit test statistics designed to evaluate various parametric assumptions involved in common factor models. The examples include assessing the distributional assumptions of latent variables and the functional form assumptions of individual manifest variables. The performance of the proposed statistics is examined through simulation studies and an empirical data analysis. Our findings suggest that generalized residuals are promising tools for detecting misfit in measurement models, often masked when assessed by conventional GOF testing methods.

Advances in deep learning and representation learning have transformed item factor analysis (IFA) in the item response theory (IRT) literature by enabling more efficient and accurate parameter estimation. Variational autoencoders (VAEs) are widely used to model high-dimensional latent variables in this context, but the limited expressiveness of their inference networks can still hinder performance. We introduce adversarial variational Bayes (AVB) and an importance-weighted extension (IWAVB) as more flexible inference algorithms for IFA. By combining VAEs with generative adversarial networks (GANs), AVB uses an auxiliary discriminator network to frame estimation as a two-player game and removes the restrictive standard normal assumption on the latent variables. Theoretically, AVB and IWAVB can achieve likelihoods that match or exceed those of VAEs and importance-weighted autoencoders (IWAEs). In exploratory analyses of empirical data, IWAVB attained higher likelihoods than IWAE, indicating greater expressiveness. In confirmatory simulations, IWAVB achieved comparable mean-square error in parameter recovery while consistently yielding higher likelihoods, and it clearly outperformed IWAE when the latent distribution was multimodal. These findings suggest that IWAVB can scale IFA to complex, large-scale, and potentially multimodal settings, supporting closer integration of psychometrics with modern multimodal data analysis.

Survey questionnaires are commonly used by psychologists and social scientists to measure various latent traits of study subjects. Various causal inference methods such as the potential outcome framework and structural equation models have been used to infer causal effects. However, the majority of these methods assume the knowledge of true causal structure, which is unknown for many applications in psychological and social sciences. This calls for alternative causal approaches for analyzing such questionnaire data. Bayesian networks are a promising option as they do not require causal structure to be known a priori but learn it objectively from data. Although we have seen some recent successes of using Bayesian networks to discover causality for psychological questionnaire data, their techniques tend to suffer from causal non-identifiability with observational data. In this paper, we propose the use of a state-of-the-art Bayesian network that is proven to be fully identifiable for observational ordinal data. We develop a causal structure learning algorithm based on an asymptotically justified BIC score function, a hill-climbing search strategy, and the bootstrapping technique, which is able to not only identify a unique causal structure but also quantify the associated uncertainty. Using simulation studies, we demonstrate the power of the proposed learning algorithm by comparing it with alternative Bayesian network methods. For illustration, we consider a dataset from a psychological study of the functional relationships among the symptoms of obsessive-compulsive disorder and depression. Without any prior knowledge, the proposed algorithm reveals some plausible causal relationships. This paper is accompanied by a user-friendly open-source R package OrdCD on CRAN.

Fitting propensity (FP) analysis quantifies model complexity but has been impeded in item response theory (IRT) due to the computational infeasibility of uniformly and randomly sampling multinomial item response patterns under a full-information approach. We adopt a limited-information (LI) approach, wherein we generate data only up to the lower-order margins of the complete item response patterns. We present an algorithm that builds upon classical work on sampling contingency tables with fixed margins by implementing a Sequential Importance Sampling algorithm to Quickly and Uniformly Obtain Contingency tables (SISQUOC). Theoretical justification and comprehensive validation demonstrate the effectiveness of the SISQUOC algorithm for IRT and offer insights into sampling from the complete data space defined by the lower-order margins. We highlight the efficiency and simplicity of the LI approach for generating large and uniformly random datasets of dichotomous and polytomous items. We further present an iterative proportional fitting procedure to reconstruct joint multinomial probabilities after LI-based data generation, facilitating FP evaluation using traditional estimation strategies. We illustrate the proposed approach by examining the FP of the graded response model and generalized partial credit model, with results suggesting that their functional forms express similar degrees of configural complexity.

Rating procedure is crucial in many applied fields (e.g., educational, clinical, emergency). In these contexts, a rater (e.g., teacher, doctor) scores a subject (e.g., student, doctor) on a rating scale. Given raters’ variability, several statistical methods have been proposed for assessing and improving the quality of ratings. The analysis and the estimate of the Intraclass Correlation Coefficient (ICC) are major concerns in such cases. As evidenced by the literature, ICC might differ across different subgroups of raters and might be affected by contextual factors and subject heterogeneity. Model estimation in the presence of heterogeneity has been one of the recent challenges in this research line. Consequently, several methods have been proposed to address this issue under a parametric multilevel modelling framework, in which strong distributional assumptions are made. We propose a more flexible model under the Bayesian nonparametric (BNP) framework, in which most of those assumptions are relaxed. By eliciting hierarchical discrete nonparametric priors, the model accommodates clusters among raters and subjects, naturally accounts for heterogeneity, and improves estimates’ accuracy. We propose a general BNP heteroscedastic framework to analyze continuous and coarse rating data and possible latent differences among subjects and raters. The estimated densities are used to make inferences about the rating process and the quality of the ratings. By exploiting a stick-breaking representation of the discrete nonparametric priors, a general class of ICC indices might be derived for these models. Our method allows us to independently identify latent similarities between subjects and raters and can be applied in precise education to improve personalized teaching programs or interventions. Theoretical results about the ICC are provided together with computational strategies. Simulations and a real-world application are presented, and possible future directions are discussed.

To assess country-level progress toward these educational goals it is important to monitor trends in educational outcomes over time. The purpose of this article is to demonstrate how optimally predictive growth models can be constructed to monitor the pace of progress at which countries are moving toward (or way from) the education sustainable development goals as specified by the United Nations. A number of growth curve models can be specified to estimate the pace of progress, however, choosing one model and using it for predictive purposes assumes that the chosen model is the one that generated the data, and this choice runs the risk of “over-confident inferences and decisions that are more risky than one thinks they are” (Hoeting et al., 1999). To mitigate this problem, we adapt and apply Bayesian stacking to form mixtures of predictive distributions from an ensemble of individual models specified to predict country-level pace of progress. We demonstrate Bayesian stacking using country-level data from the Program on International Student Assessment. Our results show that Bayesian stacking yields better predictive accuracy than any single model as measured by the Kullback–Leibler divergence. Issues of Bayesian model identification and estimation for growth models are also discussed.

Numerous applications in psychological research require that a data set is partitioned via the inverse of a clustering criterion. This anticlustering seeks for high similarity between groups (maximum diversity) or high pairwise dissimilarity within groups (maximum dispersion). Brusco et al. (2020) proposed a bicriterion heuristic (BILS) that simultaneously seeks for maximum diversity and dispersion, introducing the bicriterion approach for anticlustering. We investigate if the bicriterion approach can be improved using exact algorithms that guarantee globally optimal criterion values. Despite the theoretical computational intractability of anticlustering, we present a new exact algorithm for maximum dispersion that scales to quite large data sets ($N = 1,000$). However, a fully exact bicriterion approach was only feasible for small data sets (about $N = 30$). We therefore developed hybrid approaches that maintain optimal dispersion but use heuristics to maximize diversity on top of it. In a simulation study and an example application, we compared several hybrid approaches. An adaptation of BILS that initiates each iteration with a partition having optimal dispersion (BILS-Hybrid-All) performed best across a variety of data inputs. All of the methods developed here as well as the original BILS algorithm are available via the free and open-source R package anticlust.

Multivariate analysis of psychological variables using graphical models has become a standard analysis in the psychometric literature. Most cross-sectional measures are either binary or ordinal, and the methodology for inferring the structure of networks of binary and ordinal variables is developing rapidly. In practice, however, research questions often focus on whether and how networks differ between observed groups. While Bayes factor methods for inferring network structure are well established, a similar methodology for assessing group differences in networks of binary or ordinal variables is currently lacking. In this article, we extend the Bayesian framework for the analysis of ordinal Markov random fields, a network model for binary and ordinal variables, and develop Bayes factor tests for assessing parameter differences in the networks of two independent groups. The proposed methods are implemented in the R package bgms, and we use numerical illustrations to show that the implemented methods work correctly and how well the methods work compared to existing methods in situations resembling empirical research.

This article introduces item response models for rating relational data. The relational data are obtained via ratings of senders and receivers in a directed network. The proposed models allow comparisons of senders and receivers on a one-dimensional latent scale while accounting for unobserved homophilic relationships. We show that the approach effectively captures reciprocity and clustering phenomena in the relational data. We estimate model parameters using a Bayesian specification and utilize Markov Chain Monte Carlo methods to approximate the full conditional posterior distributions. Simulation studies demonstrate that model parameters can be recovered satisfactorily even when the dimensionality of the network is small. We also present an extensive empirical application to illustrate the usefulness of the proposed models for complete and incomplete networks.

The latent Markov model (LMM) has been increasingly used to analyze log data from computer-interactive assessments. An important consideration in applying the LMM to assessment data is measurement effects of items. In educational and psychological assessment, items exhibit distinct psychometric qualities and induce systematic variance to assessment outcome data. The current development in LMM, however, assumes that items have uniform effects and do not contribute to the variance of measurement outcomes. In this study, we propose a refinement of LMM that relaxes the measurement invariance constraint and examine empirical performance of the new framework through numerical experimentation. We modify the LMM for noninvariant measurements and refine the inferential scheme to accommodate the event-specific measurement effects. Numerical experiments are conducted to validate the proposed inference methods and evaluate the performance of the new framework. Results suggest that the proposed inferential scheme performs adequately well in retrieving the model parameters and state profiles. The new LMM framework demonstrated reliable and stable performance in modeling latent processes while appropriately accounting for items’ measurement effects. Compared with the traditional scheme, the refined framework demonstrated greater relevance to real assessment data and yielded more robust inference results when the model was ill-specified. The findings from the empirical evaluations suggest that the new framework has potential for serving large-scale assessment data that exhibit distinct measurement effects.

Extended redundancy analysis (ERA) is a statistical approach to component-based multivariate regression modeling that explores interrelationships among multiple sets of while incorporating regression with a data-reduction technique. The extant models that utilize ERA have assumed the outcome variables with the same data type. Also, ERA models focused on estimating direct pathways only without explicitly addressing mediation effects. In this paper, ERA is extended to handle multiple mediators and mixed types of outcome variables by adopting a Bayesian framework, taking into account correlation structure among all of the outcome variables. The proposed method develops an algorithm that derives the joint posterior distribution of parameters using a Markov chain Monte Carlo algorithm. Simulations and an empirical dataset are provided to illustrate the usefulness of the proposed method.

For classroom teaching and learning, classifying students’ skills into more than two categories (e.g., no, basic, and advanced masteries) is more instructionally relevant. Such classifications require polytomous attributes, for which most existing cognitive diagnosis models (CDMs) are inapplicable. This paper proposes the saturated polytomous cognitive diagnosis model (sp-CDM), a general model that subsumes existing CDMs for polytomous attributes as special cases. The generalization is shown by mathematically illustrating the relationships between the proposed and existing CDMs. Moreover, algorithms to estimate the proposed model is proposed. A simulation study is conducted to evaluate the parameter recovery of the sp-CDM using the proposed estimation algorithms, as well as to illustrate the consequences of improperly fitting constrained or unnecessarily complex polytomous-attribute CDMs. A real-data example involving polytomous attributes is presented to demonstrate the practical utility of the proposed model.

Factor score indeterminacy is a characteristic property of factor analysis (FA) models. This research introduces a novel procedure, regression-based factor score exploration (RFE), which uniquely determines factor scores and simultaneously estimates other parameters of the FA model. RFE uniquely determines factor scores by minimizing a loss function that balances FA and multivariate regression, regulated by a tuning parameter. Theoretical aspects of RFE, including the uniqueness of factor scores, the relationship between observed and latent variables, and rotational indeterminacy, are examined. Additionally, clustering-based factor exploration (CFE) is presented as a variant of RFE, derived by generalizing the penalty term to enable the clustering of factor scores. It is demonstrated that CFE creates cluster structures more accurately than the existing method. A simulation study shows that the proposed procedures accurately recover true parameter matrices even in the presence of error-contaminated data, with lower computational demand compared to existing methods. Real data examples illustrate that the proposed procedures provide interpretable results, demonstrating high relevance to the factor scores obtained by existing methods.

With increasingly available computer-based or online assessments, researchers have shown keen interest in analyzing log data to improve our understanding of test takers’ problem-solving processes. In this article, we propose a multi-state survival model (MSM) to action sequence data from log files, focusing on modeling test takers’ reaction times between actions, in order to investigate which factors and how they influence test takers’ transition speed between actions. We specifically identify the key actions that differentiate correct and incorrect answers, compare transition probabilities between these groups, and analyze their distinct problem-solving patterns. Through simulation studies and sensitivity analyses, we evaluate the robustness of our proposed model. We demonstrate the proposed approach using problem-solving items from the Programme for the International Assessment of Adult Competencies (PIAAC).

We develop a fully nonparametric, easy-to-use, and powerful test for the missing completely at random (MCAR) assumption on the missingness mechanism of a dataset. The test compares distributions of different missing patterns on random projections in the variable space of the data. The distributional differences are measured with the Kullback-Leibler Divergence, using probability Random Forests (Malley et al., 2011). We thus refer to it as “Projected Kullback–Leibler MCAR” (PKLM) test. The use of random projections makes it applicable even if very few or no fully observed observations are available or if the number of dimensions is large. An efficient permutation approach guarantees the level for any finite sample size, resolving a major shortcoming of most other available tests. Moreover, the test can be used on both discrete and continuous data. We show empirically on a range of simulated data distributions and real datasets that our test has consistently high power and is able to avoid inflated type-I errors. Finally, we provide an R-package PKLMtest with an implementation of our test.

The basic local independence model (BLIM) is a probabilistic model developed in knowledge space theory (KST). Recently, Stefanutti, de Chiusole, et al. (2020, Psychometrika 85, 684–715) proposed the polytomous local independence model (PoLIM), which is an extension of the BLIM to items with more than two response alternatives (polytomous items). In a Commentary to this paper, Chiu et al. (2023, Psychometrika 88, 656–671) claimed that (i) the BLIM is just a deterministic input noisy AND-gate (DINA) model where every item has a single skill and, as a consequence of this, (ii) the “PoLIM is simply a paraphrase of a DINA model in cognitive diagnosis (CD) for polytomous items” (p. 656). This rejoinder shows that such statements are invalid and totally misleading. Its aim is to clarify the nature of the relationship between the BLIM and the DINA, as well as that between the PoLIM and the Polytomous DINA. It builds upon formal results by Heller, et al. (2015, Psychometrika 80(4), 995–1019) on the intimate relation between KST and CD notions, and shows that the BLIM/PoLIM may be conceived as marginal models for whole classes of CD models.

We propose an optimally regularized Bayesian estimator of multilevel latent variable models that aims to outperform traditional maximum likelihood (ML) estimation in mean squared error (MSE) performance. We focus on the between-group slope in a two-level model with a latent covariate. Our estimator combines prior information with data-driven insights for optimal parameter estimation. We present a “proof of concept” by computer simulations, involving varying numbers of groups, group sizes, and intraclass correlations (ICCs), which we conducted to compare the newly proposed estimator with ML. Additionally, we provide a step-by-step tutorial on applying the regularized Bayesian estimator to real-world data using our MultiLevelOptimalBayes package.

Encouragingly, our results show that our estimator offers improved MSE performance, especially in small samples with low ICCs. These findings suggest that the estimator can be an effective means for enhancing estimation accuracy.