# R/itemfit.R In irtoys: A Collection of Functions Related to Item Response Theory (IRT)

#### Documented in apiitf

grp = function(theta, bins=9, breaks=NULL, equal="count", type="means") {
if (!is.null(dim(theta))) theta = theta[,1]
if (is.null(breaks)) {
if (equal=="count") {
qu = seq(0, 1, 1/bins)
breaks=quantile(theta, probs=qu, names=FALSE, na.rm=TRUE)
} else {
r = range(theta)
st = (r-r)/bins
breaks = seq(r, r, by=st)
}
}
grmemb = findInterval(theta, breaks, rightmost.closed=TRUE)
counts = table(grmemb)
ref = switch(type,
mids = {
dd = diff(breaks,1)
md = mean(dd)/2
c(breaks-md, breaks[-length(breaks)]+dd/2, breaks[length(breaks)]+md)
},
meds = tapply(theta, grmemb, median),
means= tapply(theta, grmemb, mean)
)
return(list(breaks=breaks,counts=counts,ref=ref,grmemb=grmemb))
}
# compute an item fit statistic with df and pval.
# Optionally plot the IRF with residuals shown

#' Test item fit
#'
#' Returns a statistic of item fit together with its degrees of freedom and
#' p-value. Optionally produces a plot.
#'
#' Given a long test, say 20 items or more, a large-test statistic of item fit
#' could be constructed by dividing examinees into groups of similar ability,
#' and comparing the observed proportion of correct answers in each group with
#' the expected proportion under the proposed model. Different statistics have
#' been proposed for this purpose.
#'
#' The chi-squared statistic
#' \deqn{X^2=\sum_g(N_g\frac{(p_g-\pi_g)^2}{\pi_g(1-\pi_g)},} where \eqn{N_g}
#' is the number of examinees in group \eqn{g}, \eqn{p_g=r_g/N_g}, \eqn{r_g} is
#' the number of correct responses to the item in group \eqn{g}, and
#' \eqn{\pi_g} is the IRF of the proposed model for the median ability in group
#' \eqn{g}, is attributed by Embretson & Reise to R. D. Bock, although the
#' article they cite does not actually mention it. The statistic is the sum of
#' the squares of quantities that are often called "Pearson residuals" in the
#' literature on categorical data analysis.
#'
#' BILOG uses the likelihood-ratio statistic
#' \deqn{X^2=2\sum_g\left[r_g\log\frac{p_g}{\pi_g} +
#' (N_g-r_g)\log\frac{(1-p_g)}{(1-\pi_g)}\right],} where \eqn{\pi_g} is now the
#' IRF for the mean ability in group \eqn{g}, and all other symbols are as
#' above.
#'
#' Both statistics are assumed to follow the chi-squared distribution with
#' degrees of freedom equal to the number of groups minus the number of
#' parameters of the model (eg 2 in the case of the 2PL model). The first
#' statistic is obtained in \code{itf} with \code{stat="chi"}, and the second
#' with \code{stat="lr"} (or not specifying \code{stat} at all).
#'
#' In the real world we can only work with estimates of ability, not with
#' ability itself. \code{irtoys} allows use of any suitable ability measure
#' via the argument \code{theta}. If \code{theta} is not specified, \code{itf}
#' will compute EAP estimates of ability, group them in 9 groups having
#' approximately the same number of cases, and use the means of the ability
#' eatimates in each group. This is the approximate behaviour of BILOG.
#'
#' If the test has less than 20 items, \code{itf} will issue a warning.
#' For tests of 10 items or less, BILOG has a special statistic of fit, which
#' can be found in the BILOG output. Also of interest is the fit in 2- and
#' 3-way marginal tables in package \code{ltm}.
#'
#' @param resp A matrix of responses: persons as rows, items as columns,
#' entries are either 0 or 1, no missing data
#' @param ip Item parameters: the object returned by \eqn{est},
#' or the equivalent of its first part.
#' @param item A single number pointing to the item (column of \code{resp}, row
#' of \code{ip}), for which fit is to be tested
#' @param stat The statistic to be computed, either \code{"chi"} or
#' \code{"lr"}. Default is \code{"lr"}. See details below.
#' @param theta A vector containing some viable estimate of the latent variable
#' for the same persons whose responses are given in \code{resp}. If not given
#' (and \code{group} is also missing), WLE estimates will be computed from
#' \code{resp} and \code{ip}.
#' @param standardize Standardize the distribution of ability estimates?
#' @param mu Mean of the standardized distribution of ability estimates
#' @param sigma Standard deviation of the standardized distribution of ability
#' estimates
#' @param bins Desired number of bins (default is 9)
#' @param breaks A vector of cutpoints. Overrides \code{bins} if present.
#' @param equal Either \code{"width"} for bins of equal width, or
#' \code{"count"} for bins with roughly counts of observations. Default is
#' \code{"quant"}
#' @param type The points at which \code{itf} will evaluate the IRF. One of
#' \code{"mids"} (the mid-point of each bin), \code{"meds"} (the median of the
#' values in the bin), or \code{"means"} (the mean of the values in the bin).
#' Default is \code{"means"}.
#' @param do.plot Whether to do a plot
#' @param main The title of the plot if one is desired
#' @return A vector of three numbers: \item{Statistic}{The value of the statistic of item fit}
#' \item{DF}{The degrees of freedom} \item{P-value}{The p-value}
#' @author Ivailo Partchev
#' @references S. E. Embretson and S. P. Reise (2000), Item Response Theory for
#' Psychologists, Lawrence Erlbaum Associates, Mahwah, NJ
#'
#' M. F. Zimowski, E. Muraki, R. J. Mislevy and R. D. Bock (1996), BILOG--MG.
#' Multiple-Group IRT Analysis and Test Maintenance for Binary Items, SSI
#' Scientific Software International, Chicago, IL
#' @keywords models
#'
#' @examples
#'
#' fit   <- itf(resp=Scored, ip=Scored2pl, item=7)
#'
itf = function(resp, ip, item, stat = "lr", theta,
standardize=TRUE, mu=0, sigma=1,
bins = 9, breaks = NULL, equal = "count", type = "means",
do.plot=TRUE, main="Item fit") {
if (is.list(ip)) ip = ip$est if (!(item %in% 1:nrow(ip))) stop("bad item number") if (nrow(ip)<20) warning("item fit statistic computed for a test of less than 20 items") if (missing(theta)) theta = wle(resp, ip) if (!is.null(dim(theta))) theta = theta[, 1] pa = ip[item,] iresp = resp[,item] oo = !is.na(iresp) iresp = iresp[oo] theta = theta[oo] if (standardize) theta = scale(theta)*sigma + mu groups = grp(theta, bins=bins, breaks=breaks, equal=equal, type=type) ep = as.vector(irf(ip=pa, x=groups$ref)$f) nn = groups$counts
gr = groups$grmemb * iresp rg = tabulate(gr[gr>0], nbins=length(nn)) op = rg / nn npar = 3 - (all(ip[,3]==0)) - (var(ip[,1])==0) if (stat=="chi") tst = sum((op-ep)^2/(ep*(1-ep))*nn) else { cmp = rg*log(op/ep)+(nn-rg)*log((1-op)/(1-ep)) tst = 2*sum(cmp[is.finite(cmp)]) } dfr = sum(nn>0) - npar pval=pchisq(tst,dfr,lower.tail=FALSE) if (do.plot) { ifit = paste("Q = ",round(tst,2)," d.f. = ",dfr," P = ",round(pval,4),sep="") plot(c(-4,4),c(0,1),main="",xlab="Ability",ylab="Proportion right",type="n") title(main=main,mtext(ifit,line=0.5,cex=0.7)) points(groups$ref, op)
}
result = c(tst,dfr,pval)
names(result) = c("Statistic","DF","P-value")
result
}

