Statistical Methods
Reference

Classical Test Theory

CTT characterizes items and tests using observed-score statistics from an \(N\times K\) response matrix (\(N\) respondents, \(K\) items). Items may be dichotomous (0/1) or polytomous (non-negative integers). All variance estimates use the sample \((n-1)\) denominator throughout. Missing data are handled pairwise (default) or listwise.

Item-Level Descriptive Statistics

Mean, SD, and Standard Error

Let \(\mathcal{V}_j\) be the set of respondents with non-missing responses to item \(j\) and \(n_j=|\mathcal{V}_j|\). Then

Relative Difficulty (p-value)

Items with \(p_j<0.15\) or \(p_j>0.85\) are flagged. For dichotomous items this equals the classical p-value.

Discrimination Indices

Adjusted Point-Biserial (Item-Rest Correlation)

Item \(j\) is correlated with the sum of the remaining \(K-1\) items, eliminating the inflation from item \(j\) appearing in its own total:

This is Pearson's \(r\), computed pairwise over non-missing observations. Items with \(r_{pb,j}<0.20\) are flagged.

Biserial and Polyserial Correlation

These indices estimate the latent linear association assuming the true variable is continuous and normal but observed only in coarsened categories. For an item with \(M\) categories, let \(\hat{F}_m\) be the cumulative proportion through category \(m\). The internal threshold is \(z_m=\Phi^{-1}(1-\hat{F}_m)\) and its standard normal ordinate is \(\varphi(z_m)\):

The last category is excluded because its threshold is \(+\infty\). For \(M=2\) this reduces to the classical biserial formula. \(\Phi^{-1}\) is computed via the Beasley-Springer-Moro rational approximation (accuracy \(\approx10^{-9}\)).

Reliability and Standard Error of Measurement

Cronbach's \(\alpha\)

where \(s_j^2\) is the sample variance of item \(j\) and \(s_T^2\) is the sample variance of total scores (Cronbach, 1951). Both use the \((n-1)\) denominator.

Alpha-if-Item-Deleted

For each item \(j\), \(\alpha\) is recomputed on the \((K-1)\)-item submatrix formed by removing column \(j\). An increase signals that item \(j\) is reducing internal consistency.

Standard Error of Measurement

Item Response Models

PACER supports binary and polytomous items in any mixture within a single test. The scaling constant \(D\) defaults to 1.7 (close to the normal-ogive metric; \(D=1\) gives the pure logistic metric). All probabilities are clamped to \([\varepsilon,\,1-\varepsilon]\) with \(\varepsilon\approx1.49\times10^{-8}\) to prevent \(\log 0\) in the likelihood.

Binary Item Response Models

Three-Parameter Logistic (3PL)

2PL and 1PL / Rasch

Setting \(c=0\) gives the 2PL. The 1PL additionally fixes \(a=1\) and estimates a free population standard deviation \(\sigma\) that scales the quadrature nodes: \(\theta_k^*=\sigma z_k\), where \(\{z_k\}\) are the standard GH nodes.

Graded Response Model (GRM)

The GRM (Samejima, 1969) models items with \(K+1\) ordered categories \((0,1,\ldots,K)\). The boundary response function gives the cumulative probability of responding in category \(k\) or higher:

with ordered thresholds \(d_1 < d_2 < \cdots < d_K\) and sentinels \(P^*(X\ge0)\equiv1\), \(P^*(X\ge K+1)\equiv0\). Category probabilities are differences of adjacent boundaries:

Boundary probabilities use the sign convention \(p_{kq}=\sigma(Da(d_k-\theta_q))\), so larger \(d_k\) means a harder threshold. Category probabilities are formed by differencing consecutive rows of the \(K\times Q\) boundary matrix.

Generalized Partial Credit Model (GPCM) and PCM

The GPCM (Muraki, 1992) uses adjacent-category transitions. For an item with \(K+1\) categories \((0,\ldots,K)\), define the log-numerator of category \(k\):

with \(\eta_0\equiv0\) and a leading \(d_0=0\) prepended to the step-parameter vector. Category probabilities follow the multinomial logit (softmax):

Evaluated via log-sum-exp: \(P(X=k)=\exp(\eta_k-M)/\sum_j\exp(\eta_j-M)\) where \(M=\max_j\eta_j\), preventing overflow at extreme \(\theta\). The PCM (Masters, 1982) is the GPCM with \(a\equiv1\).

Rating Scale Model (RSM)

The RSM (Andrich, 1978) constrains the PCM so that all items share identical step-difficulty offsets \(\tau_1,\ldots,\tau_K\) (\(\tau_0=0\) prepended by convention). Each item has a unique location \(\delta_i\), and its step parameters are

