1. Overview
The tirt package provides a light-version user-friendly comprehensive toolkit for Item Response Theory (IRT) and Testlet Response Theory (TRT) analysis. It supports:
Data Simulation: Generating data for simple IRT and Testlet models.
Estimation: Marginal Maximum Likelihood (MML) and Joint Maximum Likelihood (JML) estimation.
Fixed Calibration: Estimating item parameters with fixed persons, or person parameters with fixed items.
Equating: Linking scales using Mean-Mean, Mean-Sigma, and Stocking-Lord methods.
This vignette walks through the primary workflows.
2. Data Simulation
Before estimating models, let’s generate some synthetic data using sim_irt (standard IRT) and sim_trt (testlet structure).
2.1 Standard IRT Simulation
We generate responses for 500 examinees on 15 items using the 2-Parameter Logistic (2PL) model.
set.seed(123)
# Define item structure: 1 block of 15 items using 2PL
sim_data <- sim_irt(
n_people = 500,
item_structure = list(
list(model = "2PL", n_items = 15)
),
theta_mean = 0,
theta_sd = 1
)
#> ----------------------------------------------------------------
#> Starting Simulation for N = 500 examinees...
#> >> Ability (Theta): Generated from N(mean=0.00, sd=1.00).
#> >> Block 1: 15 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: 15 items, 500 examinees.
#> ----------------------------------------------------------------
# Access responses and true parameters
head(sim_data$resp)
#> item_1 item_2 item_3 item_4 item_5 item_6 item_7 item_8 item_9 item_10
#> 1 1 0 1 0 1 1 1 1 1 1
#> 2 0 0 1 0 0 1 0 1 0 0
#> 3 1 1 1 1 1 1 1 1 1 1
#> 4 0 1 1 1 0 1 0 0 0 1
#> 5 0 0 1 0 1 1 1 1 1 0
#> 6 1 1 1 1 1 1 1 1 1 1
#> item_11 item_12 item_13 item_14 item_15
#> 1 0 0 0 1 0
#> 2 1 0 0 0 0
#> 3 1 1 1 0 1
#> 4 0 1 0 0 1
#> 5 1 0 1 0 0
#> 6 1 1 1 1 1
head(sim_data$true_params)
#> item_id block model categories discrimination guessing difficulty
#> 1 item_1 1 2PL 2 1 0 0
#> 2 item_2 1 2PL 2 1 0 0
#> 3 item_3 1 2PL 2 1 0 0
#> 4 item_4 1 2PL 2 1 0 0
#> 5 item_5 1 2PL 2 1 0 0
#> 6 item_6 1 2PL 2 1 0 02.2 Testlet Simulation
For Testlet Response Theory, we define blocks that share a common local dependency (Gamma).
# Define 2 testlets with 5 items each
trt_data <- sim_trt(
n_people = 500,
item_structure = list(
list(model = "2PLT", n_items = 5, testlet_id = "T1"),
list(model = "2PLT", n_items = 5, testlet_id = "T2")
)
)
#> ================================================================
#> STARTING TRT SIMULATION (N = 500)
#> ================================================================
#> >> Ability (Theta): Generated from N(mean=0.00, sd=1.00).
#> >> Testlet Effects (Gamma):
#> - Testlet 'T1': Generated Gamma ~ N(0, 0.50) [Default].
#> - Testlet 'T2': Generated Gamma ~ N(0, 0.50) [Default].
#> >> Block 1: 5 items (Model: 2PLT, Testlet: T1)
#> - Discrimination a: Default used (1.00).
#> - Difficulty/Loc b: Default used (0.00).
#> >> Block 2: 5 items (Model: 2PLT, Testlet: T2)
#> - Discrimination a: Default used (1.00).
#> - Difficulty/Loc b: Default used (0.00).
#> ================================================================
#> Constructing final dataframes...
#> Simulation Complete.
#> Summary: 10 items, 500 examinees.
#> ================================================================
head(trt_data$resp)
#> item_1 item_2 item_3 item_4 item_5 item_6 item_7 item_8 item_9 item_10
#> 1 0 1 1 0 0 1 0 0 0 0
#> 2 1 1 1 1 1 1 0 1 0 1
#> 3 1 0 1 1 1 0 0 0 0 0
#> 4 1 0 0 0 1 1 0 0 0 0
#> 5 0 1 1 0 0 0 0 1 0 1
#> 6 0 0 1 0 1 0 1 1 1 13. IRT Estimation
We can estimate parameters using binary_irt (for dichotomous data), polytomous_irt (for graded data), or mixed_irt (for both).
Here, we estimate the parameters for the data we just simulated.
# Estimate 2PL model
# Note: max_iter set low for vignette speed; use default (100) in practice
fit_2pl <- binary_irt(
data = sim_data$resp,
model = "2PL",
method = "EM",
control = list(max_iter = 20, verbose = FALSE)
)
# View estimated item parameters
print(head(fit_2pl$item_params))
#> item discrimination discrimination_se difficulty difficulty_se number
#> 1 item_1 0.926 0.114 -0.063 0.105 500
#> 2 item_2 0.986 0.117 0.018 0.100 500
#> 3 item_3 0.840 0.110 -0.101 0.115 500
#> 4 item_4 0.983 0.117 -0.061 0.100 500
#> 5 item_5 1.036 0.119 0.111 0.096 500
#> 6 item_6 0.914 0.113 0.029 0.106 500
#> pvalue
#> 1 0.512
#> 2 0.496
#> 3 0.518
#> 4 0.512
#> 5 0.476
#> 6 0.494
# View model fit indices
print(fit_2pl$model_fit)
#> Index Value
#> 1 LogLikelihood -4219.624
#> 2 AIC 8499.247
#> 3 BIC 8625.685
#> 4 Iterations 13.0004. Testlet (TRT) Estimation
When items are grouped into testlets (e.g., reading passages), standard IRT assumptions are violated. tirt allows you to model this dependency using Bi-factor or Testlet models. This package allows estimations for binary and polytomous testlet models
We use trt_binary for this example. We must define the group structure.
# Items 1-5 are Testlet 1, Items 6-10 are Testlet 2
testlet_structure <- list(
c(1, 2, 3, 4, 5), # Indices for Testlet 1
c(6, 7, 8, 9, 10) # Indices for Testlet 2
)
fit_trt <- trt_binary(
data = trt_data$resp,
group = testlet_structure,
model = "2PLT",
method = "EM",
control = list(max_iter = 15, verbose = FALSE)
)
# The output includes item parameters
head(fit_trt$item_params)
#> item discrimination discrimination_se difficulty difficulty_se guessing
#> 1 item_1 1.887640 0.1768678 0.18400878 0.06118852 0
#> 2 item_2 1.030662 0.1142394 -0.06652782 0.09738936 0
#> 3 item_3 1.824778 0.1713166 0.02173113 0.06217262 0
#> 4 item_4 1.775003 0.1675142 -0.12431961 0.06403315 0
#> 5 item_5 1.823515 0.1713748 0.10875628 0.06233405 0
#> 6 item_6 1.636682 0.1541821 0.03051786 0.06822841 0
#> guessing_se
#> 1 0
#> 2 0
#> 3 0
#> 4 0
#> 5 0
#> 6 0
# And person parameters including the testlet specific effects (Gamma)
head(fit_trt$person_params)
#> person ability ability_se testlet_1 testlet_1_se testlet_2
#> 1 1 -0.4179722 0.6295140 0.001012444 0.6714548 -0.3640896
#> 2 2 0.5696989 0.6508913 0.894816372 0.7296370 -0.4251104
#> 3 3 -0.1633864 0.6509442 0.951443765 0.6901421 -1.1210697
#> 4 4 -0.3173443 0.6253186 0.146364703 0.6644582 -0.4262548
#> 5 5 -0.3581111 0.6258050 -0.040190281 0.6693177 -0.2694689
#> 6 6 0.2154850 0.6239477 -0.249288095 0.6626316 0.4437022
#> testlet_2_se
#> 1 0.6957473
#> 2 0.6885031
#> 3 0.7306446
#> 4 0.6925863
#> 5 0.6879875
#> 6 0.69111335. Universal Estimation (Mixed standaline & testlet items)
Real-world exams often contain a mix of independent items (i.e., using standard IRT) and item bundles/testlets (i.e., using TRT). The irt_trt function is a universal estimator that handles both simultaneously. It can also handle mixed data types (binary and polytomous) in the same run.
To use irt_trt, you must provide an item_spec data frame that maps each item to its model and testlet ID (if any).
5.1. Step 1: Simulate Mixed Structure Data
We will generate data for 1000 examinees with:
10 Independent Binary Items (2PL Model)
5 Independent Polytomous Items (GRM Model)
2 Testlet Items belonging to Testlet “T1” (2PLT Model)
3 Testlet Items belonging to Testlet “T2” (GPCMT Model)
cat("\n\n--- Generating Simulated Data ---\n")
#>
#>
#> --- Generating Simulated Data ---
set.seed(2025)
N <- 1000
# 20 Items Total:
# 1-10: Binary (2PL) - Independent
# 11-15: Polytomous (GRM, 3 cats) - Independent
# 16-17: Binary (2PLT) - Testlet 1
# 18-20: Polytomous (GPCT, 3 cats) - Testlet 2
theta <- rnorm(N, 0, 1)
gamma_1 <- rnorm(N, 0, 0.5) # Testlet 1 effect
gamma_2 <- rnorm(N, 0, 0.6) # Testlet 2 effect
resp_matrix <- matrix(NA, N, 20)
item_names <- paste0("Item_", 1:20)
colnames(resp_matrix) <- item_names
# Define Item Parameters (True Values)
a_true <- runif(20, 0.8, 1.5)
b_true <- seq(-1.5, 1.5, length.out = 20)
# 1. Simulate Binary Independent (1-10)
for(j in 1:10) {
p <- 1 / (1 + exp(-a_true[j] * (theta - b_true[j])))
resp_matrix[,j] <- rbinom(N, 1, p)
}
# 2. Simulate Poly Independent (11-15) - GR
for(j in 11:15) {
b_k <- sort(c(b_true[j] - 0.7, b_true[j] + 0.7))
p1 <- 1 / (1 + exp(-a_true[j] * (theta - b_k[1])))
p2 <- 1 / (1 + exp(-a_true[j] * (theta - b_k[2])))
# GRM Probabilities
prob_0 <- 1 - p1
prob_1 <- p1 - p2
prob_2 <- p2
# Sample
r <- runif(N)
resp_matrix[,j] <- ifelse(r < prob_0, 0, ifelse(r < prob_0 + prob_1, 1, 2))
}
# 3. Simulate Binary Testlet 1 (16-17)
for(j in 16:17) {
eff_theta <- theta + gamma_1
p <- 1 / (1 + exp(-a_true[j] * (eff_theta - b_true[j])))
resp_matrix[,j] <- rbinom(N, 1, p)
}
# 4. Simulate Poly Testlet 2 (18-20) - GPCMT
for(j in 18:20) {
eff_theta <- theta + gamma_2
# GPCM steps (centered at b_true)
b_vec <- c(b_true[j] - 0.5, b_true[j] + 0.5)
numer <- matrix(0, N, 3)
numer[,1] <- 0
numer[,2] <- a_true[j] * (eff_theta - b_vec[1])
numer[,3] <- a_true[j] * (eff_theta - b_vec[1]) + a_true[j] * (eff_theta - b_vec[2])
max_n <- apply(numer, 1, max)
exps <- exp(numer - max_n)
probs <- exps / rowSums(exps)
r <- runif(N)
resp_matrix[,j] <- apply(cbind(r, probs), 1, function(x) {
if(x[1] < x[2]) return(0)
if(x[1] < x[2]+x[3]) return(1)
return(2)
})
}
sim_data <- as.data.frame(resp_matrix)
#view the data
print(head(sim_data))
#> Item_1 Item_2 Item_3 Item_4 Item_5 Item_6 Item_7 Item_8 Item_9 Item_10
#> 1 1 1 1 1 1 1 1 1 1 0
#> 2 1 1 1 1 0 0 1 1 0 1
#> 3 1 1 1 1 1 1 1 1 0 1
#> 4 1 1 1 1 1 1 1 1 1 1
#> 5 1 1 1 1 1 0 0 1 0 0
#> 6 1 1 1 1 1 0 0 0 1 0
#> Item_11 Item_12 Item_13 Item_14 Item_15 Item_16 Item_17 Item_18 Item_19
#> 1 1 2 1 2 1 0 1 1 0
#> 2 1 2 2 0 2 0 0 2 1
#> 3 1 1 0 0 0 1 0 1 0
#> 4 1 2 2 2 1 0 0 2 2
#> 5 0 1 0 2 0 0 0 1 1
#> 6 2 1 1 0 0 0 0 0 0
#> Item_20
#> 1 1
#> 2 1
#> 3 1
#> 4 2
#> 5 1
#> 6 05.2. Step 2: Create the Item Specification
The item_spec data frame tells the function how to treat each item.
item: Must match column names in data.
model: The model to use (e.g., “2PL”, “2PLT”, “GPCMT”, “RaschT”, etc.).
testlet: The ID of the testlet. Use NA for independent items.
# Create Item Specification
spec <- data.frame(
item = item_names,
model = c(rep("2PL", 10), rep("GRM", 5), rep("2PLT", 2), rep("GPCMT", 3)),
testlet = c(rep(NA, 15), rep("testlet1", 2), rep("testlet2", 3)),
stringsAsFactors = FALSE
)
#view the spec
print(spec)
#> item 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 2PL <NA>
#> 8 Item_8 2PL <NA>
#> 9 Item_9 2PL <NA>
#> 10 Item_10 2PL <NA>
#> 11 Item_11 GRM <NA>
#> 12 Item_12 GRM <NA>
#> 13 Item_13 GRM <NA>
#> 14 Item_14 GRM <NA>
#> 15 Item_15 GRM <NA>
#> 16 Item_16 2PLT testlet1
#> 17 Item_17 2PLT testlet1
#> 18 Item_18 GPCMT testlet2
#> 19 Item_19 GPCMT testlet2
#> 20 Item_20 GPCMT testlet25.3. Step 3: Run Estimation
We run the estimation. The function automatically detects that “T1” items share a nuisance dimension (Gamma) while the others depend only on Theta.
# Run the Universal Function
results <- irt_trt(sim_data, spec, method = "EM", control = list(max_iter=50, verbose=TRUE))
#>
#> ================================================
#> UNIVERSAL IRT/TRT ESTIMATION
#> ================================================
#> Data Summary:
#> - Total Items: 20
#> - Binary Items: 12
#> - Polytomous Items: 8
#>
#> Structure:
#> - Independent Items: 15
#> - Testlet Items: 5 (grouped into 2 testlets)
#>
#> NOTE on Testlets:
#> Any item with a value in the 'testlet' string within the item_spec statement will be estimated
#> as a Testlet item (Theta + Gamma), regardless of the model string.
#> (e.g., model='GPCM' + testlet='T1' -> Testlet GPCT).
#> ------------------------------------------------
#> Data: 1000 Persons, 20 Items, 2 Testlets
#> ------------------------------------------------
#> Starting Estimation Loop...
#> Iter 1 | LogLik: -14918.96 | Max Change: 1.25812 Iter 2 | LogLik: -14004.14 | Max Change: 0.14006 Iter 3 | LogLik: -13988.50 | Max Change: 0.05962 Iter 4 | LogLik: -13979.68 | Max Change: 0.03897 Iter 5 | LogLik: -13973.19 | Max Change: 0.02922 Iter 6 | LogLik: -13968.36 | Max Change: 0.02305 Iter 7 | LogLik: -13964.75 | Max Change: 0.01837 Iter 8 | LogLik: -13962.04 | Max Change: 0.01476 Iter 9 | LogLik: -13959.98 | Max Change: 0.01190 Iter 10 | LogLik: -13958.41 | Max Change: 0.00966 Iter 11 | LogLik: -13957.20 | Max Change: 0.00791 Iter 12 | LogLik: -13956.26 | Max Change: 0.00648 Iter 13 | LogLik: -13955.53 | Max Change: 0.00534 Iter 14 | LogLik: -13954.95 | Max Change: 0.00440 Iter 15 | LogLik: -13954.50 | Max Change: 0.00361 Iter 16 | LogLik: -13954.13 | Max Change: 0.00298 Iter 17 | LogLik: -13953.84 | Max Change: 0.00245 Iter 18 | LogLik: -13953.61 | Max Change: 0.00202 Iter 19 | LogLik: -13953.42 | Max Change: 0.00168 Iter 20 | LogLik: -13953.27 | Max Change: 0.00139 Iter 21 | LogLik: -13953.15 | Max Change: 0.00114 Iter 22 | LogLik: -13953.04 | Max Change: 0.00095 Iter 23 | LogLik: -13952.96 | Max Change: 0.00076 Iter 24 | LogLik: -13952.90 | Max Change: 0.00065 Iter 25 | LogLik: -13952.84 | Max Change: 0.00054 Iter 26 | LogLik: -13952.80 | Max Change: 0.00045 Iter 27 | LogLik: -13952.76 | Max Change: 0.00038 Iter 28 | LogLik: -13952.73 | Max Change: 0.00030 Iter 29 | LogLik: -13952.71 | Max Change: 0.00025 Iter 30 | LogLik: -13952.69 | Max Change: 0.00021 Iter 31 | LogLik: -13952.67 | Max Change: 0.00018 Iter 32 | LogLik: -13952.66 | Max Change: 0.00018 Iter 33 | LogLik: -13952.64 | Max Change: 0.00015 Iter 34 | LogLik: -13952.64 | Max Change: 0.00013 Iter 35 | LogLik: -13952.63 | Max Change: 0.00010
#>
#> >>> Convergence Confirmation: Model Converged!
#> ------------------------------------------------
#> Estimating Final Person Parameters (Theta)...
#> Using Theta range: [-4.0, 4.0]. Adjust via control=list(theta_range=c(...)).
#> Estimating Testlet Effects (Gamma) per person...
#> All parameters estimated. Formatting results...
#> All completed. You can view the results now.
#> ================================================
# Save Results
item_par=results$item_params
fit=results$model_fit
person_par=results$person_params6. Fixed Parameter Calibration
In equating scenarios or computer-adaptive testing, we often need to estimate one set of parameters while holding the other fixed.
6.1 Fixed Item Calibration (Estimating Person Ability)
If we know the item parameters (e.g., from an item bank), we can estimate person ability. Below we first simulate items, including some items have “known parameters” and some items are “new items”.
# --- 1. Simulation Helper Functions ---
# Function to simulate Binary 2PL
sim_2pl <- function(theta, a, b) {
D <- 1.7
prob <- 1 / (1 + exp(-D * a * (theta - b)))
resp <- ifelse(runif(length(theta)) < prob, 1, 0)
return(resp)
}
# Function to simulate Polytomous (GPCM structure to match fixed_item logic)
sim_poly <- function(theta, a, steps) {
n_cat <- length(steps) + 1
n_stud <- length(theta)
# Calculate numerators for categories 0, 1, 2...
# P(k) proportional to exp( sum(a*(theta - step_j)) )
probs <- matrix(0, nrow = n_stud, ncol = n_cat)
probs[, 1] <- 1 # Reference category (unnormalized exp is 1 or e^0)
current_sum <- 0
for (k in 2:n_cat) {
# Step parameters correspond to boundaries
current_sum <- current_sum + a * (theta - steps[k-1])
probs[, k] <- exp(current_sum)
}
# Normalize
prob_sum <- rowSums(probs)
probs_final <- probs / prob_sum
# Random draw
resp <- numeric(n_stud)
for (i in 1:n_stud) {
resp[i] <- sample(0:(n_cat-1), 1, prob = probs_final[i, ])
}
return(resp)
}
# --- 2. Generate Parameters for 50 Items ---
n_students <- 5497
true_theta <- rnorm(n_students, mean = 0.4, sd = 1.7)
all_items_meta <- list()
response_matrix <- matrix(NA, nrow = n_students, ncol = 50)
colnames(response_matrix) <- paste0("Item_", 1:50)
# --- Set 1: Known Dichotomous (Items 1-10) ---
# Model: 2PL
for (i in 1:10) {
a <- runif(1, 0.8, 1.4)
b <- rnorm(1, 0, 1)
resp <- sim_2pl(true_theta, a, b)
response_matrix[, i] <- resp
all_items_meta[[i]] <- list(item = colnames(response_matrix)[i], type = "Known", model = "2PL", a = a, b = b)
}
# --- Set 2: Known Polytomous (Items 11-18) ---
# Model: Poly (3 categories: 0, 1, 2). Requires 2 steps (d1, d2).
for (i in 11:18) {
a <- runif(1, 0.7, 1.2)
d1 <- rnorm(1, -0.5, 0.3)
d2 <- rnorm(1, 0.5, 0.3) # Steps usually ordered
resp <- sim_poly(true_theta, a, c(d1, d2))
response_matrix[, i] <- resp
all_items_meta[[i]] <- list(item = colnames(response_matrix)[i], type = "Known", model = "Poly", a = a, d1 = d1, d2 = d2)
}
# --- Set 3: Unknown Dichotomous (Items 19-40) ---
# 22 Items
for (i in 19:40) {
a <- runif(1, 0.8, 1.4)
b <- rnorm(1, 0, 1)
resp <- sim_2pl(true_theta, a, b)
response_matrix[, i] <- resp
all_items_meta[[i]] <- list(item = colnames(response_matrix)[i], type = "Unknown", model = "2PL", a = a, b = b)
}
# --- Set 4: Unknown Polytomous (Items 41-50) ---
# 10 Items (Varying categories: some 3 cat, some 4 cat)
for (i in 41:50) {
a <- runif(1, 0.7, 1.2)
# Mix of 3 categories (2 steps) and 4 categories (3 steps)
if (i %% 2 == 0) {
steps <- c(-0.8, 0, 0.8) # 4 cats
} else {
steps <- c(-0.5, 0.5) # 3 cats
}
resp <- sim_poly(true_theta, a, steps)
response_matrix[, i] <- resp
all_items_meta[[i]] <- list(item = colnames(response_matrix)[i], type = "Unknown", model = "Poly", a = a, steps = steps)
}
#check
response_df <- as.data.frame(response_matrix)
head(response_df[, c(1:3, 11:13, 41:43)]) # Peak at data
#> Item_1 Item_2 Item_3 Item_11 Item_12 Item_13 Item_41 Item_42 Item_43
#> 1 1 0 1 0 2 1 0 0 0
#> 2 0 0 0 0 0 0 1 0 1
#> 3 1 0 1 1 1 1 2 1 2
#> 4 0 0 0 0 0 0 1 0 0
#> 5 1 0 0 1 2 2 1 1 1
#> 6 1 0 0 0 1 2 1 1 0Then we create a necessary data structure for the next calibration step. In this step, we need to specify known item with corresponding parameters.
# Create empty dataframe structure
known_params_df <- data.frame(
item = character(),
model = character(),
a = numeric(),
b = numeric(),
d1 = numeric(),
d2 = numeric(),
stringsAsFactors = FALSE
)
# Populate with Known Items (1-18)
for (i in 1:18) {
meta <- all_items_meta[[i]]
if (meta$model == "2PL") {
# Add row for Binary
known_params_df <- rbind(known_params_df, data.frame(
item = meta$item,
model = "2PL",
a = meta$a,
b = meta$b,
d1 = NA, d2 = NA
))
} else {
# Add row for Poly
# Note: Using "GPCM" model string as our function uses GPCM math
known_params_df <- rbind(known_params_df, data.frame(
item = meta$item,
model = "GPCM",
a = meta$a,
b = NA,
d1 = meta$d1,
d2 = meta$d2
))
}
}
# Display input to verify
print("Known Item Parameters Input:")
#> [1] "Known Item Parameters Input:"
print(known_params_df)
#> item model a b d1 d2
#> 1 Item_1 2PL 1.0287464 -0.9197240 NA NA
#> 2 Item_2 2PL 1.0048075 1.2354626 NA NA
#> 3 Item_3 2PL 0.9183354 0.6649479 NA NA
#> 4 Item_4 2PL 1.2716232 -1.1443195 NA NA
#> 5 Item_5 2PL 0.9924154 0.7640520 NA NA
#> 6 Item_6 2PL 1.1279951 0.9675910 NA NA
#> 7 Item_7 2PL 1.2603260 0.2995644 NA NA
#> 8 Item_8 2PL 1.2943065 1.2325874 NA NA
#> 9 Item_9 2PL 1.3738451 -1.7190245 NA NA
#> 10 Item_10 2PL 1.0130638 -1.1422731 NA NA
#> 11 Item_11 GPCM 1.1950806 NA 0.07845972 0.74138298
#> 12 Item_12 GPCM 1.1169482 NA -0.50576871 0.32892965
#> 13 Item_13 GPCM 1.1917066 NA -0.27650908 0.51124636
#> 14 Item_14 GPCM 0.9595820 NA -0.37947396 0.22564038
#> 15 Item_15 GPCM 0.9219541 NA -0.14945831 0.44633887
#> 16 Item_16 GPCM 0.9264782 NA -0.59867093 0.03966326
#> 17 Item_17 GPCM 0.9225556 NA -0.70413230 0.71069682
#> 18 Item_18 GPCM 1.0801744 NA -0.09153471 1.07061113Finally, we can conduct a fiex-item calibration and check the final results.
# Run the calibration
# We set model_default = "2PL" (this applies to unknown binary items).
# The function automatically detects poly items in the unknown set based on response categories.
results <- fixed_item(
response_df = response_df,
item_params_df = known_params_df,
control = list(max_iter = 100)
)
#> ---------------------------------------------------------
#> Starting Fixed Item Calibration...
#> ---------------------------------------------------------
#> Inferred Unknown Binary: 2PL
#> Inferred Unknown Poly: GPCM
#> Cycle 1: LogLik = -144716.33, Change = Inf
#> Cycle 2: LogLik = -142757.17, Change = 1959.1556
#> Cycle 3: LogLik = -142546.08, Change = 211.0919
#> Cycle 4: LogLik = -142456.85, Change = 89.2311
#> Cycle 5: LogLik = -142416.96, Change = 39.8897
#> Cycle 6: LogLik = -142398.92, Change = 18.0392
#> Cycle 7: LogLik = -142390.71, Change = 8.2095
#> Cycle 8: LogLik = -142386.96, Change = 3.7538
#> Cycle 9: LogLik = -142385.24, Change = 1.7196
#> Cycle 10: LogLik = -142384.45, Change = 0.7895
#> Cycle 11: LogLik = -142384.08, Change = 0.3636
#> Cycle 12: LogLik = -142383.92, Change = 0.1672
#> Cycle 13: LogLik = -142383.84, Change = 0.0768
#> Cycle 14: LogLik = -142383.80, Change = 0.0357
#> Cycle 15: LogLik = -142383.79, Change = 0.0165
#> Cycle 16: LogLik = -142383.78, Change = 0.0075
#> Cycle 17: LogLik = -142383.78, Change = 0.0035
#> Estimating Person Parameters...
#> ---------------------------------------------------------
#> Process Complete. You can view results now.
# Check that Items 1-18 are "Fixed" and match inputs.
# Check that Items 19-50 are "Estimated".
# Note how Poly items have step_1, step_2, etc.
item=results$item_params
person=results$person_params
fit=results$model_fit6.2 Fixed Person Calibration (Estimating Item Difficulty)
Conversely, if we have known person abilities, we can calibrate new items. In this estimation, you can also choose to model a covariate if you believe that is releated to your latent variables. If not specified, then there will be no covariate in the modeling process.
The final output also contains some classic item statistics such as p-values and mean of item scores.
# --- Example: With Package Data ---
data("ela1", package = "tirt")
# Select Item Responses (Cols 1-30)
df_real <- ela1[, 1:30]
fixed_theta <- ela1$THETA
fixed_cov <- ela1$COVARIATE
real_res <- fix_person(df = df_real,
theta = fixed_theta,
model = "2PL",
covariate = fixed_cov)
head(real_res)
#> item discrimination discrimination_se discrimination_zvalue
#> 1 ITEM1 0.08298835 0.004848405 17.11663
#> 2 ITEM2 0.08620777 0.005016701 17.18416
#> 3 ITEM3 0.09418605 0.004894558 19.24301
#> 4 ITEM4 0.18985935 0.009693893 19.58546
#> 5 ITEM5 0.21869224 0.011528756 18.96928
#> 6 ITEM6 0.03696713 0.006730610 5.49239
#> discrimination_Pr_z difficulty difficulty_se difficulty_zvalue
#> 1 1.115533e-65 -0.8287487 0.1058117 -7.832297
#> 2 3.489877e-66 -7.4745653 0.3911856 -19.107464
#> 3 1.615391e-82 -3.9376999 0.1807591 -21.784239
#> 4 2.057399e-85 -13.2792707 0.6164795 -21.540489
#> 5 3.060479e-80 -13.1473999 0.6265239 -20.984675
#> 6 3.965314e-08 -25.0506752 4.4253213 -5.660759
#> difficulty_Pr_z covariate covariate_se covariate_zvalue covariate_Pr_z
#> 1 4.884981e-15 -0.001330397 0.0003265024 -4.074693 4.60751e-05
#> 2 0.000000e+00 -0.001330397 0.0003265024 -4.074693 4.60751e-05
#> 3 0.000000e+00 -0.001330397 0.0003265024 -4.074693 4.60751e-05
#> 4 0.000000e+00 -0.001330397 0.0003265024 -4.074693 4.60751e-05
#> 5 0.000000e+00 -0.001330397 0.0003265024 -4.074693 4.60751e-05
#> 6 1.507053e-08 -0.001330397 0.0003265024 -4.074693 4.60751e-05
#> pvalue number correlation mean_theta mean_covariate
#> 1 0.4983116 52417 0.07441120 -0.9305780 -1.196020
#> 2 0.6370262 52417 0.07473806 -0.9305780 -1.196020
#> 3 0.5702539 52417 0.08380568 -0.9305780 -1.196020
#> 4 0.9100940 39063 0.09939991 -0.8452539 -1.197829
#> 5 0.9331212 36454 0.09975399 -0.8022383 -1.185759
#> 6 0.7106403 34642 0.02892981 -0.7638062 -1.1679997. Equating / Linking
equate_irt links two scales (Form X and Form Y) using common anchor items.
# Create dummy parameters for Form X (Base)
base_params <- data.frame(
item = c("Item1", "Item2", "Item3", "Item4", "Item5"),
model = "2PL",
a = c(1.0, 1.2, 0.9, 1.1, 1.0),
b = c(-1.0, -0.5, 0.0, 0.5, 1.0),
stringsAsFactors = FALSE
)
# Create dummy parameters for Form Y (New)
# Suppose Form Y is harder and has different discrimination scaling
new_params <- data.frame(
item = c("Item1", "Item2", "Item3", "Item6", "Item7"), # Items 1-3 are anchors
model = "2PL",
a = c(1.0, 1.2, 0.9, 1.1, 1.0) / 0.9, # Scale shift
b = (c(-1.0, -0.5, 0.0, 0.8, 1.2) * 0.9) + 0.2, # Location shift
stringsAsFactors = FALSE
)
# Perform Equating
linked <- equate_irt(
base_params = base_params,
new_params = new_params,
methods = c("Stocking-Lord", "Mean-Mean")
)
#> ---------------------------------------------------------
#> Starting IRT Equating / Linking
#> ---------------------------------------------------------
#> Methods selected: Stocking-Lord, Mean-Mean
#> Models detected: 2PL
#> Items detected: Base = 5 | New = 5
#> Anchor items detected: 3
#> Anchor IDs: Item1, Item2, Item3
#> Estimating constants: Mean-Mean...
#> Estimating constants: Stocking-Lord (TCC Optimization)...
#> Transforming parameters...
#> Equating completed.
#> ---------------------------------------------------------
# View Linking Constants (A and B)
print(linked$linking_constants)
#> Method A B
#> 1 Mean-Mean 0.900000 -0.2750000
#> 2 Stocking-Lord 1.111111 -0.2222221
# View Transformed Parameters for Form Y (put on Form X scale)
print(linked$transformed_item_params)
#> item anchor_status method model discrimination discrimination_se
#> 1 Item1 Anchor Mean-Mean 2PL 1.2345679 NA
#> 2 Item2 Anchor Mean-Mean 2PL 1.4814815 NA
#> 3 Item3 Anchor Mean-Mean 2PL 1.1111111 NA
#> 4 Item6 New Mean-Mean 2PL 1.3580247 NA
#> 5 Item7 New Mean-Mean 2PL 1.2345679 NA
#> 6 Item1 Anchor Stocking-Lord 2PL 0.9999997 NA
#> 7 Item2 Anchor Stocking-Lord 2PL 1.1999997 NA
#> 8 Item3 Anchor Stocking-Lord 2PL 0.8999998 NA
#> 9 Item6 New Stocking-Lord 2PL 1.0999997 NA
#> 10 Item7 New Stocking-Lord 2PL 0.9999997 NA
#> difficulty difficulty_se guessing guessing_se
#> 1 -9.050000e-01 NA 0 NA
#> 2 -5.000000e-01 NA 0 NA
#> 3 -9.500000e-02 NA 0 NA
#> 4 5.530000e-01 NA 0 NA
#> 5 8.770000e-01 NA 0 NA
#> 6 -1.000000e+00 NA 0 NA
#> 7 -5.000000e-01 NA 0 NA
#> 8 1.380493e-07 NA 0 NA
#> 9 8.000003e-01 NA 0 NA
#> 10 1.200000e+00 NA 0 NA8. Final Comment
The tirt package streamlines complex psychometric analyses. For detailed documentation on specific functions, please use the help command, e.g., ?binary_irt.