rasch.mml2: Estimation of the Generalized Logistic Item Response Model, Ramsay's Quotient Model and a Nonparametric Item Response Model

Description

This function employs marginal maximum likelihood estimation of three item response models for dichotomous data. First, the Rasch type model (generalized item response model) can be estimated. The generalized logistic link function (Stukel, 1988) can be estimated or fixed for conducting IRT with different link functions than the logistic one. The Four-Parameter logistic item response model is a special case of this model (Loken & Rulison, 2010). Second, Ramsay's quotient model (Ramsay, 1989) can be estimated by specifying irtmodel="ramsay.qm". Third, quite general item response functions can be estimated in a nonparametric framework (Rossi, Wang & Ramsay, 2002).

Usage

rasch.mml2( dat , theta.k=seq(-6,6,len=21) , group=NULL , weights=NULL ,  
   constraints=NULL , glob.conv=10^(-5) , parm.conv=10^(-4) , mitermax=1 , 
   mmliter=1000 , progress=TRUE ,  fixed.a=rep(1,ncol(dat)) , 
   fixed.c=rep(0,ncol(dat)) , fixed.d=rep(1,ncol(dat)) , 
   fixed.K=rep(3,ncol(dat)) , b.init=NULL , est.a=NULL , est.b=NULL , 
   est.c=NULL , est.d=NULL , min.c=0 , max.c=1 , min.d=0 , max.d=1 ,
   est.K=NULL , min.K=1 , max.K=20 , pid=1:(nrow(dat)) , trait.weights=NULL ,
   center.trait=TRUE , alpha1=0 , alpha2=0 ,est.alpha=FALSE , equal.alpha=FALSE , 
   designmatrix=NULL , alpha.conv=parm.conv , numdiff.parm=0.001 ,
   numdiff.alpha.parm= numdiff.parm , distribution.trait="normal" , Qmatrix=NULL , 
   variance.fixed=NULL , mu.fixed=cbind(seq(1,ncol(Qmatrix)),rep(0,ncol(Qmatrix))) ,
   irtmodel="raschtype" , npformula=NULL , npirt.monotone=TRUE ,
   use.freqpatt = is.null(group)  , ... )
   
## S3 method for class 'rasch.mml':
summary(object,...)

## S3 method for class 'rasch.mml':
plot(x,items=NULL, xlim=NULL, main=NULL, ...)

Arguments

dat

An $N \times I$ data frame of dichotomous item responses

theta.k

Optional vector of discretized theta values. For multidimensional IRT models with $D$ dimensions, it is a matrix with $D$ columns.

group

Vector of integers with group identifiers in multiple group estimation.

weights

Optional vector of person weights (sample weights).

constraints

Constraints on b parameters (item difficulties). It must be a matrix with two columns: the first column contains item names, the second column fixed parameter values.

glob.conv

Convergence criterion for deviance

parm.conv

Convergence criterion for item parameters

mitermax

Maximum number of iterations in M step. This argument does only apply for the estimation of the $b$ parameters.

mmliter

Maximum number of iterations

progress

Should progress be displayed at the console?

fixed.a

Fixed or initial $a$ parameters

fixed.c

Fixed or initial $c$ parameters

fixed.d

Fixed or initial $d$ parameters

fixed.K

Fixed or initial $K$ parameters in Ramsay's quotient model.

b.init

Initial $b$ parameters

est.a

Vector of integers which indicate which $a$ parameters should be estimated. Equal integers correspond to the same estimated parameters.

est.b

Vector of integers which indicate which $b$ parameters should be estimated. Equal integers correspond to the same estimated parameters.

est.c

Vector of integers which indicate which $c$ parameters should be estimated. Equal integers correspond to the same estimated parameters.

est.d

Vector of integers which indicate which $d$ parameters should be estimated. Equal integers correspond to the same estimated parameters.

min.c

Minimal $c$ parameter to be estimated

max.c

Maximal $c$ parameter to be estimated

min.d

Minimal $d$ parameter to be estimated

max.d

Maximal $d$ parameter to be estimated

est.K