#' The api appropriateness index
#'
#' Computes the api appropriateness index, a measure of person fit in IRT models
#'
#'
#' @param resp A matrix of responses: persons as rows, items as columns,
#' entries are either 0 or 1, no missing data
#' @param ip Item parameters: the object returned by \eqn{est}, or its first part.
#' @param theta A measure of ability, typically produced with \code{mlebme},
#' \code{wle} etc. If missing, ML estimates will be computed automatically.
#' @return A vector of length equal to the number of rows in \code{resp},
#' containing the appropriateness indices
#' @author Ivailo Partchev
#' @references Drasgow, F., Levine, M. V., & Williams, E. A. (1985).
#' Appropriateness measurement with polychotomous item response models and
#' standardized indices. British Journal of Mathematical and Statistical
#' Psychology, 38, 67--80
#' @keywords models
#'
#' @examples
#'
#' api(Scored, Scored2pl)
#'
api = function(resp, ip, theta){
if (is.list(ip)) ip = ip$est if(missing(theta)) theta = wle(resp, ip) if (!is.null(dim(theta))) theta = theta[, 1] p = irf(ip=ip, x=theta)$f
p = pmax(p, .00001); pr = pmin(p, .99999)
if(is.null(dim(p))) p=matrix(p,ncol=length(p))
q = 1 - p
lp= log(p)
lq= log(q)
l = apply(lp*resp + lq*(1-resp), 1, sum)
m = apply(p*lp+q*lq, 1, sum)
v = apply(p*q*(log(p/q))^2, 1, sum)
(l-m)/sqrt(v)
}


## Try the irtoys package in your browser

Any scripts or data that you put into this service are public.

irtoys documentation built on May 2, 2019, 3:03 a.m.