Response probabilities are evaluated with the PCM formula using these compound parameters. The free parameter vector is \([\delta_1,\ldots,\delta_I,\,\tau_1,\ldots,\tau_K]\).

Item Calibration: Marginal Maximum Likelihood via EM

PACER estimates item parameters by marginal maximum likelihood (MML), integrating the latent trait \(\theta\) out of the likelihood. Because \(\theta\) is unknown, direct maximization of the item parameters is not feasible; instead PACER forms the marginal likelihood and maximizes it using the Bock-Aitkin (1981) EM algorithm. Unique response patterns are identified before the E-step to eliminate redundant likelihood evaluations.

The Marginal Likelihood

The goal is to estimate the full vector of item parameters \(\bm{\gamma}\). Since the person parameters \(\theta_i\) are also unknown, they are integrated out under an assumed population distribution \(g(\theta\mid\mu,\sigma)\):

The integral has no closed form and is evaluated by Gauss-Hermite quadrature over \(Q\) nodes:

In practice, logs of individual pattern likelihoods are accumulated rather than computing products, ensuring numerical stability. When the data contain both binary and polytomous items (\(J\) binary, \(M\) polytomous), the pattern likelihood factors as a product over both item types:

where \(P_j\) is the binary ICC (1PL/2PL/3PL) and \(P_l\) is the polytomous category probability (GPCM or GRM). Both item types are estimated within the same EM loop.

Model Identification

For the 2PL and 3PL, the population distribution is normally fixed at \(N(0,1)\) to identify the latent scale. For the 1PL (Rasch model), only the mean need be fixed; PACER therefore estimates \(N(0,\sigma)\), with \(\sigma\) as a free parameter interpretable as the common discrimination across all items.

When anchor items are present in a single-group calibration, the scale is externally defined by the anchor constraints and the \(N(0,1)\) identification constraint is no longer required. In this case PACER automatically estimates \(\mu\) and \(\sigma\) for the single group using the same closed-form M-step update as the multigroup focal groups (see §Multigroup below). This behavior can be overridden by explicitly setting EstimatePopulation = false in the request.

Starting Values

Because the EM algorithm is iterative, starting values are required. Three strategies are available. Random draws each parameter from a uniform distribution over the admissible parameter space — useful when the likelihood may be multimodal; running the analysis multiple times with different random starts is a recommended diagnostic. Fixed sets \(a_j=1\), \(b_j=0\), \(c_j=0.25\) for all binary items and uses ordered random starting values for polytomous step parameters. Logistic uses logistic regression intercepts to initialize the difficulty parameters (and proportional-odds intercepts for polytomous step parameters), with \(a\) and \(c\) initialized as in the fixed option — this is the default and generally performs best.

E-Step: Posterior Weights and Sufficient Statistics

Given current item parameters, the likelihood of response pattern \(\bm{x}_i\) at node \(\theta_k\) is

The posterior weight is

The marginal log-likelihood used to monitor convergence is

Binary Sufficient Statistics

For each binary item \(j\) and node \(k\), summed over respondents \(\mathcal{C}_j\) who answered item \(j\):

\(r_{jk}\) is the expected number of correct responses at node \(k\); \(f_{jk}\) is the expected total count.

Polytomous Sufficient Statistics

For polytomous item \(j\), the sufficient statistic at node \(k\) for category \(s\) is

In the pattern-based implementation, the posterior weight for unique pattern \(p\) at node \(k\) is \(r_{p,k}=c_p\cdot w_k\cdot L_k/m_p\) where \(c_p\) is the pattern count and \(m_p=\sum_k L_k w_k\).

M-Step: Binary Items (2PL / 1PL / 3PL)

The M-step maximizes the expected complete-data log-likelihood for item \(j\):

PACER computes analytical gradients for the 2PL. The gradients with respect to \(a_j\) and \(b_j\) are:

In the multigroup path, the gradient accumulates across all groups \(g\) using group-scaled nodes \(\theta_{gk}=\mu_g+\sigma_g z_k\):

For the 3PL, \(c_j\) is parameterized on the logit scale internally; its gradient is computed by forward differences.

M-Step: GPCM and PCM

The complete-data log-likelihood for GPCM item \(j\) is the multinomial log-likelihood in expected category counts:

The analytical gradient with respect to \(a_j\) uses the conditional mean of the log-numerator:

For step parameter \(d_j\) (\(j=1,\ldots,K\)), let \(P_k^{\ge j}=\sum_{s\ge j}P(X=s\mid\theta_k)\) be the exceedance probability. Then:

This is a score-residual: model-expected exceedance minus the observed indicator. For the PCM, \(a_j\equiv1\) and the \(a\)-gradient is not computed.

M-Step: GRM

The complete-data log-likelihood for GRM item \(j\) is:

