Esto es una replicación del tutorial de Philipp K. Masur sobre IRT en R, disponible en: https://philippmasur.de/2022/05/13/how-to-run-irt-analyses-in-r/, pero únicamente del modelo 3pl.
En este breve tutorial introductorio, nos centramos en tres de los modelos IRT más utilizados: el modelo 3PL, el modelo 2PL y, finalmente, el modelo 1PL o Rasch. Estos modelos se nombran por el número de parámetros de elementos utilizados en la función que modela la relación entre θ y la respuesta del ítem (0/1). Cada modelo tiene propiedades únicas, pero todos ellos son adecuados para estimar variables latentes a partir de elementos binarios (por ejemplo, pruebas de conocimiento), que trataremos en este tutorial.
library(readxl)
library(tidyverse) #esto incluye 8 paquetes: ggplot, tibble, tidyr, readr, purr, dplyr, stringr y forcats
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.2 ──
## ✔ ggplot2 3.4.1 ✔ purrr 1.0.1
## ✔ tibble 3.2.1 ✔ dplyr 1.1.1
## ✔ tidyr 1.3.0 ✔ stringr 1.5.0
## ✔ readr 2.1.3 ✔ forcats 1.0.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
library(sjPlot)
## Learn more about sjPlot with 'browseVignettes("sjPlot")'.
library(psych)
##
## Attaching package: 'psych'
##
## The following objects are masked from 'package:ggplot2':
##
## %+%, alpha
library(mirt)
## Loading required package: stats4
## Loading required package: lattice
# El paqueteggmirt incluye una función conveniente para simular datos para análisis IRT. Vamos a crear rápidamente un conjunto de datos con 500 observaciones y 10 elementos que deberían ajustarse a un modelo 3PL, 2PL y quizás incluso 1PL.
library(ggmirt)
library(gtable)
#definir la base de datos
#esta se puede cambiar
set.seed(42)
d <- sim_irt(500, 10, discrimination = .25, seed = 42)
## Warning: The `x` argument of `as_tibble.matrix()` must have unique column names if
## `.name_repair` is omitted as of tibble 2.0.0.
## ℹ Using compatibility `.name_repair`.
## ℹ The deprecated feature was likely used in the ggmirt package.
## Please report the issue to the authors.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
head(d)
## # A tibble: 6 × 10
## V1 V2 V3 V4 V5 V6 V7 V8 V9 V10
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0 1 0 1 0 1 1 1 0 1
## 2 0 1 1 0 0 0 0 0 0 1
## 3 0 1 0 1 0 1 0 1 0 0
## 4 0 0 1 0 0 1 0 0 0 1
## 5 0 1 1 1 0 0 1 1 1 0
## 6 0 0 1 0 0 0 0 1 0 1
#modelo 3pl
unimodel <- 'F1 = 1-10'
fit3PL <- mirt(data = d,
model = 1, # alternatively, we could also just specify model = 1 in this case
itemtype = "3PL",
verbose = FALSE)
# en IRT, generalmente estamos más interesados en los parámetros reales de IRT como se discutió anteriormente (discriminación, dificultad y probabilidad de adivinar). Se pueden extraer de la siguiente manera:
params3PL <- coef(fit3PL, IRTpars = TRUE, simplify = TRUE)
round(params3PL$items, 2) # g = c = guessing parameter
## a b g u
## V1 1.70 1.17 0.00 1
## V2 1.02 -0.58 0.00 1
## V3 1.08 0.35 0.00 1
## V4 1.80 0.81 0.09 1
## V5 1.19 0.52 0.00 1
## V6 1.17 -0.14 0.00 1
## V7 1.29 1.87 0.00 1
## V8 0.84 0.03 0.00 1
## V9 1.76 2.11 0.00 1
## V10 0.93 -0.05 0.01 1
# Los valores de los parámetros de la pendiente (parámetros a), son una medida de qué tan bien un elemento diferencia a los individuos con diferentes niveles theta. Los valores más grandes o las pendientes más pronunciadas son mejores para diferenciar a las personas. Una pendiente también puede interpretarse como un indicador de la fuerza de una relación entre un elemento y un rasgo latente, con valores de pendiente más altos que corresponden a relaciones más fuertes. Los parámetros de ubicación o dificultad (parámetro b) también se enumeran para cada elemento. Los parámetros de ubicación se interpretan como el valor de theta que corresponde a una probabilidad de .50 de responder correctamente en o por encima de esa ubicación en un elemento. Los parámetros de ubicación muestran que los elementos cubren una amplia gama del rasgo latente.
# Los parámetros de ubicación o dificultad (parámetro b) también se enumeran para cada ítem. Los parámetros de ubicación se interpretan como el valor de theta que corresponde a una probabilidad de .50 de responder correctamente en o por encima de esa ubicación en un elemento. Los parámetros de ubicación muestran que los elementos cubren una amplia gama del rasgo latente.
#Ahora calculamos el índice M2, que está diseñado para evaluar el ajuste de los modelos de respuesta al ítem.
M2(fit3PL)
## M2 df p RMSEA RMSEA_5 RMSEA_95 SRMSR TLI CFI
## stats 17.77133 25 0.8519424 0 0 0.02036034 0.03371296 1.018438 1
# Como podemos ver, la estadística M2 es comparativamente baja y no significativa. Por lo tanto, no hay diferencias preocupantes entre el modelo y los datos. Esto está respaldado además por un RMSEA muy bajo y un CFA y TLI de 1. En la TRI, sin embargo, por lo general estamos más interesados en los índices de ajuste de ítems y personas. IRT nos permite evaluar qué tan bien se ajusta cada elemento al modelo y si los patrones de respuesta individuales se alinean con el modelo.
#Comencemos con el ajuste de los elementos: al igual que en muchas otras áreas, se han propuesto diferentes índices para evaluar el ajuste de los elementos. Podemos usar la función itemfit() para obtener una variedad de ellos. Por defecto, recibimos el S_X2 de Orlando y Thissen (2000) y los valores dfs, RMSEA y p correspondientes. Esta prueba no debe ser significativa para indicar un buen ajuste del ítem. Como podemos ver aquí, solo el artículo V9 muestra un menor ajuste con el modelo.
#En términos generales, los valores no estandarizados deben estar entre 0,5 y 1,5 para que no sean degradantes. En el paquete ggmirt, también podemos usar la función ´itemfitPlot()` para inspeccionar esto visualmente.
itemfit(fit3PL)
## item S_X2 df.S_X2 RMSEA.S_X2 p.S_X2
## 1 V1 4.709 4 0.019 0.319
## 2 V2 2.503 5 0.000 0.776
## 3 V3 6.762 5 0.027 0.239
## 4 V4 3.701 5 0.000 0.593
## 5 V5 4.093 5 0.000 0.536
## 6 V6 2.824 5 0.000 0.727
## 7 V7 8.648 5 0.038 0.124
## 8 V8 2.901 5 0.000 0.715
## 9 V9 11.610 4 0.062 0.020
## 10 V10 2.694 5 0.000 0.747
itemfit(fit3PL, fit_stats = "infit") # typical for Rasch modeling
## item outfit z.outfit infit z.infit
## 1 V1 0.639 -2.365 0.874 -1.736
## 2 V2 0.848 -3.059 0.871 -3.713
## 3 V3 0.832 -3.792 0.879 -3.444
## 4 V4 0.836 -2.533 0.873 -2.566
## 5 V5 0.804 -3.634 0.873 -3.186
## 6 V6 0.810 -4.159 0.843 -4.604
## 7 V7 0.731 -1.594 0.994 -0.033
## 8 V8 0.893 -3.420 0.911 -3.296
## 9 V9 0.526 -1.359 1.053 0.383
## 10 V10 0.871 -3.689 0.895 -3.644
itemfitPlot(fit3PL)
## `geom_line()`: Each group consists of only one observation.
## ℹ Do you need to adjust the group aesthetic?
## `geom_line()`: Each group consists of only one observation.
## ℹ Do you need to adjust the group aesthetic?
# Evaluando ajuste personal
#Técnicamente, podemos producir exactamente las mismas medidas para cada persona para evaluar qué tan bien se alinean los patrones de respuesta de cada persona con este modelo. Piénselo de esta manera: si las respuestas de una persona con un theta alto (es decir, una capacidad latente alta) son incorrectas para un ítem fácil, esta persona no encaja en el modelo. Por el contrario, si una persona con poca capacidad responde correctamente a una pregunta muy difícil, tampoco encaja en el modelo. En la práctica, lo más probable es que haya algunas personas que no encajen bien en el modelo. Pero mientras el número de encuestados no aptos sea bajo, estamos bien. En su mayoría, volvemos a mirar las estadísticas de infit y outfit. Si menos del 5% de los encuestados tienen valores de infit y outfit superiores o inferiores a 1,96 y -1,96, estamos bien.
head(personfit(fit3PL))
## outfit z.outfit infit z.infit Zh
## 1 1.0161479 0.17495415 1.0183924 0.15840392 -0.07323077
## 2 0.6331922 -0.13275101 0.8411215 -0.49761279 0.56159341
## 3 0.7456650 -0.19057356 0.8912849 -0.36118534 0.43906570
## 4 0.7209540 -0.03928867 0.9571507 -0.06374702 0.26012606
## 5 2.1363674 2.36953932 1.7713727 2.18326557 -2.69486327
## 6 0.7538809 0.06104147 1.0221023 0.17247403 0.08435612
personfit(fit3PL) %>%
summarize(infit.outside = prop.table(table(z.infit > 1.96 | z.infit < -1.96)),
outfit.outside = prop.table(table(z.outfit > 1.96 | z.outfit < -1.96))) # lower row = non-fitting people
## Warning: Returning more (or less) than 1 row per `summarise()` group was deprecated in
## dplyr 1.1.0.
## ℹ Please use `reframe()` instead.
## ℹ When switching from `summarise()` to `reframe()`, remember that `reframe()`
## always returns an ungrouped data frame and adjust accordingly.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## infit.outside outfit.outside
## 1 0.958 0.98
## 2 0.042 0.02
personfitPlot(fit3PL)
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
# Item Characteristics Curves (Trace Plots)
# Las curvas características de los elementos o los gráficos de trazas visualizan los tres parámetros IRT para todos los elementos. Esta visualización es útil para comprender mejor las propiedades únicas de cada elemento. También puede ser útil identificar lagunas en la evaluación, así como diferencias en la pendiente.
tracePlot(fit3PL)
tracePlot(fit3PL, facet = F, legend = T) + scale_color_brewer(palette = "Set3")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
tracePlot(fit3PL,items = c(1:3), facet = F, legend = T) + scale_color_brewer(palette = "Set2")
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
# Item Information Curves
# Otra forma de ver la calidad de cada elemento es trazando las llamadas curvas de información del elemento. La información es un concepto estadístico que se refiere a la capacidad de un elemento para estimar con precisión las puntuaciones en theta. La información a nivel de elemento aclara qué tan bien contribuye cada elemento a la precisión de la estimación de puntaje con niveles más altos de información que conducen a estimaciones de puntaje más precisas.
itemInfoPlot(fit3PL) + scale_color_brewer(palette = "Set3")
## Warning: The `<scale>` argument of `guides()` cannot be `FALSE`. Use "none" instead as
## of ggplot2 3.3.4.
## ℹ The deprecated feature was likely used in the ggmirt package.
## Please report the issue to the authors.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Scale for colour is already present.
## Adding another scale for colour, which will replace the existing scale.
itemInfoPlot(fit3PL, facet = T)
# Test Information Curves
#El concepto de “información” también se puede aplicar a toda la escala. Aquí, vemos que la escala es muy buena para estimar puntajes theta entre -2 y 3, pero tiene menos precisión para estimar puntajes theta en los extremos.
testInfoPlot(fit3PL, adj_factor = 2)
#Scale Characteristic Curves
#Una última cualidad de un modelo IRT puede ser que la puntuación correcta del número simple (puntuación de la suma de las respuestas correctas) sea una estimación suficientemente buena del rasgo subyacente. Un gráfico de la llamada curva característica de escala permite evaluar esto visualmente al graficar la relación entre theta y el puntaje numérico correcto.
scaleCharPlot(fit3PL)
# Esta curva generalmente toma la forma de una S ya que la relación es más fuerte en el rango medio de theta y peor en los extremos (como ya se vio en la curva de información de la prueba). Por supuesto, también podemos probar esto con una correlación simple. Primero, extraemos la puntuación IRT latente usando la función fscores(). Luego lo correlacionamos con el simple número-puntaje-correcto.
score <- fscores(fit3PL)
sumscore <- rowSums(d)
cor.test(score, sumscore)
##
## Pearson's product-moment correlation
##
## data: score and sumscore
## t = 181.54, df = 498, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9910997 0.9937302
## sample estimates:
## cor
## 0.9925295
#Como podemos ver, se correlacionan casi a la perfección.