The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.

Psychometric Tools in tirt: Information, Scoring, Fit, DIF, and Mixtures

1. Overview

Once item and person parameters have been estimated, a psychometric study usually continues with a battery of follow-up analyses: evaluating measurement precision, building score-conversion tables, checking person and item fit, screening for local dependence and differential item functioning, summarizing reliability, and sometimes probing for latent subpopulations. This vignette walks through the post-estimation tools in tirt. Every function accepts the item-parameter data frame produced by binary_irt(), polytomous_irt(), or mixed_irt(), so the workflow is uniform across models.

library(tirt)

We start from a small simulated 2PL data set and a quick calibration.

set.seed(2025)
sim <- sim_irt(
  n_people = 600,
  item_structure = list(list(model = "2PL", n_items = 12))
)
#> ----------------------------------------------------------------
#> Starting Simulation for N = 600 examinees...
#> >> Ability (Theta): Generated from N(mean=0.00, sd=1.00).
#> >> Block 1: 12 items using 2PL
#>    - Discrimination (a): Default values used (Fixed at 1).
#>    - Difficulty (b): Default values used (Fixed at 0).
#> ----------------------------------------------------------------
#> Constructing final data frames...
#> Simulation Complete.
#> Summary: 12 items, 600 examinees.
#> ----------------------------------------------------------------

fit <- binary_irt(sim$resp, model = "2PL", method = "EM",
                  control = list(max_iter = 20, verbose = FALSE))
head(fit$item_params)
#>     item discrimination discrimination_se difficulty difficulty_se number
#> 1 item_1          1.010             0.108     -0.057         0.089    600
#> 2 item_2          0.946             0.105     -0.001         0.094    600
#> 3 item_3          0.974             0.106      0.073         0.092    600
#> 4 item_4          0.866             0.102     -0.099         0.102    600
#> 5 item_5          0.906             0.103     -0.070         0.098    600
#> 6 item_6          0.970             0.106      0.016         0.092    600
#>   pvalue
#> 1  0.512
#> 2  0.500
#> 3  0.485
#> 4  0.518
#> 5  0.513
#> 6  0.497

2. Item and Test Information

item_info() returns the Fisher information of each item across a grid of ability values, and test_info() sums it into the test information function together with the conditional standard error of measurement (SEM).

theta_grid <- seq(-3, 3, by = 0.5)

# Item information (items in rows, theta in columns)
info <- item_info(fit$item_params, theta = theta_grid)
round(info[1:3, ], 3)
#>           -3  -2.5    -2  -1.5    -1  -0.5     0   0.5     1   1.5     2   2.5
#> item_1 0.047 0.074 0.110 0.156 0.205 0.243 0.255 0.236 0.194 0.145 0.101 0.067
#> item_2 0.047 0.070 0.102 0.140 0.180 0.212 0.224 0.212 0.180 0.140 0.102 0.070
#> item_3 0.043 0.066 0.098 0.139 0.183 0.220 0.237 0.227 0.195 0.151 0.109 0.075
#>            3
#> item_1 0.043
#> item_2 0.047
#> item_3 0.049

# Test information function and conditional SEM
tif <- test_info(fit$item_params, theta = theta_grid)
tif
#>    theta test_info       sem reliability
#> 1   -3.0 0.5479329 1.3509408   0.3539772
#> 2   -2.5 0.8379575 1.0924184   0.4559178
#> 3   -2.0 1.2398965 0.8980640   0.5535508
#> 4   -1.5 1.7458850 0.7568193   0.6358187
#> 5   -1.0 2.2902138 0.6607878   0.6960684
#> 6   -0.5 2.7337536 0.6048119   0.7321730
#> 7    0.0 2.9109843 0.5861113   0.7443099
#> 8    0.5 2.7414786 0.6039592   0.7327260
#> 9    1.0 2.2991850 0.6594973   0.6968948
#> 10   1.5 1.7505310 0.7558143   0.6364338
#> 11   2.0 1.2393432 0.8982645   0.5534405
#> 12   2.5 0.8341527 1.0949069   0.4547891
#> 13   3.0 0.5430402 1.3570131   0.3519287

