if (FALSE) {
name <- 'old2PL'
par <- c(a = .5, b = -2)
est <- c(TRUE, TRUE)
P.old2PL <- function(par,Theta, ncat){
a <- par[1]
b <- par[2]
P1 <- 1 / (1 + exp(-1*a*(Theta - b)))
cbind(1-P1, P1)
}
x <- createItem(name, par=par, est=est, P=P.old2PL)
# So, let's estimate it!
dat <- expand.table(LSAT7)
sv <- mirt(dat, 1, c(rep('2PL',4), 'old2PL'), customItems=list(old2PL=x), pars = 'values')
tail(sv) #looks good
mod <- mirt(dat, 1, c(rep('2PL',4), 'old2PL'), customItems=list(old2PL=x))
coef(mod)
mod2 <- mirt(dat, 1, c(rep('2PL',4), 'old2PL'), customItems=list(old2PL=x), method = 'MHRM')
coef(mod2)
# same definition as above, but using symbolic derivative computations
# (can be more accurate/stable)
xs <- createItem(name, par=par, est=est, P=P.old2PL, derivType = 'symbolic')
mod <- mirt(dat, 1, c(rep('2PL',4), 'old2PL'), customItems=list(old2PL=xs))
coef(mod, simplify=TRUE)
# several secondary functions supported
M2(mod, calcNull=FALSE)
itemfit(mod)
fscores(mod, full.scores=FALSE)
plot(mod)
# fit the same model, but specify gradient function explicitly (use of a browser() may be helpful)
gr <- function(x, Theta){
# browser()
a <- x@par[1]
b <- x@par[2]
P <- probtrace(x, Theta)
PQ <- apply(P, 1, prod)
r_P <- x@dat / P
grad <- numeric(2)
grad[2] <- sum(-a * PQ * (r_P[,2] - r_P[,1]))
grad[1] <- sum((Theta - b) * PQ * (r_P[,2] - r_P[,1]))
## check with internal numerical form to be safe
# numerical_deriv(x@par[x@est], mirt:::EML, obj=x, Theta=Theta)
grad
}
x <- createItem(name, par=par, est=est, P=P.old2PL, gr=gr)
mod <- mirt(dat, 1, c(rep('2PL',4), 'old2PL'), customItems=list(old2PL=x))
coef(mod, simplify=TRUE)
### non-linear
name <- 'nonlin'
par <- c(a1 = .5, a2 = .1, d = 0)
est <- c(TRUE, TRUE, TRUE)
P.nonlin <- function(par,Theta, ncat=2){
a1 <- par[1]
a2 <- par[2]
d <- par[3]
P1 <- 1 / (1 + exp(-1*(a1*Theta + a2*Theta^2 + d)))
cbind(1-P1, P1)
}
x2 <- createItem(name, par=par, est=est, P=P.nonlin)
mod <- mirt(dat, 1, c(rep('2PL',4), 'nonlin'), customItems=list(nonlin=x2))
coef(mod)
### nominal response model (Bock 1972 version)
Tnom.dev <- function(ncat) {
T <- matrix(1/ncat, ncat, ncat - 1)
diag(T[-1, ]) <- diag(T[-1, ]) - 1
return(T)
}
name <- 'nom'
par <- c(alp=c(3,0,-3),gam=rep(.4,3))
est <- rep(TRUE, length(par))
P.nom <- function(par, Theta, ncat){
alp <- par[1:(ncat-1)]
gam <- par[ncat:length(par)]
a <- Tnom.dev(ncat) %*% alp
c <- Tnom.dev(ncat) %*% gam
z <- matrix(0, nrow(Theta), ncat)
for(i in 1:ncat)
z[,i] <- a[i] * Theta + c[i]
P <- exp(z) / rowSums(exp(z))
P
}
nom1 <- createItem(name, par=par, est=est, P=P.nom)
nommod <- mirt(Science, 1, 'nom1', customItems=list(nom1=nom1))
coef(nommod)
Tnom.dev(4) %*% coef(nommod)[[1]][1:3] #a
Tnom.dev(4) %*% coef(nommod)[[1]][4:6] #d
}
Run the code above in your browser using DataLab