Let \(f_k^{(m)}=P^*(X\ge m\mid\theta_k)\cdot[1-P^*(X\ge m\mid\theta_k)]\) be the logistic slope of the \(m\)-th boundary. The gradient with respect to \(a_j\) depends on which category was observed:

The term \(\Delta_k(s)\) depends on the observed category \(s\):

The gradient with respect to threshold \(d_m\) is nonzero only for the two categories adjacent to it:

GRM item gradients are parallelized over items at the M-step.

M-Step: RSM

The RSM free parameter vector is \(\bm{\xi}=[\delta_1,\ldots,\delta_I,\,\tau_1,\ldots,\tau_K]\). Within the optimizer, item \(i\)'s step vector is \(\bm{d}_i=[0,\,\delta_i+\tau_1,\ldots,\delta_i+\tau_K]\) and probabilities are evaluated with the PCM formula (\(a=1\)). The gradients follow from the GPCM step gradient with \(a=1\):

The \(\tau_t\) gradient is the sum of per-item step gradients across all \(I\) items, since \(\partial d_{it}/\partial\tau_t=1\) for all \(i\). Item gradients are computed in parallel; \(\tau\) gradients are accumulated sequentially to avoid data races.

Parameter Priors

Priors are incorporated by augmenting the marginal log-likelihood with log-prior terms, one per item parameter:

Priors on \(a\) and \(b\) are optional. A prior on \(c\) is applied by default when estimating the 3PL and can be disabled for very large samples.

Beta Prior for the Guessing Parameter

The beta distribution is parameterized by shape parameters \(\alpha\) and \(\beta\), but it is more natural to specify a prior mean and standard deviation. PACER converts these to shape parameters via the method of moments. Starting from

solving for \(\alpha\) and \(\beta\) gives

The defaults \(\mu_c=0.20\) and \(\sigma_c=0.087\) correspond to \(\mathcal{B}(4,16)\), which places most mass below 0.35 and is appropriate for items with three to five response options. The logit-normal prior is described in the next section.

Logit-Normal Prior for the Guessing Parameter

The logit-normal prior places a normal distribution directly on the internal logit-scale parameter \(\lambda_j = \operatorname{logit}(c_j) = \ln(c_j/(1-c_j))\), which is the parameterization PACER uses internally during optimization. Because the optimizer works on \(\lambda_j\), no Jacobian correction is required in the objective function — the log-prior contribution is simply the log-density of a normal:

The gradient of the log-prior with respect to \(\lambda_j\) is therefore

and the corresponding contribution to the negative log-likelihood gradient used by the optimizer is \((\lambda_j - \mu_\lambda)/\sigma_\lambda^2\).

The default hyperparameters \(\mu_\lambda = -1.386\) and \(\sigma_\lambda = 0.544\) are chosen so that the prior mean on the probability scale is \(\operatorname{logit}^{-1}(-1.386) \approx 0.20\), matching the beta prior default mean. The prior is specified in the UI by entering the mean and standard deviation on the logit scale directly.

The logit-normal is preferable to the beta prior when software interoperability requires a prior specified on the logit scale (for example, mirt's parprior = list(c(i, "norm", mu, sd)) applies a normal prior to the internally stored logit of \(c\), which is exactly equivalent to PACER's logit-normal prior). It may also be preferred when the guessing parameter is expected to cluster near zero, since the logit transform spreads mass more evenly in the tails of the probability scale.

Per-Item Prior Overrides

In addition to the global priors described above, PACER supports per-item prior overrides that selectively replace or suppress the global prior for any individual item–parameter pair. The resolution rule is applied independently for each combination of item \(j\) and parameter \(\psi\in\{a,b,c,d\}\):

1. If an explicit per-item prior is specified for item \(j\), parameter \(\psi\): use that prior family and hyperparameters.
2. If an explicit None override is set for item \(j\), parameter \(\psi\): apply no prior for that parameter (pure MLE).
3. Otherwise: inherit the global prior for \(\psi\) (if the global prior is enabled).

The table below lists the parameters that accept per-item priors and their associated prior families:

Augmented Log-Posterior with Per-Item Priors

Replacing the global prior family with the resolved per-item family \(f_j^\psi\) for each item \(j\) and parameter \(\psi\), the augmented log-posterior becomes

where \(f_j^\psi\) evaluates to the per-item prior if one is set, the global prior if the item has no override, or 1 (no contribution) if an explicit None override is in effect for that pair. The gradient contributions to the optimizer's negative log-posterior objective are identical in form to the global-prior case — only the hyperparameters differ:

The hyperparameters \((\mu_a,\sigma_a)\), \((\mu_b,\sigma_b)\), \((\alpha,\beta)\), or \((\mu_\lambda,\sigma_\lambda)\) are taken from the per-item override when present, otherwise from the global prior setting.

