Description Usage Arguments Details Value Author(s) References See Also Examples
mixedmirt
fits MIRT models using FIML estimation to dichotomous and polytomous
IRT models conditional on fixed and random effect of person and item level covariates.
This can also be understood as 'explanatory IRT' if only fixed effects are modeled, or
multilevel/mixed IRT if random and fixed effects are included. The method uses the MH-RM
algorithm exclusively. Additionally, computation of the log-likelihood can be sped up by
using parallel estimation via mirtCluster
.
1 2 3 4 5 |
data |
a |
covdata |
a |
model |
an object returned from, or a string to be passed to, |
fixed |
a right sided R formula for specifying the fixed effect (aka 'explanatory')
predictors from |
random |
a right sided formula or list of formulas containing crossed random effects
of the form |
itemtype |
same as itemtype in |
lr.fixed |
an R formula (or list of formulas) to specify regression
effects in the latent variables from the variables in |
lr.random |
a list of random effect terms for modeling variability in the
latent trait scores, where the syntax uses the same style as in the |
itemdesign |
a |
constrain |
a list indicating parameter equality constrains. See |
pars |
used for parameter starting values. See |
return.design |
logical; return the design matrices before they have (potentially) been reassigned? |
SE |
logical; compute the standard errors by approximating the information matrix using the MHRM algorithm? Default is TRUE |
internal_constraints |
logical; use the internally defined constraints for constraining effects across persons and items? Default is TRUE. Setting this to FALSE runs the risk of under-identification |
technical |
the technical list passed to the MH-RM estimation engine, with the
SEtol default increased to .0001. Additionally, the argument |
... |
additional arguments to be passed to the MH-RM estimation engine. See
|
For dichotomous response models, mixedmirt
follows the general form
P(x = 1|θ, ψ) = g + \frac{(u - g)}{1 + exp(-1 * [θ a + X β + Z δ])}
where X is a design matrix with associated β fixed effect intercept coefficients,
and Z is a design matrix with associated δ random effects for the intercepts.
For simplicity and easier interpretation, the unique item intercept values typically found in
X β
are extracted and reassigned within mirt's 'intercept' parameters (e.g., 'd'
).
To observe how the design matrices are structured prior to reassignment and estimation pass
the argument return.design = TRUE
.
Polytomous IRT models follow a similar format except the item intercepts are automatically
estimated internally, rendering the items
argument in the fixed formula redundant and
therefore must be omitted from the specification. If there are a mixture of dichotomous and
polytomous items the intercepts for the dichotomous models are also estimated for consistency.
The decomposition of the θ parameters is also possible to form
latent regression and multilevel IRT models by using the lr.fixed
and lr.random
inputs. These effects decompose θ such that
θ = V Γ + W ζ + ε
where V and W are fixed and random effects design matrices for the associated coefficients.
To simulate maximum a posteriori estimates for the random effect terms
use the randef
function.
function returns an object of class MixedClass
(MixedClass-class).
Phil Chalmers rphilip.chalmers@gmail.com
Chalmers, R., P. (2012). mirt: A Multidimensional Item Response Theory Package for the R Environment. Journal of Statistical Software, 48(6), 1-29. doi: 10.18637/jss.v048.i06
Chalmers, R. P. (2015). Extended Mixed-Effects Item Response Models with the MH-RM Algorithm. Journal of Educational Measurement, 52, 200-222. doi: 10.1111/jedm.12072
mirt
, randef
, fixef
, boot.mirt
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 | ## Not run:
#make some data
set.seed(1234)
N <- 750
a <- matrix(rlnorm(10,.3,1),10,1)
d <- matrix(rnorm(10), 10)
Theta <- matrix(sort(rnorm(N)))
pseudoIQ <- Theta * 5 + 100 + rnorm(N, 0 , 5)
pseudoIQ <- (pseudoIQ - mean(pseudoIQ))/10 #rescale variable for numerical stability
group <- factor(rep(c('G1','G2','G3'), each = N/3))
data <- simdata(a,d,N, itemtype = rep('2PL',10), Theta=Theta)
covdata <- data.frame(group, pseudoIQ)
#use parallel computing
mirtCluster()
#specify IRT model
model <- 'Theta = 1-10'
#model with no person predictors
mod0 <- mirt(data, model, itemtype = 'Rasch')
#group as a fixed effect predictor (aka, uniform dif)
mod1 <- mixedmirt(data, covdata, model, fixed = ~ 0 + group + items)
anova(mod0, mod1)
summary(mod1)
coef(mod1)
#same model as above in lme4
wide <- data.frame(id=1:nrow(data),data,covdata)
long <- reshape2::melt(wide, id.vars = c('id', 'group', 'pseudoIQ'))
library(lme4)
lmod0 <- glmer(value ~ 0 + variable + (1|id), long, family = binomial)
lmod1 <- glmer(value ~ 0 + group + variable + (1|id), long, family = binomial)
anova(lmod0, lmod1)
#model using 2PL items instead of Rasch
mod1b <- mixedmirt(data, covdata, model, fixed = ~ 0 + group + items, itemtype = '2PL')
anova(mod1, mod1b) #better with 2PL models using all criteria (as expected, given simdata pars)
#continuous predictor with group
mod2 <- mixedmirt(data, covdata, model, fixed = ~ 0 + group + items + pseudoIQ)
summary(mod2)
anova(mod1b, mod2)
#view fixed design matrix with and without unique item level intercepts
withint <- mixedmirt(data, covdata, model, fixed = ~ 0 + items + group, return.design = TRUE)
withoutint <- mixedmirt(data, covdata, model, fixed = ~ 0 + group, return.design = TRUE)
#notice that in result above, the intercepts 'items1 to items 10' were reassigned to 'd'
head(withint$X)
tail(withint$X)
head(withoutint$X) #no intercepts design here to be reassigned into item intercepts
tail(withoutint$X)
###################################################
### random effects
#make the number of groups much larger
covdata$group <- factor(rep(paste0('G',1:50), each = N/50))
#random groups
rmod1 <- mixedmirt(data, covdata, 1, fixed = ~ 0 + items, random = ~ 1|group)
summary(rmod1)
coef(rmod1)
#random groups and random items
rmod2 <- mixedmirt(data, covdata, 1, random = list(~ 1|group, ~ 1|items))
summary(rmod2)
eff <- randef(rmod2) #estimate random effects
#random slopes with fixed intercepts (suppressed correlation)
rmod3 <- mixedmirt(data, covdata, 1, fixed = ~ 0 + items, random = ~ -1 + pseudoIQ|group)
summary(rmod3)
eff <- randef(rmod3)
str(eff)
###################################################
##LLTM, and 2PL version of LLTM
data(SAT12)
data <- key2binary(SAT12,
key = c(1,4,5,2,3,1,2,1,3,1,2,4,2,1,5,3,4,4,1,4,3,3,4,1,3,5,1,3,1,5,4,5))
model <- 'Theta = 1-32'
# Suppose that the first 16 items were suspected to be easier than the last 16 items,
# and we wish to test this item structure hypothesis (more intercept designs are possible
# by including more columns).
itemdesign <- data.frame(itemorder = factor(c(rep('easier', 16), rep('harder', 16))))
#notice that the 'fixed = ~ ... + items' argument is omitted
LLTM <- mixedmirt(data, model = model, fixed = ~ 0 + itemorder, itemdesign = itemdesign,
SE = TRUE) # SE argument ensures that the information matrix is computed accurately
summary(LLTM)
coef(LLTM)
wald(LLTM)
L <- matrix(c(-1, 1, 0), 1)
wald(LLTM, L) #first half different from second
#compare to items with estimated slopes (2PL)
twoPL <- mixedmirt(data, model = model, fixed = ~ 0 + itemorder, itemtype = '2PL',
itemdesign = itemdesign)
#twoPL not mixing too well (AR should be between .2 and .5), decrease MHcand
twoPL <- mixedmirt(data, model = model, fixed = ~ 0 + itemorder, itemtype = '2PL',
itemdesign = itemdesign, technical = list(MHcand = 0.8))
anova(twoPL, LLTM) #much better fit
summary(twoPL)
coef(twoPL)
wald(twoPL)
L <- matrix(0, 1, 34)
L[1, 1] <- 1
L[1, 2] <- -1
wald(twoPL, L) #n.s., which is the correct conclusion. Rasch approach gave wrong inference
##LLTM with item error term
LLTMwithError <- mixedmirt(data, model = model, fixed = ~ 0 + itemorder, random = ~ 1|items,
itemdesign = itemdesign)
summary(LLTMwithError)
#large item level variance after itemorder is regressed; not a great predictor of item difficulty
coef(LLTMwithError)
###################################################
### Polytomous example
#make an arbitrary group difference
covdat <- data.frame(group = rep(c('m', 'f'), nrow(Science)/2))
#partial credit model
mod <- mixedmirt(Science, covdat, model=1, fixed = ~ 0 + group)
coef(mod)
#gpcm to estimate slopes
mod2 <- mixedmirt(Science, covdat, model=1, fixed = ~ 0 + group,
itemtype = 'gpcm')
summary(mod2)
anova(mod, mod2)
#graded model
mod3 <- mixedmirt(Science, covdat, model=1, fixed = ~ 0 + group,
itemtype = 'graded')
coef(mod3)
###################################################
# latent regression with Rasch and 2PL models
set.seed(1)
n <- 300
a <- matrix(1, 10)
d <- matrix(rnorm(10))
Theta <- matrix(c(rnorm(n, 0), rnorm(n, 1), rnorm(n, 2)))
covdata <- data.frame(group=rep(c('g1','g2','g3'), each=n))
dat <- simdata(a, d, N=n*3, Theta=Theta, itemtype = '2PL')
#had we known the latent abilities, we could have computed the regression coefs
summary(lm(Theta ~ covdata$group))
#but all we have is observed test data. Latent regression helps to recover these coefs
#Rasch model approach (and mirt equivalent)
rmod0 <- mirt(dat, 1, 'Rasch') # unconditional
# these two models are equivalent
rmod1a <- mirt(dat, 1, 'Rasch', covdata = covdata, formula = ~ group)
rmod1b <- mixedmirt(dat, covdata, 1, fixed = ~ 0 + items + group)
anova(rmod0, rmod1b)
coef(rmod1a, simplify=TRUE)
summary(rmod1b)
# 2PL, requires different input to allow Theta variance to remain fixed
mod0 <- mirt(dat, 1) # unconditional
mod1a <- mirt(dat, 1, covdata = covdata, formula = ~ group, itemtype = '2PL')
mod1b <- mixedmirt(dat, covdata, 1, fixed = ~ 0 + items, lr.fixed = ~group, itemtype = '2PL')
anova(mod0, mod1b)
coef(mod1a)$lr.betas
summary(mod1b)
# specifying specific regression effects is accomplished by passing a list of formula
model <- 'F1 = 1-5
F2 = 6-10'
covdata$contvar <- rnorm(nrow(covdata))
mod2 <- mirt(dat, model, itemtype = 'Rasch', covdata=covdata,
formula = list(F1 = ~ group + contvar, F2 = ~ group))
coef(mod2)[11:12]
mod2b <- mixedmirt(dat, covdata, model, fixed = ~ 0 + items,
lr.fixed = list(F1 = ~ group + contvar, F2 = ~ group))
summary(mod2b)
####################################################
## Simulated Multilevel Rasch Model
set.seed(1)
N <- 2000
a <- matrix(rep(1,10),10,1)
d <- matrix(rnorm(10))
cluster = 100
random_intercept = rnorm(cluster,0,1)
Theta = numeric()
for (i in 1:cluster)
Theta <- c(Theta, rnorm(N/cluster,0,1) + random_intercept[i])
group = factor(rep(paste0('G',1:cluster), each = N/cluster))
covdata <- data.frame(group)
dat <- simdata(a,d,N, itemtype = rep('2PL',10), Theta=matrix(Theta))
# null model
mod1 <- mixedmirt(dat, covdata, 1, fixed = ~ 0 + items, random = ~ 1|group)
summary(mod1)
# include level 2 predictor for 'group' variance
covdata$group_pred <- rep(random_intercept, each = N/cluster)
mod2 <- mixedmirt(dat, covdata, 1, fixed = ~ 0 + items + group_pred, random = ~ 1|group)
# including group means predicts nearly all variability in 'group'
summary(mod2)
anova(mod1, mod2)
# can also be fit for Rasch/non-Rasch models with the lr.random input
mod1b <- mixedmirt(dat, covdata, 1, fixed = ~ 0 + items, lr.random = ~ 1|group)
summary(mod1b)
mod2b <- mixedmirt(dat, covdata, 1, fixed = ~ 0 + items + group_pred, lr.random = ~ 1|group)
summary(mod2b)
anova(mod1b, mod2b)
mod3 <- mixedmirt(dat, covdata, 1, fixed = ~ 0 + items, lr.random = ~ 1|group, itemtype = '2PL')
summary(mod3)
anova(mod1b, mod3)
head(cbind(randef(mod3)$group, random_intercept))
## End(Not run)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.