Vector of integers which indicate which $K$ parameters should be estimated. Equal integers correspond to the same estimated parameters.

min.K

Minimal $K$ parameter to be estimated

max.K

Maximal $K$ parameter to be estimated

pid

Optional vector of person identifiers

trait.weights

Optional vector of trait weights for a fixing the trait distribution.

center.trait

Should the trait distribution be centered

alpha1

Fixed or initial $\alpha_1$ parameter

alpha2

Fixed or initial $\alpha_2$ parameter

est.alpha

Should $\alpha$ parameters be estimated?

equal.alpha

Estimate $\alpha$ parameters under the assumption $\alpha_1 = \alpha_2$?

designmatrix

Design matrix for item difficulties $b$ to estimate linear logistic test models

alpha.conv

Convergence criterion for $\alpha$ parameter

numdiff.parm

Parameter for numerical differentiation

numdiff.alpha.parm

Parameter for numerical differentiation for $\alpha$ parameter

distribution.trait

Assumed trait distribution. The default is the normal distribution ("normal"). Log-linear smoothing of the trait distribution is also possible ("smooth2", "smooth3" or "smooth4" for smoothing up to

Qmatrix

The Q-matrix

variance.fixed

Matrix for fixing covariance matrix (See Examples)

mu.fixed

Matrix for fixing mean vector (See Examples)

irtmodel

Specify estimable IRT models: raschtype (Rasch type model), ramsay.qm (Ramsay's quotient model), npirt (Nonparametric item response model). If npirt is used as the argument for irtmodel,

npformula

A string or a vector which contains Rformula objects for specifying the item response function. For example, "y~theta" is the specification of the 2PL model (see Details). If irtmodel="npirt" and npformula is not

npirt.monotone

Should nonparametrically estimated item response functions be monotone? The default is TRUE. This function applies only to irtmodel='npirt' and npformula=NULL.

use.freqpatt

A logical if frequencies of pattern should be used or not. The default is is.null(group). This means that for single group analyses, frequency patterns are used but not for multiple groups. If data processing times are large, then

object

Object of class rasch.mml

items

Vector of integer or item names which should be plotted

xlim

Specification for xlim in plot

main

Title of the plot

...

Further arguments to be passed

Value

A list with following entries
datOriginal data frame
itemEstimated item parameters in the generalized item response model
item2Estimted item parameters for Ramsay's quotient model
trait.distrDiscretized ability distribution points and probabilities
mean.traitEstimated mean vector
sd.traitEstimated standard deviations
skewness.traitEstimated skewnesses
devianceDeviance
pjkEstimated probabilities of item correct evaluated at theta.k
personPerson parameter estimates: mode (MAP) and mean (EAP) of the posterior distribution
pidPerson identifier
ability.est.patternResponse pattern estimates
f.qk.yiIndividual posterior distribution
f.yi.qkIndividual likelihood
fixed.aEstimated $a$ parameters
fixed.cEstimated $c$ parameters
GNumber of groups
alpha1Estimated $\alpha_1$ parameter in generalized logistic item response model
alpha2Estimated $\alpha_2$ parameter in generalized logistic item response model
se.bStandard error of $b$ parameter in generalized logistic model or Ramsay's quotient model
se.aStandard error of $a$ parameter in generalized logistic model
se.cStandard error of $c$ parameter in generalized logistic model
se.dStandard error of $d$ parameter in generalized logistic model
se.alphaStandard error of $\alpha$ parameter in generalized logistic model
se.KStandard error of $K$ parameter in Ramsay's quotient model
iterNumber of iterations
reliabilityEAP reliability
irtmodelType of estimated item response model
DNumber of dimensions
muMean vector (for multdimensional models)
Sigma.covCovariance matrix (for multdimensional models)
theta.kGrid of discretized ability distributions
trait.weightsFixed vector of probabilities for the ability distribution
icInformation criteria
...

Details

The item response function of the generalized item response model (irtmodel="raschtype"; Stukel, 1988) can be written as $$P( X_{pi} = 1 | \theta_{pd} ) = c_i + (d_i - c_i ) g_{\alpha_1 , \alpha_2} [ a_i ( \theta_{pd} - b_i ) ]$$ where $g$ is the generalized logistic link function depending on parameters $\alpha_1$ and $\alpha_2$. The logistic link function corresponds to $\alpha_1=0$ and $\alpha_2=0$. A $D$-dimensional model can also be specified but only allows for between item dimensionality (one item loads on one and only dimension). Setting $c_i=0$, $d_i=1$ and $a_i=1$ for all items $i$, an additive item response model $$P( X_{pi} = 1 | \theta_p ) = g_{\alpha_1 , \alpha_2} ( \theta_p - b_i )$$ is estimated. Ramsay's quotient model (irtmodel="qm.ramsay") uses the item response function $$P( X_{pi} = 1 | \theta_p ) = \frac{ \exp(\theta_p / b_i)} { K_i + \exp (\theta_p / b_i )}$$ Quite general unidimensional item response models can be estimated in a nonparametric framework (irtmodel="npirt"). The response functions are a linear combination of transformed $\theta$ values $$logit[ P( X_{pi} = 1 | \theta_p ) ] = Y_\theta \beta$$ Where $Y_\theta$ is a design matrix of $\theta$ and $\beta$ are item parameters to be estimated. The formula $Y_\theta \beta$ can be specified in the Rformula framework.

References

Loken, E., & Rulison, K. L. (2010). Estimation of a four-parameter item response theory model. British Journal of Mathematical and Statistical Psychology, 63, 509-525. Ramsay, J. O. (1989). A comparison of three simple test theory models. Psychometrika, 54, 487-499. Rossi, N., Wang, X., & Ramsay, J. O. (2002). Nonparametric item response function estimates with the EM algorithm. Journal of Educational and Behavioral Statistics, 27, 291-317. Stukel, T. A. (1988). Generalized logistic models. Journal of the American Statistical Association, 83, 426-431. van der Maas, H. J. L., Molenaar, D., Maris, G., Kievit, R. A., & Borsboom, D. (2011). Cognitive psychology meets psychometric theory: On the relation between process models for decision making and latent variable models for individual differences. Psychological Review, 318, 339-356.

Examples

Run this code

#############################################################################
# EXAMPLE 1: Reading dataset
#############################################################################

data(data.read)
dat <- data.read
I <- ncol(dat) # number of items

# Rasch model
mod1 <- rasch.mml2( dat )
summary(mod1)
plot( mod1 )    # plot all items
# title 'Rasch model', display curves from -3 to 3 only for items 1, 5 and 8
plot(mod1, main="Rasch model Items 1, 5 and 8", xlim=c(-3,3) , items=c(1,5,8) )

# Rasch model with constraints on item difficulties
# set item parameters of A1 and C3 equal to -2
constraints <- data.frame( c("A1","C3") , c(-2,-2) )
mod1a <- rasch.mml2( dat , constraints=constraints)
summary(mod1a)

# estimate equal item parameters for 1st and 11th item
est.b <- 1:I
est.b[11] <- 1
mod1b <- rasch.mml2( dat , est.b = est.b )
summary(mod1b)

# estimate Rasch model with skew trait distribution
mod1c <- rasch.mml2( dat , distribution.trait="smooth3")
summary(mod1c)

# 2PL model
mod2 <- rasch.mml2( dat , est.a = 1:I )
summary(mod2)
plot(mod2)	# plot 2PL item response curves

# 3PL model
mod3 <- rasch.mml2( dat , est.a = 1:I , est.c = 1:I , 
        mmliter = 400 # maximal 400 iterations
                 )
summary(mod3)

# 3PL model with fixed guessing paramters of 1/4
# equal slopes
mod4 <- rasch.mml2( dat , fixed.c = rep(1/4 , I )   )
summary(mod4)

# 3PL model with equal guessing paramters for all items
mod5 <- rasch.mml2( dat , est.c = rep(1, I )   )
summary(mod5)

# difficulty + guessing model
mod6 <- rasch.mml2( dat , est.c = 1:I   )
summary(mod6)

# 4PL model
mod7 <- rasch.mml2( dat , est.a = 1:I , est.c=1:I , est.d = 1:I  ,
        min.d = .95 , max.c = .25)  
        # set minimal d and maximal c parameter to .95 and .25
summary(mod7)

# constrained 4PL model
# equal slope, guessing and slipping parameters
mod8 <- rasch.mml2( dat ,est.c=rep(1,I) , est.d = rep(1,I) )
summary(mod8)

# estimation of an item response model with an
# uniform theta distribution
theta.k <- seq( 0.01 , .99 , len=20 )
trait.weights <- rep( 1/length(theta.k) , length(theta.k) )
mod9 <- rasch.mml2( dat , theta.k=theta.k , trait.weights = trait.weights ,
            normal.trait=FALSE , est.a = 1:12  )
summary(mod9)

#############################################################################
# EXAMPLE 2: Longitudinal data
#############################################################################

data(data.long)
dat <- data.long[,-1]

# define Q loading matrix
Qmatrix <- matrix( 0 , 12 , 2 )
Qmatrix[1:6,1] <- 1 # T1 items
Qmatrix[7:12,2] <- 1    # T2 items

# define restrictions on item difficulties
est.b <- c(1,2,3,4,5,6,   3,4,5,6,7,8) 
mu.fixed <- cbind(1,0) 
    # set first mean to 0 for identification reasons

# Model 1: 2-dimensional Rasch model
mod1 <- rasch.mml2( dat , Qmatrix=Qmatrix , miterstep=4,
    est.b = est.b   ,  mu.fixed = mu.fixed , mmliter=30 )
summary(mod1)
plot(mod1)
##     Plot function is only applicable for unidimensional models

#############################################################################
# SIMULATED EXAMPLE 3
#  one group, estimation of alpha parameter 
#		in generalized logistic model
#############################################################################

# simulate theta values
set.seed(786)
N <- 1000                  # number of persons
theta <- rnorm( N , sd =1.5 ) # N persons with SD 1.5
b <- seq( -2 , 2 , len=15)

# simulate data
dat <- sim.raschtype( theta = theta , b = b , alpha1 = 0 , alpha2 = -0.3 )

#  estimating alpha parameters
mod1 <- rasch.mml2( dat , est.alpha = TRUE , mmliter=30 )
summary(mod1)
plot(mod1)

# fixed alpha parameters
mod1b <- rasch.mml2( dat , est.alpha = FALSE , alpha1=0 , alpha2=-.3 )
summary(mod1b)

# estimation with equal alpha parameters
mod1c <- rasch.mml2( dat , est.alpha = TRUE , equal.alpha=TRUE )
summary(mod1c)

# Ramsay QM
mod2a <- rasch.mml2( dat , irtmodel ="ramsay.qm" )
summary(mod2a)

# Ramsay QM with estimated K parameters
mod2b <- rasch.mml2( dat , irtmodel ="ramsay.qm" , est.K=1:15 , mmliter=30)
summary(mod2b)
plot(mod2b)

# nonparametric estimation of monotone item response curves
mod3a <- rasch.mml2( dat , irtmodel ="npirt" , mmliter =100 ,
            theta.k = seq( -3 , 3 , len=10) ) # evaluations at 10 theta grid points
# nonparametric ICC of first 4 items
round( t(mod3a$pjk)[1:4,] , 3 ) 
summary(mod3a)
plot(mod3a)

# nonparametric IRT estimation without monotonicity assumption
mod3b <- rasch.mml2( dat , irtmodel ="npirt" , mmliter =10 ,
            theta.k = seq( -3 , 3 , len=10) , npirt.monotone=FALSE) 
plot(mod3b)
            
# B-Spline estimation of ICCs
library(splines)
mod3c <- rasch.mml2( dat , irtmodel ="npirt" , 
        npformula = "y~bs(theta,df=3)" , theta.k = seq(-3,3,len=15) )
summary(mod3c)
round( t(mod3c$pjk)[1:6,] , 3 )
plot(mod3c)

# nonparametric estimation of the link function
# In this example, a bit complicated ICC is estimated:
# "theta + I( theta^2)" is a quadratic item response function
mod3d <- rasch.mml2( dat , irtmodel ="npirt" , 
        npformula = "y~theta + I(theta^2)" )
summary(mod3d)
plot(mod3d)

# estimation of a stepwise ICC function
# ICCs are constant on the theta domains: [-Inf,-1], [-1,1], [1,Inf]
mod3e <- rasch.mml2( dat , irtmodel ="npirt" , 
        npformula = "y~I(theta>-1 )+I(theta>1)" )
summary(mod3e)
plot(mod3e , xlim=c(-2.5,2.5) )

# 2PL model
mod4 <- rasch.mml2( dat ,  est.a=1:15)
summary(mod4)

#############################################################################
# SIMULATED EXAMPLE 4
#  two groups, estimation of generalized logistic model
#############################################################################

# simulate generalized logistic Rasch model in two groups
set.seed(8765)
N1 <- 1000     # N1=1000 persons in group 1
N2 <- 500      # N2= 500 persons in group 2
dat1 <- sim.raschtype( theta = rnorm( N1 , sd = 1.5 )  , b = b , 
			alpha1 = -0.3 , alpha2=0)
dat2 <- sim.raschtype( theta = rnorm( N2 , mean=-.5 , sd =.75)  , 
			b = b , alpha1 = -0.3 , alpha2=0)
dat1 <- rbind( dat1 , dat2 )
group <- c( rep(1,N1) , rep(2,N2))

mod1 <-  rasch.mml2( dat1 , parm.conv=.0001 , group=group , est.alpha = TRUE )
summary(mod1)

#############################################################################
# SIMULATED EXAMPLE 5: Multidimensional model
#############################################################################

#***
# (1) simulate data
set.seed(785)
library(mvtnorm)
N <- 500
theta <- rmvnorm( N,mean=c(0,0), sigma=matrix( c(1.45,.5,.5,1.7) , 2 , 2 ))
I <- 10
# 10 items load on the first dimension
p1 <- plogis( outer( theta[,1] , seq( -2 , 2 , len=I ) , "-" ) )
resp1 <- 1 * ( p1 > matrix( runif( N*I ) , nrow=N , ncol=I ) )
# 10 items load on the second dimension
p1 <- plogis( outer( theta[,2] , seq( -2 , 2 , len=I ) , "-" ) )
resp2 <- 1 * ( p1 > matrix( runif( N*I ) , nrow=N , ncol=I ) )
#Combine the two sets of items into one response matrix
resp <- cbind(resp1,resp2)
colnames(resp) <- paste("I" , 1:(2*I), sep="")
dat <- resp

# define Q-matrix
Qmatrix <- matrix( 0 , 2*I , 2 )
Qmatrix[1:I,1] <- 1
Qmatrix[1:I+I,2] <- 1

#***
# (2) estimation of models
# 2-dimensional Rasch model
mod1 <- rasch.mml2( dat , Qmatrix=Qmatrix )
summary(mod1)

# 2-dimensional 2PL model
mod2 <- rasch.mml2( dat , Qmatrix=Qmatrix , est.a = 1:(2*I) )               
summary(mod2)

# estimation with some fixed variances and covariances
# set variance of 1st dimension to 1 and 
#  covariance to zero
variance.fixed <- matrix( cbind(c(1,1) , c(1,2) , c(1,0)) , 
             byrow=FALSE , ncol= 3 )
mod3 <- rasch.mml2( dat , Qmatrix=Qmatrix , variance.fixed = variance.fixed )
summary(mod3)

# constraints on item difficulties
#  useful for example in longitudinal linking
est.b <- c( 1:I , 1:I )     
    # equal indices correspond to equally estimated item parameters
mu.fixed <- cbind( 1 , 0 )
mod4 <- rasch.mml2( dat, Qmatrix=Qmatrix, est.b = est.b , mu.fixed = mu.fixed )
summary(mod4)

Run the code above in your browser using DataLab