# 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
# structural parameters
# see ?coef
coef(mod1) # item probabilities of success for each latent group
coef(mod1, withSE = TRUE) # item probabilities of success & standard errors
coef(mod1, what = "delta") # delta parameters
coef(mod1, what = "delta",withSE=TRUE) # delta parameters
coef(mod1, what = "gs") # guessing and slip parameters
coef(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
coef(mod2, what = "gs") # guess and slip parameters
coef(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
coef(mod3, what = "gs") # guess and slip parameters
coef(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
coef(mod11d,"delta")
coef(mod11d,"rrum")
####################################
# Example 3a. #
# 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(nquad=31,model = "2PL"))
personparm(mod12,"HO") # higher-order ability
# structural parameters
# first column is slope and the second column is intercept
coef(mod12,"lambda")
# --- Higher order DINA model ---#
mod22 <- GDINA(dat = dat, Q = Q, model = "DINA", att.dist="higher.order",
higher.order=list(model = "2PL",Prior=TRUE))
####################################
# Example 3b. #
# Model estimations #
# With log-linear att structure #
####################################
# --- DINA model with loglinear smoothed attribute space ---#
dat <- sim10GDINA$simdat
Q <- sim10GDINA$simQ
mod23 <- GDINA(dat = dat, Q = Q, model = "DINA",att.dist="loglinear",loglinear=1)
coef(mod23,"lambda") # intercept and three main effects
####################################
# Example 3c. #
# Model estimations #
# With independent att structure #
####################################
# --- GDINA model with independent attribute space ---#
dat <- sim10GDINA$simdat
Q <- sim10GDINA$simQ
mod33 <- GDINA(dat = dat, Q = Q, att.dist="independent")
coef(mod33,"lambda") # mastery probability for each attribute
####################################
# 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
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
# prior dist.
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
# estimation- only 20 iterations for illustration purposes
mod.initial <- GDINA(dat,Q,control = list(maxitr=20))
par <- coef(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),control = list(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, control = list(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)
coef(sDINA) # processing function
coef(sDINA,"itemprob") # success probabilities for each item
coef(sDINA,"LCprob") # success probabilities for each category for all latent classes
####################################
# Example 10a. #
# Multiple-Group G-DINA model #
####################################
Q <- sim10GDINA$simQ
K <- ncol(Q)
# parameter simulation
# Group 1 - female
N1 <- 3000
gs1 <- matrix(rep(0.1,2*nrow(Q)),ncol=2)
# Group 2 - male
N2 <- 3000
gs2 <- matrix(rep(0.2,2*nrow(Q)),ncol=2)
# data simulation for each group
sim1 <- simGDINA(N1,Q,gs.parm = gs1,model = "DINA",att.dist = "higher.order",
higher.order.parm = list(theta = rnorm(N1),
lambda = data.frame(a=rep(1.5,K),b=seq(-1,1,length.out=K))))
sim2 <- simGDINA(N2,Q,gs.parm = gs2,model = "DINO",att.dist = "higher.order",
higher.order.parm = list(theta = rnorm(N2),
lambda = data.frame(a=rep(1,K),b=seq(-2,2,length.out=K))))
# combine data
# see ?bdiagMatrix
dat <- bdiagMatrix(list(extract(sim1,"dat"),extract(sim2,"dat")),fill=NA)
Q <- rbind(Q,Q)
gr <- rep(c(1,2),c(3000,3000))
# Fit G-DINA model
mg.est <- GDINA(dat = dat,Q = Q,group = gr,att.dist="higher.order",
higher.order=list(model = "1PL",Prior=TRUE))
summary(mg.est)
extract(mg.est,"posterior.prob")
####################################
# Example 10b. #
# Multiple-Group G-DINA model #
####################################
Q <- sim30GDINA$simQ
K <- ncol(Q)
# parameter simulation
N1 <- 3000
gs1 <- matrix(rep(0.1,2*nrow(Q)),ncol=2)
N2 <- 3000
gs2 <- matrix(rep(0.2,2*nrow(Q)),ncol=2)
# data simulation for each group
# two groups have different theta distributions
sim1 <- simGDINA(N1,Q,gs.parm = gs1,model = "DINA",att.dist = "higher.order",
higher.order.parm = list(theta = rnorm(N1),
lambda = data.frame(a=rep(1,K),b=seq(-2,2,length.out=K))))
sim2 <- simGDINA(N2,Q,gs.parm = gs2,model = "DINO",att.dist = "higher.order",
higher.order.parm = list(theta = rnorm(N2,1,1),
lambda = data.frame(a=rep(1,K),b=seq(-2,2,length.out=K))))
# combine data
# see ?bdiagMatrix
dat <- bdiagMatrix(list(extract(sim1,"dat"),extract(sim2,"dat")),fill=NA)
Q <- rbind(Q,Q)
gr <- rep(c(1,2),c(3000,3000))
# Fit G-DINA model
mg.est <- GDINA(dat = dat,Q = Q,group = gr,att.dist="higher.order",
higher.order=list(model = "Rasch",Prior=FALSE,Type = "SameLambda"))
summary(mg.est)
coef(mg.est,"lambda")
####################################
# Example 11. #
# Bug DINO model #
####################################
set.seed(123)
Q <- sim10GDINA$simQ # 1 represents misconceptions/bugs
ip <- list(
c(0.8,0.2),
c(0.7,0.1),
c(0.9,0.2),
c(0.9,0.1,0.1,0.1),
c(0.9,0.1,0.1,0.1),
c(0.9,0.1,0.1,0.1),
c(0.9,0.1,0.1,0.1),
c(0.9,0.1,0.1,0.1),
c(0.9,0.1,0.1,0.1),
c(0.9,0.1,0.1,0.1,0.1,0.1,0.1,0.1))
sim <- simGDINA(N=1000,Q=Q,catprob.parm = ip)
dat <- extract(sim,"dat")
# use latent.var to specify a bug model
est <- GDINA(dat=dat,Q=Q,latent.var="bugs",model="DINO")
coef(est)
####################################
# Example 12. #
# Bug DINA model #
####################################
set.seed(123)
Q <- sim10GDINA$simQ # 1 represents misconceptions/bugs
ip <- list(
c(0.8,0.2),
c(0.7,0.1),
c(0.9,0.2),
c(0.9,0.9,0.9,0.1),
c(0.9,0.9,0.9,0.1),
c(0.9,0.9,0.9,0.1),
c(0.9,0.9,0.9,0.1),
c(0.9,0.9,0.9,0.1),
c(0.9,0.9,0.9,0.1),
c(0.9,0.9,0.9,0.9,0.9,0.9,0.9,0.1))
sim <- simGDINA(N=1000,Q=Q,catprob.parm = ip)
dat <- extract(sim,"dat")
# use latent.var to specify a bug model
est <- GDINA(dat=dat,Q=Q,latent.var="bugs",model="DINA")
coef(est)
####################################
# Example 13a. #
# user specified design matrix #
# LCDM (logit G-DINA) #
####################################
dat <- sim30GDINA$simdat
Q <- sim30GDINA$simQ
#find design matrix for each item => must be a list
D <- lapply(rowSums(Q),designmatrix,model="GDINA")
# for comparison, use change in -2LL as convergence criterion
# LCDM
lcdm <- GDINA(dat = dat, Q = Q, model = "UDF", design.matrix = D,
linkfunc = "logit", control=list(conv.type="neg2LL"),solver="slsqp")
# identity link GDINA
iGDINA <- GDINA(dat = dat, Q = Q, model = "GDINA",
control=list(conv.type="neg2LL"),solver="slsqp")
# compare two models => identical
anova(lcdm,iGDINA)
####################################
# Example 13b. #
# user specified design matrix #
# RRUM #
####################################
dat <- sim30GDINA$simdat
Q <- sim30GDINA$simQ
# specify design matrix for each item => must be a list
# D can be defined by the user
D <- lapply(rowSums(Q),designmatrix,model="ACDM")
# for comparison, use change in -2LL as convergence criterion
# RRUM
logACDM <- GDINA(dat = dat, Q = Q, model = "UDF", design.matrix = D,
linkfunc = "log", control=list(conv.type="neg2LL"),solver="slsqp")
# identity link GDINA
RRUM <- GDINA(dat = dat, Q = Q, model = "RRUM",
control=list(conv.type="neg2LL"),solver="slsqp")
# compare two models => identical
anova(logACDM,RRUM)
####################################
# Example 14. #
# Multiple-strategy DINA model #
####################################
Q <- matrix(c(1,1,1,1,0,
1,2,0,1,1,
2,1,1,0,0,
3,1,0,1,0,
4,1,0,0,1,
5,1,1,0,0,
5,2,0,0,1),ncol = 5,byrow = TRUE)
d <- list(
item1=c(0.2,0.7),
item2=c(0.1,0.6),
item3=c(0.2,0.6),
item4=c(0.2,0.7),
item5=c(0.1,0.8))
set.seed(12345)
sim <- simGDINA(N=1000,Q = Q, delta.parm = d,
model = c("MSDINA","MSDINA","DINA",
"DINA","DINA","MSDINA","MSDINA"))
# simulated data
dat <- extract(sim,what = "dat")
# estimation
# MSDINA need to be specified for each strategy
est <- GDINA(dat,Q,model = c("MSDINA","MSDINA","DINA",
"DINA","DINA","MSDINA","MSDINA"))
coef(est,"delta")
# }
# NOT RUN {
# }
Run the code above in your browser using DataLab