## File Name: tam_mml_ic.R
## File Version: 9.208
#--- information criteria
tam_mml_ic <- function( nstud, deviance, xsi, xsi.fixed,
beta, beta.fixed, ndim, variance.fixed, G, irtmodel,
B_orig=NULL, B.fixed, E, est.variance, resp,
est.slopegroups=NULL, variance.Npars=NULL, group, penalty_xsi=0,
AXsi=NULL, pweights=NULL, resp.ind=NULL, B=NULL)
{
#--- log likelihood and log prior
deviance <- deviance - penalty_xsi
loglike <- - deviance / 2
logprior <- - penalty_xsi / 2
logpost <- loglike + logprior
#***Model parameters
ic <- data.frame("n"=nstud, "deviance"=deviance )
ic$loglike <- loglike
ic$logprior <- logprior
ic$logpost <- logpost
dev <- deviance
# xsi parameters
ic$Nparsxsi <- length(xsi)
if ( ! is.null( xsi.fixed) ){
ic$Nparsxsi <- ic$Nparsxsi - nrow(xsi.fixed )
}
#-- B slopes
if (!is.null(AXsi)){
maxKi <- rowSums( ! is.na(AXsi) ) - 1
NB <- dim(B_orig)[3]
NparsB <- 0
for (dd in 1:NB){
NparsB <- NparsB + sum( maxKi*(rowSums(B_orig[,,dd]) > 0) )
}
}
ic$NparsB <- 0
if ( irtmodel=="2PL" ){
ic$NparsB <- sum( B_orig !=0 )
}
if ( irtmodel=="GPCM" ){
B1 <- B[,2,]
NparsB <- sum(B1!=0)
ic$NparsB <- NparsB
}
if ( irtmodel=="GPCM.groups" ){
ic$NparsB <- length( unique( setdiff(est.slopegroups,0)))
# not correctly implemented in general setting
}
if ( irtmodel=="GPCM.design" ){
ic$NparsB <- ncol(E)
}
if ( irtmodel=="2PL.groups" ){
ic$NparsB <- length( unique( setdiff(est.slopegroups,0 ) ))
# This is not yet correct for multiple dimensions and multiple
# categories
}
if ( ! is.null( B.fixed ) ){
nB <- nrow(B.fixed)
if ( irtmodel=="GPCM" ){
nB <- length(B.fixed[ B.fixed[,2]==2,1])
}
ic$NparsB <- max(ic$NparsB - nB, 0)
}
# beta regression parameters
ic$Nparsbeta <- dim(beta)[1] * dim(beta)[2]
if ( ! is.null( beta.fixed) ){
ic$Nparsbeta <- ic$Nparsbeta - nrow(beta.fixed )
}
# variance/covariance matrix
ic$Nparscov <- ndim + ndim*(ndim-1)/2
if ( ! est.variance ){
ic$Nparscov <- ic$Nparscov - ndim
}
if ( ! is.null(variance.fixed) ){
ic$Nparscov <- max(0, ic$Nparscov - nrow(variance.fixed ) )
}
if ( ! is.null(variance.Npars) ){
ic$Nparscov <- variance.Npars
}
if ( ! is.null(group) ){
ic$Nparscov <- ic$Nparscov + length( unique(group) ) - 1
}
# total number of parameters
ic$Npars <- ic$np <- ic$Nparsxsi + ic$NparsB + ic$Nparsbeta + ic$Nparscov
#- compute total number of observations
ic$ghp_obs <- tam_ghp_number_informations(pweights=pweights, resp.ind=resp.ind)
#--- calculate all criteria
ic <- tam_mml_ic_criteria(ic=ic)
#--- OUTPUT
return(ic)
}
.TAM.ic <- tam_mml_ic
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.