# Where does the test measure most precisely?
tif$theta[which.max(tif$test_info)]
#> [1] 0

3. Summed-Score to Theta Conversion Table

Operational programs report scale scores from raw (summed) scores. score_table() produces the conversion using the Lord-Wingersky recursion, supporting expected a posteriori ("EAP"), weighted likelihood ("WLE"), and maximum likelihood ("MLE") scoring.

# EAP conversion table (0 to 12 correct)
score_table(fit$item_params, method = "EAP")
#>    summed_score  theta    se
#> 1             0 -1.845 0.643
#> 2             1 -1.468 0.605
#> 3             2 -1.131 0.572
#> 4             3 -0.826 0.548
#> 5             4 -0.541 0.532
#> 6             5 -0.269 0.522
#> 7             6 -0.004 0.519
#> 8             7  0.262 0.522
#> 9             8  0.534 0.531
#> 10            9  0.819 0.547
#> 11           10  1.124 0.571
#> 12           11  1.461 0.604
#> 13           12  1.840 0.643

# Maximum-likelihood conversion
score_table(fit$item_params, method = "MLE")
#>    summed_score  theta    se
#> 1             0 -4.000    NA
#> 2             1 -2.471 1.096
#> 3             2 -1.654 0.810
#> 4             3 -1.129 0.679
#> 5             4 -0.713 0.625
#> 6             5 -0.348 0.595
#> 7             6 -0.004 0.588
#> 8             7  0.339 0.595
#> 9             8  0.703 0.624
#> 10            9  1.117 0.678
#> 11           10  1.639 0.783
#> 12           11  2.449 1.054
#> 13           12  4.000    NA

4. Person Fit

person_fit() computes the lz standardized log-likelihood index to identify examinees whose response patterns are unlikely under the model.

pf <- person_fit(sim$resp, fit$item_params, fit$person_params)
head(pf)
#>   person n_items  theta loglik     lz  flag
#> 1      1      12  0.668 -7.585  0.102 FALSE
#> 2      2      12  0.325 -8.121  0.069 FALSE
#> 3      3      12  0.044 -8.410 -0.763 FALSE
#> 4      4      12  0.399 -8.259 -0.259 FALSE
#> 5      5      12  0.406 -8.186 -0.159 FALSE
#> 6      6      12 -1.591 -5.518  0.012 FALSE

# Number of examinees flagged as potentially misfitting
sum(pf$flag, na.rm = TRUE)
#> [1] 0

5. Item Fit

item_fit() returns the infit and outfit mean-square statistics and their standardized versions, for dichotomous and polytomous items alike.

item_fit(sim$resp, fit$item_params)
#>       item   n outfit outfit_z infit infit_z
#> 1   item_1 600  0.861   -3.893 0.892  -3.746
#> 2   item_2 600  0.875   -3.783 0.903  -3.518
#> 3   item_3 600  0.863   -3.989 0.901  -3.517
#> 4   item_4 600  0.894   -3.476 0.916  -3.257
#> 5   item_5 600  0.887   -3.550 0.909  -3.431
#> 6   item_6 600  0.871   -3.779 0.899  -3.621
#> 7   item_7 600  0.871   -3.849 0.902  -3.530
#> 8   item_8 600  0.848   -3.991 0.880  -3.998
#> 9   item_9 600  0.891   -3.482 0.911  -3.392
#> 10 item_10 600  0.834   -4.241 0.876  -4.036
#> 11 item_11 600  0.829   -4.278 0.875  -4.012
#> 12 item_12 600  0.869   -3.843 0.900  -3.561

6. Local Dependence (Yen’s Q3)

ld_stats() computes Yen’s Q3 residual correlations for every item pair. Large positive values signal that a pair depends on something beyond the common trait, which is the situation a testlet model is designed to absorb.

