[R] IRT Likelihood problem

From: Doran, Harold <HDoran_at_air.org>
Date: Wed, 12 Dec 2007 08:27:27 -0500

I have the following item response theory (IRT) likelihood that I want to maximize w.r.t. to theta (student ability).

L(\theta) = \prod(p(x))

Where p(x) is the 3-parameter logistic model when items are scored dichotomously (x_{ij} = 0 or 1) and p(x) is Muraki's generalized partial credit model when items are scored polytomously (x_{ij} = 0 \ldots J).

Now, I wrote the following two functions to maximize the likelihood. The first one is for the 3PL and works when all items are scored dichotmously. The second is for the GPCM and gives the MLE when all items are scored polytomously.

In the code below, this requires that I first have in hand estimates of the item parameters so that the maximization is only w.r.t. theta.

# x = response pattern
# b = vector of location parameters
# a = vector of discrimination parameters
# c = vector of guessing paramete

theta.3pl <- function(x, b, a, c){
 opt <- function(theta) -sum(dbinom(x, 1, c + ((1-c)/(1 + exp(-1.7*a*(theta - b)))), log = TRUE))  start_val <- log(sum(x)/(length(x)/sum(x)))  out <- optim(start_val , opt, method = "BFGS", hessian = TRUE)  out$par

# score = the category the student scored in for item i
# d = the item parameters for the ith item
# a = the discrimination p
# Muraki's GPCM

pcm.max <- function(score, d, a){

    pcm <- function(theta, d, score, a)

    opt <- function(theta) -sum(log(mapply(pcm, d, theta = theta, score= score )))

    start_val <- log(sum(score-1)/(length(score-1)/sum(score-1)))     out <- optim(start_val, opt, method = "BFGS", hessian = TRUE)     round(out$par, 2)

However, I have data for which there is a mixture of item types. Some are dichotomous and others are polytomous. Therefore, I need to somehow modify these functions to work in a conditional statement that first evaluates whether the item is dichotomous or not and then uses the right function to write out and then maximize the likelihood. However, I'm a bit stumped on how I might code this. Can anyone suggest how that might work?

For example, assume I have a test consisting of 3 items. The first two are dichotomous and the last is polytomous. Assume the students score on these three items is:

x <- c(1,0,3) # that is, 'right', 'wrong', 'scored in category 3'

And further assume the item parameters for these items are

Item 1 c = .11, b = 1.2, a = .58
Item 2 c = .20, b = .65, a = 1.2
Item 3 d = (0, -1.4, -.28, .95)

Now, my function pcm.max also reduces to Master's partial credit model when a = 1 for all items. And, because Master's partial credit reduces to the Rasch model when items are scored dichotmously, the only function I need is pcm.max when working within the Rasch family of models. However, Muraki's model does not reduce to the 3PL because of the guessing parameter when items are scored dichotomously, so I think the problem is slightly more complex in this scenario and requires the use of both functions.

Last, and a bit tangentially, it is possible when maximizing the 3PL that I will converge upon a local and not global maximum, a situation that does not occur with Rasch or 2PL. I only know of two methods for knowing whether I converged upon the right max. First, plot the likelihood and look at it or second, use different starting points and see if I converge to the same max. However, I am maximizing this likelihood over 80,000 students (or so) and I don't think either of these two methods are viable. If anyone has a suggestion on how I could proceed thoughtfully on that issue I welcome that as well.

Many thanks,

R-help_at_r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Received on Wed 12 Dec 2007 - 13:31:00 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Wed 12 Dec 2007 - 14:00:19 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.