Prior on Polytomous Threshold Parameters

For polytomous items (GRM, GPCM, PCM), a single normal prior \(\mathcal{N}(\mu_d,\sigma_d^2)\) is applied independently to each threshold or step parameter \(d_{jk}\). The joint log-prior contribution for item \(j\) with \(K_j\) categories is

The gradient contribution to the optimizer objective (negative log-posterior) for threshold \(d_{jk}\) is

This term is added to the EM sufficient-statistic gradient before passing to L-BFGS-B, analogously to the \(a\) and \(b\) prior gradient contributions for binary items.

Multigroup Calibration (Bock-Zimowski)

For data from multiple groups, PACER implements the Bock-Zimowski (1997) multigroup MML-EM. Each group \(g\) has distribution \(\theta\sim N(\mu_g,\sigma_g^2)\); the reference group (\(g=0\)) is fixed at \(N(0,1)\) for identification.

Group-Scaled Nodes

Weights \(\{w_k\}\) are invariant because the Jacobian and density ratio cancel exactly.

Multigroup E-Step

For respondent \(i\) in group \(g\), probabilities are evaluated at \(\theta_{gk}\):

Sufficient statistics \(r_g[j,k]\) and \(f_g[j,k]\) are accumulated per group on the same quadrature index \(k\).

Group Distribution Update (Closed Form)

Let \(W_{gk}=\sum_{i\in g}h_{ik}\). After the item M-step, the general update for a group with both \(\mu\) and \(\sigma\) free is:

For fixed-discrimination models (1PL, PCM, RSM) calibrated without anchors, the mean is identified by convention and held at \(\mu_g=0\); only the scale is free. In this case the second moment around zero is used directly:

Both updates are clamped to \(\hat{\sigma}_g\in[0.05,\,10]\) to prevent degenerate quadrature grids. The outer loop iterates until

where \(\tau_\ell=\max(\text{Tol},\,10^{-5})\). The convergence condition on the group parameters requires that both the mean and standard deviation shift by less than \(10^{-4}\) for every non-reference group.

Optimizer: L-BFGS-B

Item parameters are estimated with L-BFGS-B, a limited-memory quasi-Newton method with box constraints. The inverse Hessian approximation is built from the most recent \(m=10\) curvature pairs \(\{(\bm{s}_t,\bm{y}_t)\}\) where \(\bm{s}_t=\bm{x}_{t+1}-\bm{x}_t\) and \(\bm{y}_t=\nabla f_{t+1}-\nabla f_t\). The two-loop recursion computes the search direction \(\bm{p}_t=-H_t^{-1}\nabla f_t\) without explicitly forming \(H_t^{-1}\). Box constraints are enforced by projecting onto the active constraint set. Convergence: \(|\Delta f|/(0.1+|f|)<\tau\) (default \(\tau=10^{-8}\)). Convergence is evaluated relative to the current function value to ensure scale-invariant termination.

Ability Estimation

Given calibrated item parameters, PACER estimates each respondent's latent trait \(\hat\theta\) from their observed response vector. Four methods are available. All operate on sparse response vectors so that respondents with missing items are scored without imputation.

Maximum Likelihood Estimation (MLE)

A dense grid from \(-8\) to \(+8\) in steps of 0.05 locates the global basin; gradient descent with backtracking then refines from the grid minimum. Extreme scores: perfect and null scores produce unbounded estimates. PACER subtracts (or adds) 0.5 from the response of the item with the smallest discrimination \(a_j\). Standard error from observed information via numerical second differentiation:

Maximum A Posteriori (MAP)

The normal prior \(N(\mu,\sigma^2)\) (defaults: \(\mu=0\), \(\sigma=1\)) is added to the log-likelihood and the same minimizer is used. MAP is always finite, including for extreme scores.

Expected A Posteriori (EAP)

EAP is the posterior mean under a normal prior, approximated on a uniform quadrature grid spanning \(\mu\pm4\sigma\) in \(Q\) steps with weights proportional to \(\phi((\theta_k-\mu)/\sigma)\) (normalized to sum to one):

The posterior standard deviation serves as the standard error:

TCC Scoring

TCC scoring (Lord, 1980) inverts the Test Characteristic Curve \(\tau(\theta)=\sum_j E[X_j\mid\theta]\) to find

solved on \([-8,8]\) by Brent's root-finding method. The SE is \(1/\sqrt{I(\hat\theta_{\text{TCC}})}\) from the test information function.

Boundary Condition

A unique TCC estimate exists only when the raw score exceeds the sum of the guessing parameters across all binary items:

When this condition is not met — that is, when the observed score is no higher than what could be achieved by guessing alone — TCC assigns the same lower-bound \(\hat\theta\) to all such examinees. Similarly, a perfect raw score produces the upper-bound value defined by the search interval. The default interval \([-8,8]\) covers the range used in practice; if scores fall outside this range the interval endpoints can be widened.

