## Not run: ------------------------------------
# ####################################
# # Example 1. #
# # GDINA, DINA, DINO #
# # ACDM, LLM and RRUM #
# # estimation and comparison #
# # #
# ####################################
#
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
#
# #--------GDINA model --------#
#
# mod1 <- GDINA(dat = dat, Q = Q, model = "GDINA")
# mod1
# # summary information
# summary(mod1)
#
# AIC(mod1) #AIC
# BIC(mod1) #BIC
# logLik(mod1) #log-likelihood value
# deviance(mod1) # deviance: -2 log-likelihood
# npar(mod1) # number of parameters
#
#
# head(indlogLik(mod1)) # individual log-likelihood
# head(indlogPost(mod1)) # individual log-posterior
#
# # item parameters
# # see ?itemparm
# itemparm(mod1) # item probabilities of success for each latent group
# itemparm(mod1, withSE = TRUE) # item probabilities of success & standard errors
# itemparm(mod1, what = "delta") # delta parameters
# itemparm(mod1, what = "delta",withSE=TRUE) # delta parameters
# itemparm(mod1, what = "gs") # guessing and slip parameters
# itemparm(mod1, what = "gs",withSE = TRUE) # guessing and slip parameters & standard errors
#
# # person parameters
# # see ?personparm
# personparm(mod1) # EAP estimates of attribute profiles
# personparm(mod1, what = "MAP") # MAP estimates of attribute profiles
# personparm(mod1, what = "MLE") # MLE estimates of attribute profiles
#
# #plot item response functions for item 10
# plotIRF(mod1,item = 10)
# plotIRF(mod1,item = 10,errorbar = TRUE) # with error bars
# plotIRF(mod1,item = c(6,10))
#
# # Use extract function to extract more components
# # See ?extract
#
# # ------- DINA model --------#
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# mod2 <- GDINA(dat = dat, Q = Q, model = "DINA")
# mod2
# itemparm(mod2, what = "gs") # guess and slip parameters
# itemparm(mod2, what = "gs",withSE = TRUE) # guess and slip parameters and standard errors
#
# # Model comparison at the test level via likelihood ratio test
# anova(mod1,mod2)
#
# # -------- DINO model -------#
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# mod3 <- GDINA(dat = dat, Q = Q, model = "DINO")
# #slip and guessing
# itemparm(mod3, what = "gs") # guess and slip parameters
# itemparm(mod3, what = "gs",withSE = TRUE) # guess and slip parameters + standard errors
#
# # Model comparison at test level via likelihood ratio test
# anova(mod1,mod2,mod3)
#
# # --------- ACDM model -------#
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# mod4 <- GDINA(dat = dat, Q = Q, model = "ACDM")
# mod4
# # --------- LLM model -------#
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# mod4b <- GDINA(dat = dat, Q = Q, model = "LLM")
# mod4b
# # --------- RRUM model -------#
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# mod4c <- GDINA(dat = dat, Q = Q, model = "RRUM")
# mod4c
#
# # --- Different CDMs for different items --- #
#
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# models <- c(rep("GDINA",3),"LLM","DINA","DINO","ACDM","RRUM","LLM","RRUM")
# mod5 <- GDINA(dat = dat, Q = Q, model = models)
# anova(mod1,mod2,mod3,mod4,mod4b,mod4c,mod5)
#
#
# ####################################
# # Example 2. #
# # Model estimations #
# # With monotonocity constraints #
# ####################################
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# # for item 10 only
# mod11 <- GDINA(dat = dat, Q = Q, model = "GDINA",mono.constraint = c(rep(FALSE,9),TRUE))
# mod11
# mod11a <- GDINA(dat = dat, Q = Q, model = "DINA",mono.constraint = TRUE)
# mod11a
# mod11b <- GDINA(dat = dat, Q = Q, model = "ACDM",mono.constraint = TRUE)
# mod11b
# mod11c <- GDINA(dat = dat, Q = Q, model = "LLM",mono.constraint = TRUE)
# mod11c
# mod11d <- GDINA(dat = dat, Q = Q, model = "RRUM",mono.constraint = TRUE)
# mod11d
# itemparm(mod11d,"delta")
# itemparm(mod11d,"rrum")
#
# ####################################
# # Example 3. #
# # Model estimations #
# # With Higher-order att structure #
# ####################################
#
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# # --- Higher order G-DINA model ---#
# mod12 <- GDINA(dat = dat, Q = Q, model = "DINA",
# att.dist="higher.order",higher.order=list(method="MMLE",nquad=31))
# hoest=hoparm(mod12) # extract higher-order parameters
# hoest$theta # ability
# hoest$lambda # structural parameters
# # --- Higher order DINA model ---#
# mod22 <- GDINA(dat = dat, Q = Q, model = "DINA",
# att.dist="higher.order",higher.order=list(method="BMLE"))
# # --- Higher order DINO model ---#
# mod23 <- GDINA(dat = dat, Q = Q, model = "DINO",att.dist="higher.order")
# # --- Higher order ACDM model ---#
# mod24 <- GDINA(dat = dat, Q = Q, model = "ACDM",
# att.dist="higher.order",higher.order=list(model="1PL"))
# # --- Higher order LLM model ---#
# mod25 <- GDINA(dat = dat, Q = Q, model = "LLM",att.dist="higher.order")
# # --- Higher order RRUM model ---#
# mod26 <- GDINA(dat = dat, Q = Q, model = "RRUM",att.dist="higher.order")
#
# ####################################
# # Example 4. #
# # Model estimations #
# # With user-specified att structure#
# ####################################
#
# # --- User-specified attribute priors ----#
# # prior distribution is fixed during calibration
# # Assume each of 000,100,010 and 001 has probability of 0.1
# # and each of 110, 101,011 and 111 has probability of 0.15
# # Note that the sum is equal to 1
# #
# prior <- c(0.1,0.1,0.1,0.1,0.15,0.15,0.15,0.15)
# # fit GDINA model with fixed prior dist.
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# modp1 <- GDINA(dat = dat, Q = Q, att.prior = prior, att.dist = "fixed")
# # See the posterior weights
# extract(modp1,what = "posterior.prob")
# extract(modp1,what = "att.prior")
# # ----Linear structure of attributes -----#
# # Assuming A1 -> A2 -> A3
# Q <- matrix(c(1,0,0,
# 1,0,0,
# 1,1,0,
# 1,1,0,
# 1,1,1,
# 1,1,1,
# 1,0,0,
# 1,0,0,
# 1,1,0,
# 1,1,0,
# 1,1,1,
# 1,1,1),ncol=3,byrow=TRUE)
# # item parameters for DINA model (guessing and slip)
# gs <- matrix(rep(0.1,24),ncol=2)
# N <- 5000
# # attribute simulation
# att <- rbind(matrix(0,nrow=500,ncol=3),
# matrix(rep(c(1,0,0),1000),ncol=3,byrow=TRUE),
# matrix(rep(c(1,1,0),1000),ncol=3,byrow=TRUE),
# matrix(rep(c(1,1,1),2500),ncol=3,byrow=TRUE))
# # data simulation
# simD <- simGDINA(N,Q,gs.parm = gs, model = "DINA",attribute = att)
# dat <- simD$dat
# # setting structure: A1 -> A2 -> A3
# # note: latent classes with prior 0 are assumed impossible
# prior <- c(0.1,0.2,0,0,0.2,0,0,0.5)
# out <- GDINA(dat, Q, att.prior = prior,att.str = TRUE, att.dist = "fixed", model = "DINA")
# # check posterior dist.
# extract(out,what = "posterior.prob")
# extract(out,what = "att.prior")
#
# out2 <- GDINA(dat, Q, att.prior = prior,att.str = TRUE, att.dist = "saturated",model = "DINA")
# # check posterior dist.
# extract(out2,what = "posterior.prob")
# extract(out2,what = "att.prior")
#
# ####################################
# # Example 5. #
# # Model estimations #
# # With user-specified att structure#
# ####################################
#
# # --- User-specified attribute structure ----#
# Q <- sim30GDINA$simQ
# K <- ncol(Q)
# # divergent structure A1->A2->A3;A1->A4->A5;A1->A4->A6
# diverg <- list(c(1,2),
# c(2,3),
# c(1,4),
# c(4,5))
# struc <- att.structure(diverg,K)
#
# # data simulation
# N <- 1000
# true.lc <- sample(c(1:2^K),N,replace=TRUE,prob=struc$att.prob)
# table(true.lc) #check the sample
# true.att <- attributepattern(K)[true.lc,]
# gs <- matrix(rep(0.1,2*nrow(Q)),ncol=2)
# # data simulation
# simD <- simGDINA(N,Q,gs.parm = gs,
# model = "DINA",attribute = true.att)
# dat <- extract(simD,"dat")
#
# modp1 <- GDINA(dat = dat, Q = Q, att.prior = struc$att.prob, att.str = TRUE, att.dist = "saturated")
# modp1
# # Note that fixed priors were used for all iterations
# extract(modp1,what = "att.prior")
# # Posterior weights were slightly different
# extract(modp1,what = "posterior.prob")
# modp2 <- GDINA(dat = dat, Q = Q, att.prior = struc$att.prob, att.str = TRUE, att.dist = "fixed")
# modp2
# extract(modp2,what = "att.prior")
# extract(modp2,what = "posterior.prob")
#
#
# ####################################
# # Example 6. #
# # Model estimations #
# # With user-specified initial pars #
# ####################################
#
# # check initials to see the format for initial item parameters
# initials <- sim10GDINA$simItempar
# dat <- sim10GDINA$simdat
# Q <- sim10GDINA$simQ
# mod.ini <- GDINA(dat,Q,catprob.parm = initials)
# extract(mod.ini,"initial.catprob")
# ####################################
# # Example 7. #
# # Model estimation #
# # Without M-step #
# ####################################
#
# # -----------Fix User-specified item parameters
# # Item parameters are not estimated
# # Only person attributes are estimated
# # attribute prior distribution matters if interested in the marginalized likelihood
# dat <- frac20$dat
# Q <- frac20$Q
# mod.initial <- GDINA(dat,Q,maxit=20) # estimation- only 10 iterations for illustration purposes
# par <- itemparm(mod.initial,digits=8)
# weights <- extract(mod.initial,"posterior.prob",digits=8) #posterior weights
# # use the weights as the priors
# mod.fix <- GDINA(dat,Q,catprob.parm = par,att.prior=c(weights),maxitr=0) # re-estimation
# anova(mod.initial,mod.fix) # very similar - good approximation most of time
# # prior used for the likelihood calculation for the last step
# priors <- extract(mod.initial,"att.prior")
# # use the priors as the priors
# mod.fix2 <- GDINA(dat,Q,catprob.parm = par,att.prior=priors,maxitr=0) # re-estimation
# anova(mod.initial,mod.fix2) # identical results
#
# ####################################
# # Example 8. #
# # polytomous attribute #
# # model estimation #
# # see Chen, de la Torre 2013 #
# ####################################
#
#
# # --- polytomous attribute G-DINA model --- #
# dat <- sim30pGDINA$simdat
# Q <- sim30pGDINA$simQ
# #polytomous G-DINA model
# pout <- GDINA(dat,Q)
#
# # ----- polymous DINA model --------#
# pout2 <- GDINA(dat,Q,model="DINA")
# anova(pout,pout2)
#
# ####################################
# # Example 9. #
# # Sequential G-DINA model #
# # see Ma, & de la Torre 2016 #
# ####################################
#
# # --- polytomous attribute G-DINA model --- #
# dat <- sim20seqGDINA$simdat
# Q <- sim20seqGDINA$simQ
# Q
# # Item Cat A1 A2 A3 A4 A5
# # 1 1 1 0 0 0 0
# # 1 2 0 1 0 0 0
# # 2 1 0 0 1 0 0
# # 2 2 0 0 0 1 0
# # 3 1 0 0 0 0 1
# # 3 2 1 0 0 0 0
# # 4 1 0 0 0 0 1
# # ...
#
# #sequential G-DINA model
# sGDINA <- GDINA(dat,Q,sequential = TRUE)
# sDINA <- GDINA(dat,Q,sequential = TRUE,model = "DINA")
# anova(sGDINA,sDINA)
# itemparm(sDINA) # processing function
# itemparm(sDINA,"itemprob") # success probabilities for each item
# itemparm(sDINA,"LCprob") # success probabilities for each category for all latent classes
#
# ####################################
# # Example 10. #
# # Multiple-Group G-DINA model #
# ####################################
# Q <- sim10GDINA$simQ
#
# # parameter simulation
# # Group 1 - female
# N1 <- 2000
# gs1 <- matrix(rep(0.1,2*nrow(Q)),ncol=2)
# # Group 2 - male
# N2 <- 2000
# gs2 <- matrix(rep(0.2,2*nrow(Q)),ncol=2)
#
# # data simulation for each group
# sim1 <- simGDINA(N1,Q,gs.parm = gs1,model = "DINA")
# sim2 <- simGDINA(N2,Q,gs.parm = gs2,model = "DINO")
#
# # combine data
# # see ?bdiagMatrix
# dat <- bdiagMatrix(list(extract(sim1,"dat"),extract(sim2,"dat")),fill=NA)
# Q <- rbind(Q,Q)
# gr <- rep(c("female","male"),each=2000)
# # Fit G-DINA model
# mg.est <- GDINA(dat = dat,Q = Q,group = gr)
# summary(mg.est)
# extract(mg.est,"posterior.prob")
#
# # Fit G-DINA model with different joint attribute dist.
# mg.est2 <- GDINA(dat = dat,Q = Q,group = gr,
# att.dist = c("saturated","fixed"))
# summary(mg.est2)
#
## ---------------------------------------------
Run the code above in your browser using DataLab