#' extract various components from model comparison
#'
#' @description
#'NULL
#'
#'
#' @param object object of class \code{modelcomp} for various S3 methods
#' @param what argument for S3 method \code{extract} indicating what to extract;
#' It can be \code{"wald"} for wald statistics, \code{"wald.p"} for associated p-values,
#' \code{"df"} for degrees of freedom,
#' and \code{"DS"} for dissimilarity between G-DINA and other CDMs.
#' @param digits How many decimal places in each number? The default is 4.
#' @param ... additional arguments
#'
#' @include modelcomp.R
#'@describeIn modelcomp extract various elements from \code{modelcomp} objects
#'@aliases extract.modelcomp
#'@export
extract.modelcomp <- function(object,
what=c("stats","pvalues","adj.pvalues","df","DS","selected.model"), digits = 4,...){
# what <- match.arg(what)
out <- switch(what,
stats = round(object$stats,digits),
pvalues = round(object$pvalues,digits),
adj.pvalues = round(object$adj.pvalues,digits),
df = object$df,
DS = {if(is.null(object$DS)) NULL else round(object$DS,digits)},
models = object$models,
decision.args = object$decision.args,
selected.model = object$selected.model,
Wald.args = object$Wald.args,
LR.args = object$LR.args,
neg2LL = object$neg2LL)
return(out)
}
#'@include itemfit.R
#' extract elements from objects of itemfit
#'
#'@description NULL
#'
#' @param object objects of class \code{itemfit} for various S3 methods
#' @param what argument for S3 method \code{extract} indicating what to extract;
#' It can be \code{"p"} for proportion correct statistics,
#' \code{"r"} for transformed correlations, \code{logOR} for log odds ratios and
#' \code{"maxitemfit"} for maximum statistics for each item.
#' @param ... additional arguments
#'
#' @describeIn itemfit extract various elements from \code{itemfit} objects
#' @aliases extract.itemfit
#' @export
extract.itemfit <- function(object,what,...){
out <- switch(what,
r = object$r,
p = object$p,
logOR = object$logOR,
maxitemfit = object$max.itemlevel.fit,
person.sim = object$options$person.sim,
p.adjust.method = object$options$p.adjust.methods,
N.resampling = object$options$N.resampling,
randomseed = object$options$randomseed,
call = object$options$call,
digits = object$options$digits,
stop(sprintf("Can not extract element \'%s\'", what), call.=FALSE))
return(out)
}
#'@title NULL
#'
#'@description NULL
#'
#' @param object \code{Qval} objects for S3 methods
#' @param what argument for S3 method \code{extract} indicating what to extract;
#' It can be \code{"sug.Q"} for suggested Q-matrix,
#' \code{"Q"} for original Q-matrix, \code{"varsigma"} for varsigma index,
#' and \code{"PVAF"} for PVAF.
#' @param ... additional arguments
#'
#'@describeIn Qval extract various elements from \code{Qval} objects
#'@aliases extract.Qval
#'@export
extract.Qval <- function(object,
what=c("sug.Q","varsigma","PVAF","eps","Q"),...){
what <- match.arg(what)
out <- switch(what,
sug.Q = object$sug.Q,
varsigma = object$varsigma,
PVAF = object$PVAF,
eps = object$eps,
Q = object$Q,
stop(sprintf("Can not show element \'%s\'", what), call.=FALSE))
return(out)
}
#' @include simGDINA.R
#' @title NULL
#'
#' @description NULL
#'
#' @param object object of class \code{simGDINA} for method \code{extract}
#' @param what argument for S3 method \code{extract} indicating what to extract
#' @param ... additional arguments
#' @rdname simGDINA
#'@export
extract.simGDINA <- function(object,
what=c("dat","Q","attribute","catprob.parm",
"delta.parm","higher.order.parm","mvnorm.parm",
"LCprob.parm"),...){
out <- switch(tolower(what),
dat = object$dat,
Q = object$Q,
attribute = object$attribute,
catprob.parm = object$catprob.parm,
delta.parm = object$delta.parm,
higher.order.parm = object$higher.order.parm,
mvnorm.parm = object$mvnorm.parm,
LCprob.parm = object$LCprob.parm,
stop(sprintf("Can not show element \'%s\'", what), call.=FALSE))
return(out)
}
#'@include GDINA-package.R GDINA.R
#'@title extract elements from objects of various classes
#'
#' @description A generic function to extract elements from objects of class \code{GDINA},
#' \code{itemfit}, \code{modelcomp}, \code{Qval} or \code{simGDINA}. This
#' page gives the elements that can be extracted from the class \code{GDINA}.
#' To see what can be extracted from \code{\link{itemfit}}, \code{\link{modelcomp}}, and
#' \code{\link{Qval}}, go to the corresponding function help page.
#'
#' Objects which can be extracted from \code{GDINA} objects include:
#'
#' \describe{
#' \item{AIC}{AIC}
#' \item{att.prior}{attribute prior weights for calculating marginalized likelihood in the last EM iteration}
#' \item{attributepattern}{all attribute patterns involved in the current calibration}
#' \item{BIC}{BIC}
#' \item{CAIC}{Consistent AIC}
#' \item{catprob.cov}{covariance matrix of item probability parameter estimates; Need to specify \code{SE.type}}
#' \item{catprob.parm}{item parameter estimates}
#' \item{catprob.se}{standard error of item probability parameter estimates; Need to specify \code{SE.type}}
#' \item{convergence}{\code{TRUE} if the calibration is converged.}
#' \item{dat}{raw data}
#' \item{del.ind}{deleted observation number}
#' \item{delta.cov}{covariance matrix of delta parameter estimates; Need to specify \code{SE.type}}
#' \item{delta.parm}{delta parameter estimates}
#' \item{delta.se}{standard error of delta parameter estimates; Need to specify \code{SE.type}}
#' \item{designmatrix}{A list of design matrices for each item/category}
#' \item{deviance}{deviance, or negative two times observed marginal log likelihood}
#' \item{discrim}{GDINA discrimination index}
#' \item{expectedCorrect}{expected # of examinees in each latent group answering item correctly}
#' \item{expectedTotal}{expected # of examinees in each latent group}
#' \item{higher.order}{higher-order model specifications}
#' \item{LCprob.parm}{success probabilities for all latent classes}
#' \item{logLik}{observed marginal log likelihood}
#' \item{linkfunc}{link functions for each item}
#' \item{initial.catprob}{initial item category probability parameters}
#' \item{natt}{number of attributes}
#' \item{ncat}{number of categories}
#' \item{ngroup}{number of groups}
#' \item{nitem}{number of items}
#' \item{nitr}{number of EM iterations}
#' \item{nobs}{number of observations, or sample size}
#' \item{nLC}{number of latent classes}
#' \item{prevalence}{prevalence of each attribute}
#' \item{posterior.prob}{posterior weights for each latent class}
#' \item{reduced.LG}{Reduced latent group for each item}
#' \item{SABIC}{Sample size Adusted BIC}
#' \item{sequential}{is a sequential model fitted?}
#' }
#' @param object objects from class \code{GDINA},\code{itemfit}, \code{modelcomp}, \code{Qval} or \code{simGDINA}
#' @param what what to extract
#' @param ... additional arguments
#'
#' @examples
#'\dontrun{
#' dat <- sim10GDINA$simdat
#' Q <- sim10GDINA$simQ
#' fit <- GDINA(dat = dat, Q = Q, model = "GDINA")
#' extract(fit,"discrim")
#' extract(fit,"designmatrix")
#' }
#'
#' @export
extract <- function (object, what, ...) {
UseMethod("extract")
}
#' @title extract elements of GDINA estimates
#' @description
#' NULL
#' @details
#' NULL
#' @describeIn GDINA extract various elements of GDINA estimates
#' @aliases extract.GDINA
#' @export
extract.GDINA <- function(object,what,SE.type = 2,...){
out <- switch(what,
AIC = object$testfit$AIC,
attributepattern = object$attributepattern,
att.dist = object$options$att.dist,
att.prior = object$options$att.prior,
att.str=object$options$att.str,
BIC = object$testfit$BIC,
CAIC = object$testfit$CAIC,
catprob.cov = {var <- OPG_p(object, SE.type = SE.type);list(cov = var$cov, index = var$ind)},
catprob.matrix = object$catprob.matrix,
catprob.parm = object$catprob.parm,
catprob.se = {
if(SE.type==3&any(extract(object,"att.dist")=="higher.order"))
stop("standard error cannot be calculated.",call. = FALSE)
Kj <- extract(object, "Kj")
se <- OPG_p(object, SE.type = SE.type, ...)$se
for (j in 1:length(se))
names(se[[j]]) <- paste0("SE[P(", apply(extract(object,"reduced.LG")[[j]], 1, paste0, collapse = ""), ")]")
names(se) <- object$options$item.names
se
},
convergence={if(all(object$control$vmaxitr==object$control$maxitr)&object$control$maxitr>0){
object$options$itr<object$control$maxitr
}else{
NA # to be determined
}},
dat = object$options$dat, #raw data
del.ind = object$technicals$del.ind,
delta.cov = {var <- OPG_d(object, SE.type = SE.type);list(cov = var$cov, index = var$ind)},
delta.parm = {format_delta(object$delta.parm,model = object$options$model,
Kj = rowSums(extract(object,"Q")),
item.names = object$options$item.names,
digits = 10)},
delta.se = {
if(SE.type==3&any(extract(object,"att.dist")=="higher.order"))stop("standard error cannot be calculated.",call. = FALSE)
format_delta(delta = OPG_d(object,SE.type = SE.type)$se,
model = object$options$model, Kj = rowSums(extract(object,"Q")),
item.names = object$options$item.names,digits = 10)
},
designmatrix = object$options$DesignMatrices,
deviance = object$testfit$Deviance,
dif.p = object$options$dif.p,
dif.prior = object$options$dif.prior,
dif.LL = object$options$dif.LL,
discrim = {
if(object$options$no.group>1) stop("Discrimination indices are not available for multiple group model.",call. = FALSE)
dj <- sapply(object$catprob.parm, function(x) x[length(x)]-x[1])
wp <- t(object$LC.prob) * c(object$posterior.prob) #L x J w*p matrix
gdi <- colSums(t((object$LC.prob - colSums(wp))^2) * c(object$posterior.prob)) # vector of length J
Discrim <- data.frame(dj = dj, GDI = gdi)
rownames(Discrim) <- object$options$item.names
colnames(Discrim) <- c("P(1)-P(0)", "GDI")
Discrim
},
eta = object$technicals$eta,
end.time = object$extra$end.time,
expectedCorrect = object$technicals$expectedCorrect,
expectedTotal = object$technicals$expectedTotal,
expectedCorrect.LC = {
if(object$options$sequential){
dat <- extract(object,"seq.dat")
}else{
dat <- extract(object,"dat")
}
out <- NgRg(as.matrix(object$technicals$logposterior.i),
as.matrix(dat),
eta(matrix(1,extract(object,"ncat"),extract(object,"natt"))),
rep(1,nrow(dat)))$Rg
row.names(out) <- object$options$item.names
out
},
expectedTotal.LC = {
if(object$options$sequential){
dat <- extract(object,"seq.dat")
}else{
dat <- extract(object,"dat")
}
out <- NgRg(as.matrix(object$technicals$logposterior.i),
as.matrix(dat),
eta(matrix(1,extract(object,"ncat"),extract(object,"natt"))),
rep(1,nrow(dat)))$Ng
row.names(out) <- object$options$item.names
out
},
group = object$options$group,
gr = object$options$gr,
group.label = object$options$group.label,
higher.order = object$options$higher.order,
initial.catprob = object$technical$initial.parm,
itemprob.parm = {
if(object$options$sequential){
itemparm <- object$catprob.matrix
p <- vector("list",extract(object,"nitem"))
Qc <- extract(object,"Qc")
for (j in 1:extract(object,"nitem")){
Qj <- Qc[which(Qc[,1]==j),-c(1:2),drop=FALSE]
itemparj <- itemparm[which(Q[,1]==j),,drop=FALSE]
colj <- which(colSums(Qj)>0)
Kj <- length(colj)
if(nrow(Qj)>1){ # polytomous items
if(length(colj)>1){
redQj <- Qj[,colj]
pj <- sj <- rbind(uP(as.matrix(eta(as.matrix(redQj))), as.matrix(itemparj)),0)
}else{
pj <- sj <- rbind(itemparj[,1:2],0)
}
for (s in 1:(nrow(sj)-1)){
pj[s,] <- apply(sj[1:s,,drop=FALSE],2,prod)*(1-sj[s+1,])
}
pj <- pj[-nrow(pj),]
}else{ #dichotomous items
pj <- matrix(itemparj[!is.na(itemparj)],nrow = 1)
}
colnames(pj) <- paste0("P(", apply(attributepattern(Kj), 1, paste0, collapse = ""), ")")
rownames(pj) <- paste("Cat",1:nrow(pj))
p[[j]] <- pj
}
names(p) <- paste("Item",1:extract(object,"nitem"))
p
}else{
object$catprob.parm
}
},
itemprob.se = { if(!object$options$sequential) extract(object,"catprob.se",SE.type = SE.type)},
item.names = object$options$item.names,
itemprob.history = object$diagnos$itemprob.matrix,
Kj = {rowSums(extract(object,"Q"))},
no.bugs = object$options$no.bugs,
LC.labels = object$technicals$LC.labels,
LCprob.parm = object$LC.prob,
LCpf.parm = {
# LCpf <- patt <- eta(as.matrix(extract(object,"Q")))
LCpf <- patt <- eta(as.matrix(extract(object,"Q")), object$attributepattern)
pf <- extract(object,"catprob.parm")
for(i in 1:nrow(LCpf)) LCpf[i,] <- pf[[i]][patt[i,]]
rownames(LCpf) <- rownames(extract(object,"catprob.matrix"))
LCpf
},# processing function
linkfunc = {c("identity","logit","log")[object$options$linkfunc]},
logLik = -0.5*object$testfit$Deviance,
logposterior.i = object$technicals$logposterior.i,
loglikelihood.i = object$technicals$loglikelihood.i,
loglinear = object$options$loglinear,
models = object$model,
models_numeric = object$options$model,
mono.constraint = object$options$mono.constraint,
npar = object$technicals$total.npar,
npar.item = object$technicals$free.item.npar, # free parameters
npar.fixeditem = object$technicals$total.item.npar - object$technicals$free.item.npar,
npar.att = object$technicals$stru.npar,
natt = ncol(extract(object,"Q")),
ncat = nrow(extract(object,"Q")),
ngroup = object$options$no.group,
nitem = ncol(object$options$dat),
nitr = object$options$itr,
nobs = nrow(object$options$dat),
nLC = dim(object$attributepattern)[1],
pf = object$pf,
posterior.prob = object$posterior.prob,
prevalence = {
Q <- extract(object,"Q")
preva <- vector("list",extract(object,"ngroup"))
pattern <- extract(object,"attributepattern")
for(g in 1:extract(object,"ngroup")) {
for (i in c(0:max(Q))) {
preva[[g]] <-
cbind(preva[[g]], t(extract(object,"posterior.prob")[g,,drop=FALSE] %*% ((pattern == i) * 1)))
}
colnames(preva[[g]]) <- paste0("Level", 0:max(Q))
rownames(preva[[g]]) <- paste0("A", seq(1, ncol(Q)))
}
names(preva) <- object$options$group.label
preva
},
Q = {
if (object$options$sequential||any(extract(object,"models")=="MSDINA")) {
out <- data.frame(object$options$Qm)
} else{
out <- data.frame(object$options$Q)
}
colnames(out) <- paste0("A", 1:ncol(out))
out
},
Qc = {if(any(extract(object,"models")=="MSDINA")){
object$options$Qcm
} else if(object$options$sequential) {
data.frame(object$options$Q)
} else{
out <- data.frame(
Item = 1:nrow(object$options$Q),
Cat = rep(1, nrow(object$options$Q)),
object$options$Q
)
colnames(out)[-c(1:2)] <-
paste0("A", 1:ncol(object$options$Q))
out
}},
reduced.LG = object$technicals$reduced.LG,
originalQ = object$options$Q,
SABIC = object$testfit$SABIC,
start.time = object$extra$start.time,
time = object$extra$timeused,
sequential = object$options$sequential,
seq.dat = object$options$seq.dat,
struc.parm = object$struc.parm,
verbose = object$options$verbose,
call = object$extra$call,
stop(sprintf("Can not extract element \'%s\'", what), call.=FALSE))
return(out)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.