q3 <- ld_stats(sim$resp, fit$item_params)
round(q3[1:5, 1:5], 3)
#>        item_1 item_2 item_3 item_4 item_5
#> item_1  1.000 -0.039 -0.004 -0.067 -0.040
#> item_2 -0.039  1.000 -0.098  0.001 -0.155
#> item_3 -0.004 -0.098  1.000 -0.056 -0.105
#> item_4 -0.067  0.001 -0.056  1.000 -0.018
#> item_5 -0.040 -0.155 -0.105 -0.018  1.000

# Largest absolute residual correlation
attr(q3, "max_abs_q3")
#> [1] 0.1841257

7. Differential Item Functioning

dif() screens dichotomous items for DIF using the Mantel-Haenszel procedure (with the ETS delta effect size and A/B/C flags) and logistic regression (which separates uniform and non-uniform DIF). Here we plant DIF in item 3.

resp_dif <- sim$resp
grp <- rep(c("Reference", "Focal"), each = 300)
flip <- grp == "Focal" & resp_dif[[3]] == 1
resp_dif[[3]][flip] <- rbinom(sum(flip), 1, 0.55)

dif(resp_dif, group = grp)[, c("item", "MH_delta", "ETS_class", "LR_p")]
#>       item MH_delta ETS_class   LR_p
#> 1   item_1    0.175         A 0.6064
#> 2   item_2   -0.235         A 0.1189
#> 3   item_3    1.910         C 0.0000
#> 4   item_4   -0.587         A 0.2569
#> 5   item_5    0.599         A 0.3319
#> 6   item_6    0.155         A 0.3628
#> 7   item_7    0.434         A 0.4699
#> 8   item_8   -0.781         A 0.1236
#> 9   item_9   -0.820         A 0.2681
#> 10 item_10   -0.144         A 0.9380
#> 11 item_11   -0.262         A 0.0020
#> 12 item_12   -0.475         A 0.7313

8. Reliability

reliability() reports the empirical (marginal) reliability from person estimates and their standard errors, the model-based marginal reliability from the test information function, and classical Cronbach’s alpha from the raw responses.

reliability(person_params = fit$person_params,
            data = sim$resp,
            item_params = fit$item_params)
#>                   Index    Value
#> 1 empirical_reliability   0.7256
#> 2  marginal_reliability   0.6886
#> 3        cronbach_alpha   0.7070
#> 4               n_items  12.0000
#> 5             n_persons 600.0000

9. Test Characteristic Curve

tcc() returns the expected score of each item and the test characteristic curve, which maps ability onto the number-correct metric.

curves <- tcc(fit$item_params, theta = seq(-3, 3, by = 1))
curves$test_curve
#>   theta expected_score
#> 1    -3      0.6204589
#> 2    -2      1.5004155
#> 3    -1      3.2870410
#> 4     0      6.0113887
#> 5     1      8.7381246
#> 6     2     10.5220784
#> 7     3     11.3938127

10. Multidimensional and Mixed Simulation

Two new simulators extend the data-generation toolkit. sim_mirt() produces compensatory multidimensional data (aligned with mirt_binary()), and sim_tirt() produces forms that mix independent items with testlets (aligned with irt_trt()).

