Skip to content

Commit

Permalink
Merge branch 'master' of github.com:philchalmers/mirt
Browse files Browse the repository at this point in the history
  • Loading branch information
@philchalmers
philchalmers committed Feb 26, 2018
2 parents 44279df + 0d93163 commit 5abbe83
Show file tree
Hide file tree
Showing 6 changed files with 42 additions and 60 deletions.
2 changes: 0 additions & 2 deletions R/01-itemclass.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -57,8 +57,6 @@ setGeneric("logLik", function(object, ...) standardGeneric("logLik"))

setGeneric('DrawValues', function(x, Theta, ...) standardGeneric("DrawValues"))

setGeneric('RandomDeriv', function(x, ...) standardGeneric("RandomDeriv"))

setGeneric('GenRandomPars', function(x) standardGeneric("GenRandomPars"))

setGeneric('CheckIntercepts', function(x) standardGeneric("CheckIntercepts"))
22 changes: 9 additions & 13 deletions R/02-item_methods.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -335,19 +335,15 @@ setMethod(
}
)

setMethod(
f = "RandomDeriv",
signature = signature(x = 'RandomPars'),
definition = function(x, estHess = TRUE){
Theta <- x@drawvals
pick <- -seq_len(ncol(Theta))
out <- .Call("dgroup", x, Theta, matrix(0L), estHess, TRUE, FALSE, FALSE)
grad <- out$grad[pick]
hess <- out$hess[pick, pick, drop=FALSE]
diag(hess) <- -abs(diag(hess)) #hack for very small clusters
list(grad=grad, hess=hess)
}
)
RandomDeriv <- function(x, estHess = TRUE){
Theta <- x@drawvals
pick <- -seq_len(ncol(Theta))
out <- .Call("dgroup", x, Theta, matrix(0L), estHess, TRUE, FALSE, FALSE)
grad <- out$grad[pick]
hess <- out$hess[pick, pick, drop=FALSE]
diag(hess) <- -abs(diag(hess)) #hack for very small clusters
list(grad=grad, hess=hess)
}

# ----------------------------------------------------------------