Post-Calibration EAP Thetas

After calibration, PACER computes a \(\hat\theta\) for every respondent using 11-point GH quadrature with a diffuse prior \(N(0,10)\) (\(\sigma=\sqrt{10}\approx3.162\)). These thetas are used for DIF matching and item-fit analysis, not score reporting.

Test Characteristic Curve, Information, and SEM

Test Characteristic Curve (TCC)

When proportional TCC is enabled, output is divided by \(X_{\max}=\sum_j\max_k k\) to give a 0–1 scale.

Binary Item Information

For 1PL and 2PL:

For 3PL:

GRM Item Information

Samejima (1969):

where \(f_k(\theta)=P^*(X\ge k\mid\theta)[1-P^*(X\ge k\mid\theta)]\) is the logistic slope of the \(k\)-th boundary, with sentinels \(f_0\equiv0\) and \(f_{K+1}\equiv0\). Concretely, consecutive logistic-slope values are differenced and each squared difference is divided by the corresponding category probability, then summed and scaled by \(D^2a_j^2\).

GPCM / PCM / RSM Item Information

Information equals \(D^2a^2\) times the conditional variance of the item score (Muraki, 1992):

The test-level SEM at each \(\theta\) is

IRT Equating and Linking

IRT linking places parameter estimates from a new test form onto the scale of an existing item bank using common anchor items. PACER supports six methods in three families. Mixed binary and polytomous items are supported; GRM uses the cumulative-logistic expected score; GPCM/PCM/RSM use the adjacent-category expected score.

The Linear Scale Transformation

All methods estimate \((A,B)\) mapping the new-form theta scale to the bank scale: \(\theta_{\text{bank}}=A\theta_{\text{form}}+B\). Item parameters transform as

The reverse transformation is \(a_{\text{rev}}=Aa\), \(b_{\text{rev}}=(b-B)/A\), \(d_{\text{rev},k}=(d_k-B)/A\).

Stocking-Lord Methods

Stocking and Lord (1983) minimize the weighted integral of squared TCC differences over \(Q\) GH nodes. Let \(\tau_1(\theta)\) be the bank TCC and \(\tau_2^*(A,B;\theta)\) the linked new-form TCC.

Haebara Methods

Haebara (1980) replaces TCC comparison with item-by-item ICC comparison. For binary items \(P_{jk}\) is the 3PL ICC; for polytomous items it is \(E[X_j\mid\theta_k]\).

All four characteristic-curve objectives are minimized by two-dimensional Nelder-Mead simplex (starting simplex side 0.5, up to 2,000 iterations, tolerance \(10^{-10}\)).

Mean/Sigma and Mean/Mean

Mean/Sigma (MS)

The location pool: \(b\)-parameters for binary items; all \(d\)-parameters pooled as a flat vector for GPCM/PCM (or averaged per item when the single-location option is selected); GRM boundary parameters; item location \(\delta_i\) only for RSM items.

Mean/Mean (MM)

dropD² Iterative Anchor Purification

After the initial equating, a per-item discrepancy statistic is computed:

For polytomous items \(P_j\) is the expected score \(E[X_j\mid\theta_k]\). Values are standardized to

The item with the largest \(|\tilde{D}^2_j|\) exceeding the cutoff (default 2.5) is removed; equating is re-run on the reduced anchor set. The process repeats until no item exceeds the cutoff or fewer than two anchors remain.

Cognitive Diagnostic Models

CDMs classify examinees into mastery profiles over \(K\) binary latent attributes \(\bm{\alpha}_i=(\alpha_{i1},\ldots,\alpha_{iK})\in\{0,1\}^K\). The Q-matrix \(\bm{Q}\in\{0,1\}^{J\times K}\) specifies which attributes are required by each item.

DINA Model

The Deterministic Inputs, Noisy And-gate (DINA) model sets the ideal response indicator as \(\xi_{ij}=\prod_{k=1}^{K}\alpha_{ik}^{q_{jk}}\). A respondent achieves \(\xi_{ij}=1\) only if they have mastered all required attributes. The response probability accommodates slipping (\(s_j\)) and guessing (\(g_j\)):

DINO Model

The Deterministic Inputs, Noisy Or-gate (DINO) model requires mastery of at least one required attribute:

EM Estimation of Item Parameters

Item parameters (\(s_j\), \(g_j\)) and class probabilities (\(\pi_c\)) are estimated jointly by marginal maximum likelihood over the \(L=2^K\) latent classes \(\bm{\alpha}^{(c)}\) with prior class probabilities \(\pi_c\ge0\), \(\sum_c\pi_c=1\). The algorithm is initialized with equal class probabilities (\(\pi_c=1/L\)) and uniform starting values for \(g_j\) and \(s_j\). The population distribution is re-estimated at each step, so the posteriors are empirical Bayes estimates.

