#############################################################################
# SIMULATED EXAMPLE 1: Simulated data according to Janssen et al. (2000, Table 2)
#############################################################################
N <- 1000
Ik <- c(4,6,8,5,9,6,8,6,5)
xi.k <- c( -.89 , -1.13 , -1.23 , .06 , -1.41 , -.66 , -1.09 , .57 , -2.44)
omega.k <- c(.98 , .91 , .76 , .74 , .71 , .80 , .79 , .82 , .54)
# select 4 attributes
K <- 4
Ik <- Ik[1:K] ; xi.k <- xi.k[1:K] ; omega.k <- omega.k[1:K]
sig2 <- 3.02
nu2 <- .09
I <- sum(Ik)
b <- rep( xi.k , Ik ) + rnorm(I , sd = sqrt(sig2) )
a <- rep( omega.k , Ik ) + rnorm(I , sd = sqrt(nu2) )
theta1 <- rnorm(N)
t1 <- rep(1,N)
p1 <- pnorm( outer(t1,a) * ( theta1 - outer(t1,b) ) )
dat <- 1 * ( p1 > runif(N*I) )
itemgroups <- rep( paste0("A" , 1:K ) , Ik )
# estimate model
mod <- mcmc.2pnoh(dat , itemgroups , burnin=200 , iter=1000 )
# summary
summary(mod)
# plot
plot(mod$mcmcobj , ask=TRUE)
# write coda files
mcmclist2coda( mod$mcmcobj , name = "simul_2pnoh" )
Run the code above in your browser using DataLab