Expand Down
8 changes: 2 additions & 6 deletions R/MHRM.utils.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -33,8 +33,6 @@ MHRM.deriv <- function(pars, gtheta, OffTerm, longpars, USE.FIXED, list, ngroups
for(i in seq_len(length(random))){
g[random[[i]]@parnum] <- 0
h[random[[i]]@parnum, random[[i]]@parnum] <- -diag(length(random[[i]]@parnum))
if(cycles == (RANDSTART - 1L)) #hack for R 3.4.0
deriv <- RandomDeriv(x=random[[i]], estHess=estHess)
}
} else {
for(i in seq_len(length(random))){
Expand Down Expand Up @@ -79,10 +77,8 @@ MHRM.deriv <- function(pars, gtheta, OffTerm, longpars, USE.FIXED, list, ngroups
}
grad <- grad[estpars & !redun_constr]
ave.h <- ave.h[estpars & !redun_constr, estpars & !redun_constr]
if(any(is.na(grad) | is.infinite(grad))){
stop('Model did not converge (unacceptable gradient caused by extreme parameter values)',
call.=FALSE)
}
if(any(is.infinite(grad) | is.nan(grad)))
stop('Inf or NaN values appeared in the gradient', call.=FALSE)
list(grad=grad, ave.h=ave.h)
}

Expand Down
1 change: 1 addition & 0 deletions R/PrepData.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@ PrepData <- function(data, model, itemtype, guess, upper, gpcm_mats, opts,
if(is.null(grsm.block)) grsm.block <- rep(1, ncol(data))
if(is.null(rsm.block)) rsm.block <- rep(1, ncol(data))
itemnames <- colnames(data)
if(any(itemnames == "")) stop("Items in data input must have valid names", call.=FALSE)
keywords <- c('COV', 'CONSTRAIN', 'CONSTRAINB', 'PRIOR', 'MEAN', 'START', 'LBOUND', 'UBOUND',
'FIXED', 'FREE', 'NEXPLORE')
data <- as.matrix(data)
Expand Down
3 changes: 2 additions & 1 deletion R/createItem.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -199,9 +199,10 @@ createItem <- function(name, par, est, P, gr=NULL, hss = NULL, gen = NULL,
dp1 <- array(x@dps(x@par, Theta, x@ncat), c(ThetaLength,x@ncat,xLength))
dp2 <- array(x@dps2(x@par, Theta, x@ncat), c(ThetaLength,x@ncat,xLength,xLength))
H <- matrix(NA,xLength,xLength)
P2 <- P^2
for (i in 1L:xLength){
for (j in i:xLength){
H[i,j] <- sum(x@dat*dp2[,,i,j]/P + x@dat*dp1[,,i]*(-dp1[,,j]/(P^2)))
H[i,j] <- sum(x@dat*dp2[,,i,j]/P + x@dat*dp1[,,i]*(-dp1[,,j]/P2))
H[j,i] <- H[i,j]
}
}
Expand Down
66 changes: 28 additions & 38 deletions tests/tests/test-07-mixedmirt.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -51,22 +51,19 @@ test_that('mixed dich', {
expect_equal(cfs, c(0.8023,1.5944,1,-1.3592,0,1,1.3988,2.9293,1,-2.4552,0,1,2.4305,4.4838,1,-3.1509,0,1,0.307,0.7794,1,-0.2878,0,1,1.4845,3.0533,1,-0.5464,0,1,1.2012,3.0349,1,-1.6528,0,1,0.9915,1.9697,1,-1.513,0,1,1.4304,2.6556,1,-2.4036,0,1,0.8535,2.0584,1,-1.8064,0,1,0.8374,2.3511,1,1.6148,0,1,0,0.1489),
tolerance = 1e-2)

mod_items <- try(mixedmirt(data, covdata, model, fixed = ~ 1, SE=FALSE, random = ~ 1|items,
verbose = FALSE, draws = 1), TRUE)
if(!is(mod_items, 'try-error')){
cfs <- c(coef(mod_items)[['GroupPars']], coef(mod_items)[['items']])
expect_equal(cfs[1:3], c(0.000, 1.083, 1.125), tolerance = 1e-2)
}
mod_items <- mixedmirt(data, covdata, model, fixed = ~ 1, SE=FALSE, random = ~ 1|items,
verbose = FALSE, draws = 1)
cfs <- c(coef(mod_items)[['GroupPars']], coef(mod_items)[['items']])
expect_equal(cfs[1:3], c(0.000, 1.083, 1.125), tolerance = 1e-2)

mod_items.group <- try(mixedmirt(data, covdata, model, fixed = ~ 1, SE=FALSE, random = ~ 1|items:group,
verbose = FALSE, draws = 1), TRUE)
if(!is(mod_items.group, 'try-error')){
cfs <- c(coef(mod_items.group)[['GroupPars']], coef(mod_items.group)[['items:group']])
expect_equal(cfs[1:3], c(0.000, 0.1431, 2.2536), tolerance = 1e-3)
set.seed(1)
bs <- boot.mirt(mod_items.group, R=2)
expect_is(bs, 'boot')
}
cfs <- c(coef(mod_items.group)[['GroupPars']], coef(mod_items.group)[['items:group']])
expect_equal(cfs[1:3], c(0.000, 0.1431, 2.2536), tolerance = 1e-3)
set.seed(1)
bs <- boot.mirt(mod_items.group, R=2)
expect_is(bs, 'boot')


#model using 2PL items instead of only Rasch, and with missing data
data[1,1] <- covdata[1,2] <- NA
Expand All @@ -81,13 +78,11 @@ test_that('mixed dich', {
covdata$group <- factor(rep(paste0('G',1:50), each = N/50))
rmod1 <- try(mixedmirt(data, covdata, 1, fixed = ~ 0 + items, random = ~ 1|group,
draws = 1, verbose = FALSE), TRUE)
if(!is(rmod1, 'try-error')){
expect_is(rmod1, 'MixedClass')
expect_equal(extract.mirt(rmod1, 'df'), 1011)
cfs <- as.numeric(do.call(c, coef(rmod1)))
expect_equal(cfs[124:129], c(0.06756121, 0.05018307, 0.08493935, 1.13460966, 0.66202405, 1.60719526),
tolerance = 1e-2)
}
expect_is(rmod1, 'MixedClass')
expect_equal(extract.mirt(rmod1, 'df'), 1011)
cfs <- as.numeric(do.call(c, coef(rmod1)))
expect_equal(cfs[124:129], c(0.06756121, 0.05018307, 0.08493935, 1.13460966, 0.66202405, 1.60719526),
tolerance = 1e-2)

#polytomous
covdat <- data.frame(group = rep(c('m', 'f'), nrow(Science)/2))
Expand Down Expand Up @@ -118,16 +113,15 @@ test_that('mixed dich', {
covdat$group <- factor(rep(paste0('G',1:20), length.out = nrow(Science)))
rmod1 <- try(mixedmirt(Science, covdat, model=model, draws=10, random = ~ 1|group,
itemtype = 'graded', verbose = FALSE, SE=FALSE), TRUE)
if(!is(rmod1, 'try-error')){
expect_is(rmod1, 'MixedClass')
expect_equal(extract.mirt(rmod1, 'df'), 238)
cfs <- as.numeric(na.omit(do.call(c, coef(rmod1))))
expect_equal(cfs, c(1.062578,4.895614,2.661662,-1.479457,1.195101,2.907836,0.897937,-2.253751,2.178545,5.083139,2.151222,-1.918151,1.08017,3.345869,0.9912499,-1.68683,0,1,0.006179096),
tolerance = 1e-4)
re <- randef(rmod1, ndraws=100)
expect_is(re, 'list')
expect_equal(length(re), 2)
}
expect_is(rmod1, 'MixedClass')
expect_equal(extract.mirt(rmod1, 'df'), 238)
cfs <- as.numeric(na.omit(do.call(c, coef(rmod1))))
expect_equal(cfs, c(1.062578,4.895614,2.661662,-1.479457,1.195101,2.907836,0.897937,-2.253751,2.178545,5.083139,2.151222,-1.918151,1.08017,3.345869,0.9912499,-1.68683,0,1,0.006179096),
tolerance = 1e-4)
re <- randef(rmod1, ndraws=100)
expect_is(re, 'list')
expect_equal(length(re), 2)


## latent regression
set.seed(1234)
Expand Down Expand Up @@ -161,16 +155,12 @@ test_that('mixed dich', {
covdata$cut <- factor(cut(Theta, breaks=2))
mod <- try(mixedmirt(dat, covdata = covdata, 1, fixed = ~ 0 + items, SE=FALSE,
random = ~ -1 + cut|group, verbose=FALSE, draws=1), TRUE)
if(!is(mod, 'try-error')){
so <- summary(mod, verbose=FALSE)
expect_equal(as.numeric(diag(so$random$group)), c(0.6155525, 0.4089267), tolerance = 1e-4)
}
so <- summary(mod, verbose=FALSE)
expect_equal(as.numeric(diag(so$random$group)), c(0.6155525, 0.4089267), tolerance = 1e-4)

mod <- try(mixedmirt(dat, covdata = covdata, 1, fixed = ~ 0 + items, SE=FALSE,
random = ~ -1 + theta|group, verbose=FALSE, draws=1), TRUE)
if(!is(mod, 'try-error')){
so <- summary(mod, verbose=FALSE)
expect_equal(as.numeric(diag(so$random$group)), c(0.2255360, 0.4568764), tolerance = 1e-4)
}
so <- summary(mod, verbose=FALSE)
expect_equal(as.numeric(diag(so$random$group)), c(0.2255360, 0.4568764), tolerance = 1e-4)

})

0 comments on commit 5abbe83

Please sign in to comment.