E-Step

Let \(\Pr(\bm{\alpha}_l\mid\bm{x}_i)\) be the posterior probability that examinee \(i\) has skill profile \(\bm{\alpha}_l\). The expected number of examinees with profile \(\bm{\alpha}_l\) is

and the expected number of those examinees who answered item \(j\) correctly is

Let \(I^0_{jl}\) be the expected count of examinees who lack at least one attribute required for item \(j\) (i.e., \(\xi_{ijl}=0\)) and \(R^0_{jl}\) be the subset of those who answered correctly. Analogously, \(I^1_{jl}\) and \(R^1_{jl}\) correspond to examinees who possess all required attributes (\(\xi_{ijl}=1\)).

M-Step

With these expected counts, the guessing and slip parameters are updated in closed form:

Class probabilities are updated as \(\hat\pi_c=N^{-1}\sum_i\tau_{ic}\). The EM loop alternates between E-step and M-step until the change in parameters across iterations falls below the convergence tolerance.

Skill Profile Estimation

Once item parameters are calibrated, each examinee's skill profile is estimated from the likelihood of their response pattern \(\bm{x}_i\) given each of the \(L=2^K\) possible profiles:

Three estimators are available.

Maximum Likelihood (MLE)

The MLE selects the profile that maximizes the likelihood across all \(L\) possible classifications:

Because skill profiles are discrete, this is computed by exhaustive evaluation over all \(L\) profiles. The MLE may be multimodal — two or more profiles can share the same maximum likelihood. PACER flags such cases and reports all profiles achieving the maximum.

Maximum A Posteriori (MAP)

The MAP incorporates a prior distribution \(\Pr(\bm{\alpha}_l)\) over skill profiles (supplied by the user), maximizing the posterior:

where \(\Pr(\bm{\alpha}_l)\) is a discrete prior such that \(\sum_l\Pr(\bm{\alpha}_l)=1\). Like the MLE, the MAP can be multimodal.

Expected A Posteriori (EAP)

The EAP uses the full posterior distribution over profiles. The posterior for profile \(\bm{\alpha}_l\) given \(\bm{x}_i\) is

Let \(\bm{\tau}_i = [\Pr(\bm{\alpha}_1\mid\bm{x}_i),\ldots,\Pr(\bm{\alpha}_L\mid\bm{x}_i)]\) be the \(1\times L\) posterior vector and \(\bm{P}\) the \(L\times K\) binary matrix of all possible skill profiles. The EAP estimate is

This is a \(1\times K\) vector where each element is the marginal posterior probability of mastery for the corresponding attribute — equivalent to \(\sum_{l:\,\alpha_k^{(l)}=1}\Pr(\bm{\alpha}_l\mid\bm{x}_i)\). Unlike the MLE and MAP, the EAP yields continuous mastery probabilities between 0 and 1 and is always unique.

Examinee Skill Profile Estimation

Once item parameters are calibrated, each examinee's mastery profile is characterized by their full posterior distribution over the \(2^K\) latent classes. Three quantities are reported.

Modal Profile (MAP Classification)

The modal-class estimate assigns the examinee to the latent class with the highest posterior probability:

The estimated attribute profile is \(\hat{\bm{\alpha}}_i = \bm{\alpha}^{(\hat{c}_i)}\), a binary vector of length \(K\) indicating mastery (1) or non-mastery (0) of each attribute. This is the primary classification output.

Marginal Attribute Mastery Probabilities

The marginal probability of mastering attribute \(k\) is obtained by summing the posterior over all latent classes in which attribute \(k\) is mastered:

These \(K\) probabilities provide a continuous measure of mastery confidence. A value near 1 indicates strong evidence of mastery; near 0 indicates strong evidence of non-mastery; intermediate values indicate classification uncertainty for that attribute.

Classification Certainty

The posterior entropy of the class distribution for examinee \(i\) quantifies overall classification certainty:

Low entropy indicates that most posterior weight is concentrated on a single class (high certainty). High entropy indicates spread across multiple classes, flagging examinees whose profiles should be interpreted cautiously. The minimum possible entropy is 0 (perfect certainty); the maximum is \(\ln(2^K)=K\ln 2\) (uniform uncertainty).

Differential Item Functioning

DIF occurs when respondents from different groups who have the same ability have systematically different probabilities of answering an item correctly. PACER implements the Mantel-Haenszel procedure (Holland & Thayer, 1988) for binary items and the Liu-Agresti (1996) generalization for polytomous items. All procedures first stratify respondents into \(K\) equal-frequency bins (default \(K=10\)) using a matching variable — either total score or IRT theta. Binning uses linear interpolation on the empirical quantile function so that equal proportions of respondents fall in each stratum; the lowest observed score is included in the first bin.