# Two correlated dimensions, simple structure
mdat <- sim_mirt(
  n_people = 400,
  dimension = 2,
  Sigma = matrix(c(1, 0.4, 0.4, 1), 2, 2),
  item_structure = list(
    list(model = "M2PL", n_items = 6, dims = 1),
    list(model = "M2PL", n_items = 6, dims = 2)
  )
)
#> ----------------------------------------------------------------
#> Starting Multidimensional Simulation (N = 400, D = 2)...
#> >> Ability (Theta): Generated from MVN with mean 0.00 and supplied Sigma.
#> >> Block 1: 6 items (M2PL), loading on dim(s) 1
#> >> Block 2: 6 items (M2PL), loading on dim(s) 2
#> ----------------------------------------------------------------
#> Constructing final data frames...
#> Simulation Complete.
#> Summary: 12 items, 400 examinees, 2 dimensions.
#> ----------------------------------------------------------------
head(mdat$true_params)
#>   item_id block model categories    a_Dim1 a_Dim2          d guessing
#> 1  item_1     1  M2PL          2 0.8192317      0  0.8491260        0
#> 2  item_2     1  M2PL          2 1.7109341      0  0.5145288        0
#> 3  item_3     1  M2PL          2 1.7031725      0 -1.0493497        0
#> 4  item_4     1  M2PL          2 0.9961762      0  1.4578654        0
#> 5  item_5     1  M2PL          2 1.5794204      0  0.5032779        0
#> 6  item_6     1  M2PL          2 1.2982607      0  1.3152575        0
# Independent items plus two testlets
tdat <- sim_tirt(
  n_people = 400,
  item_structure = list(
    list(model = "2PL",  n_items = 6),
    list(model = "2PLT", n_items = 4, testlet_id = "P1", testlet_var = 0.6),
    list(model = "GPCT", n_items = 3, categories = 3, testlet_id = "P2")
  )
)
#> ================================================================
#>    STARTING MIXED IRT/TRT SIMULATION (N = 400)
#> ================================================================
#> >> Ability (Theta): Generated from N(mean=0.00, sd=1.00).
#> >> Testlet Effects (Gamma):
#>    - Testlet 'P1': Generated Gamma ~ N(0, 0.60) [User Var].
#>    - Testlet 'P2': Generated Gamma ~ N(0, 0.50) [Default].
#> >> Block 1: 6 items (Model: 2PL, Independent)
#> >> Block 2: 4 items (Model: 2PLT, Testlet: P1)
#> >> Block 3: 3 items (Model: GPCT, Testlet: P2)
#> ================================================================
#> Constructing final data frames...
#> Simulation Complete.
#> Summary: 13 items (2 testlets), 400 examinees.
#> ================================================================
tdat$true_item_params[, c("item_id", "model", "testlet")]
#>    item_id model testlet
#> 1   item_1   2PL    <NA>
#> 2   item_2   2PL    <NA>
#> 3   item_3   2PL    <NA>
#> 4   item_4   2PL    <NA>
#> 5   item_5   2PL    <NA>
#> 6   item_6   2PL    <NA>
#> 7   item_7  2PLT      P1
#> 8   item_8  2PLT      P1
#> 9   item_9  2PLT      P1
#> 10 item_10  2PLT      P1
#> 11 item_11  GPCT      P2
#> 12 item_12  GPCT      P2
#> 13 item_13  GPCT      P2

11. Mixture (Latent-Class) IRT

mixture_irt() fits a mixture Rasch or 2PL model, in which the population is a blend of latent classes with class-specific item parameters.

set.seed(11)
N <- 300; J <- 8
b1 <- seq(-1.5, 1.5, length.out = J); b2 <- rev(b1)
theta <- rnorm(N); cls <- rep(1:2, each = N / 2)
rmat <- matrix(0, N, J)
for (i in 1:N) {
  b <- if (cls[i] == 1) b1 else b2
  rmat[i, ] <- rbinom(J, 1, 1 / (1 + exp(-(theta[i] - b))))
}
mdf <- as.data.frame(rmat); names(mdf) <- paste0("I", 1:J)

mix <- mixture_irt(mdf, n_class = 2, model = "Rasch",
                   control = list(max_iter = 40, verbose = FALSE))
mix$class_params
#>   class proportion
#> 1     1     0.4873
#> 2     2     0.5127
mix$model_fit
#>           Index     Value
#> 1 LogLikelihood -1533.969
#> 2           AIC  3101.938
#> 3           BIC  3164.902
#> 4       n_class     2.000
#> 5       entropy     0.613

12. Final Comment

Together with the estimation, calibration, and equating functions described in the companion vignette, these tools cover the routine post-estimation workflow of an operational testing program. For details on any function, use the help system, for example ?score_table.

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.