Mantel-Haenszel Procedure (Binary Items)

Within each stratum \(s\), a \(2\times2\) contingency table is formed:

The MH log-odds ratio pools across strata:

The ETS D-DIF statistic:

Negative D-DIF means the reference group has higher odds of answering correctly at the same ability — the item favors the reference group. Note: some published implementations invert this sign; PACER applies the convention stated above.

Chi-square with Yates continuity correction (1 df):

The p-value is computed from the regularized incomplete gamma function via series and continued-fraction expansions (accurate to \(\approx10^{-10}\)).

Generalized MH Procedure (Polytomous Items)

For polytomous items with score values \(v_1\lt\cdots\lt v_S\), let \(N_{g,k}(v)\) be the count of respondents in group \(g\) with score \(\le v\) in stratum \(k\). The Liu-Agresti (1996) one-df score test:

The test statistic \(T^2/\mathrm{Var}(T)\) is chi-squared with one degree of freedom.

Standardized Mean Difference (SMD)

Positive SMD means the focal group scores higher (item favors focal); negative means the item favors the reference group. Some published implementations use \(|\text{SMD}|\) for the direction determination; PACER uses the raw signed value as defined above. The standardized SMD is

ETS DIF Classification

Item Fit: Yen's Q1

After calibration, PACER assesses item fit using Yen's (1981) Q1 statistic. Q1 partitions respondents into score groups ordered by estimated ability, then compares observed item responses within each group to the values expected under the fitted IRT model. Systematic discrepancies between observed and expected values indicate that the model does not adequately account for response behavior across the ability continuum.

Person ability estimates are EAP scores computed under a diffuse \(N(0,10)\) prior, which approximates maximum likelihood estimates and avoids boundary problems for extreme scores. Respondents with missing responses on a given item are excluded from that item's computation.

Binary Items

Score Group Formation

Let \(G\) denote the number of score groups (default \(G=10\), following Yen, 1981). Respondents with valid responses to item \(j\) are sorted by their EAP estimate \(\hat\theta_i\) and divided into \(G\) equal-frequency groups. Let \(n_g\) be the number of respondents in group \(g\) and \(\bar\theta_g\) be their mean EAP estimate.

Observed and Expected Proportions

Within group \(g\), the observed proportion correct is

The expected proportion correct at the group mean ability is

where \(P(\cdot)\) is the appropriate binary IRF — 3PL for items with a nonzero guessing parameter, 2PL or 1PL otherwise.

The Q1 Statistic

The statistic follows a weighted chi-square form, with each bin's contribution scaled by the variance of a Bernoulli variable with mean \(E_g\):

where \(G^*\) is the number of non-empty bins actually formed (which may be less than \(G\) when many respondents share the same EAP estimate). Under the null hypothesis that the item fits the model, \(Q_1\) is approximately chi-square distributed with degrees of freedom

The reduction by one additional degree of freedom for the 3PL reflects the extra guessing parameter estimated from the data.

Polytomous Items

Score Group Formation

Score group formation is identical to the binary case: respondents are sorted by EAP and divided into \(G\) equal-frequency bins. Let \(K\) denote the number of response categories (0 through \(K-1\)).

Observed and Expected Category Proportions

Within group \(g\), the observed proportion in category \(k\) is

The expected proportion in category \(k\) at the group mean ability is

where \(\bm{\psi}_j\) denotes the calibrated parameters for item \(j\) and \(P(\cdot)\) is evaluated using the model-appropriate category response function. For the GRM:

with sentinels \(P^*(X \ge 0)=1\) and \(P^*(X \ge K)=0\). For the GPCM and PCM:

with the convention that the empty sum (\(k=0\)) equals zero.

The Q1 Statistic for Polytomous Items

The statistic generalizes the binary form to a sum over both bins and categories:

Under the null hypothesis of model fit, \(Q_1\) is approximately chi-square distributed with

The \(-K\) correction accounts for the \(K\) item parameters estimated from the data.

Standardization and Flagging

Standardized Q1z

To place Q1 on a common scale across items with different numbers of categories and bins, the statistic is standardized using the mean and variance of the chi-square distribution under the null:

Under correct model specification, \(Q_{1z}\) is approximately standard normal. Large positive values indicate misfit.

Practical Flagging Criterion

PACER uses Yen's (1981) sample-size-adjusted practical cutoff rather than a fixed critical value. An item is flagged for review when

where \(N\) is the number of respondents with valid responses on that item. This heuristic recognizes that with large samples even trivially small model-data discrepancies become statistically detectable. At \(N=1500\) the cutoff equals 4.0; at \(N=500\) it equals approximately 1.33.

Gauss-Hermite Quadrature

All integrals over the latent-trait distribution are approximated by Gauss-Hermite quadrature with \(Q\) nodes (default \(Q=21\)):

where \(\phi\) is the standard normal density. Nodes \(\{z_k\}\) are roots of the probabilists' Hermite polynomial \(\text{He}_Q(x)\); weights satisfy \(w_k=v_{k,1}^2\sqrt{2\pi}\) where \(v_{k,1}\) is the first component of the \(k\)-th eigenvector. For all values of \(Q\), nodes and weights are computed on demand via the symmetric QL algorithm applied to the tridiagonal Jacobi matrix whose off-diagonal elements are \(\beta_i=\sqrt{i+1}\) for \(i=0,1,\ldots,Q-2\) (the probabilists' He-polynomial recurrence). Eigenvalues of this matrix are the GH nodes; weights follow from the squared first eigenvector components scaled by \(\sqrt{2\pi}\). For multigroup calibration, nodes are shifted and scaled: \(\theta_{gk}=\mu_g+\sigma_g z_k\). Weights are invariant because the Jacobian and density ratio cancel exactly.

Optimization Routines

L-BFGS-B (Item Calibration)

Item parameters are estimated by L-BFGS-B. The two-loop recursion computes the quasi-Newton direction \(\bm{p}_t=-H_t^{-1}\nabla f_t\) from \(m=10\) curvature pairs without forming \(H_t^{-1}\) explicitly. Box constraints are enforced by projecting onto the active constraint set. Step length satisfies the strong Wolfe conditions via backtracking. Convergence: \(|\Delta f|/(0.1+|f|)<\tau\) (default \(\tau=10^{-8}\)).

Gradient Descent with Backtracking (Theta Estimation)

MLE and MAP theta estimation use gradient descent with backtracking. At each step: \(\theta_{\text{new}}=\theta-\alpha[-\ell'(\theta)]\) where \(\ell'(\theta)\approx[\ell(\theta+h)-\ell(\theta-h)]/(2h)\), \(h=10^{-6}\). Step size \(\alpha\) is increased by 1.2 on progress and halved on rejection; convergence when \(\alpha<10^{-10}\).

Brent's Method (TCC Inversion)

TCC scoring inverts \(\tau(\theta)=X_{\text{obs}}\) on \([-8,8]\). Brent's method combines bisection (guaranteed progress), secant (linear interpolation), and inverse quadratic interpolation (parabolic fit through three points). The method falls back to bisection whenever the other two strategies would not make sufficient progress, giving guaranteed \(O(\log\varepsilon^{-1})\) worst-case convergence. Tolerance: \(10^{-10}\).

Nelder-Mead Simplex (Equating Objectives)

The two-parameter SL and HB equating objectives are minimized over \((A,B)\) by Nelder-Mead simplex. The simplex of three vertices is updated by reflection, expansion, contraction, or shrinkage at each iteration. Initial simplex side: 0.5; maximum iterations: 2,000; convergence tolerance on objective range: \(10^{-10}\).

References

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

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

Bock, R. D., & Zimowski, M. F. (1997). Multiple group IRT. In W. J. van der Linden & R. K. Hambleton (Eds.), Handbook of modern item response theory. Springer.

Cronbach, L. J. (1951). Coefficient alpha and the internal structure of tests. Psychometrika, 16, 297–334.

Doran, H. (2023). A collection of numerical recipes useful for building scalable psychometric applications. Journal of Educational and Behavioral Statistics. https://doi.org/10.3102/10769986221116905

Haebara, T. (1980). Equating logistic ability scales by a weighted least squares method. Japanese Psychological Research, 22, 144–149.

Holland, P. W., & Thayer, D. T. (1988). Differential item performance and the Mantel-Haenszel procedure. In H. Wainer & H. I. Braun (Eds.), Test validity. Lawrence Erlbaum.

Liu, I. M., & Agresti, A. (1996). Mantel-Haenszel-type inference for cumulative odds ratios with a stratified ordinal response. Biometrics, 52, 1223–1234.

Lord, F. M. (1980). Applications of item response theory to practical testing problems. Erlbaum.

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

Muraki, E. (1992). A generalized partial credit model: Application of an EM algorithm. Applied Psychological Measurement, 16, 159–176.

Samejima, F. (1969). Estimation of latent ability using a response pattern of graded scores. Psychometrika Monograph Supplement, No. 17.

Yen, W. M. (1981). Using simulation results to choose a latent trait model. Applied Psychological Measurement, 5, 245–262.

Stocking, M. L., & Lord, F. M. (1983). Developing a common metric in item response theory. Applied Psychological Measurement, 7, 201–210.