sirt (version 1.14-0)

prob.guttman: Probabilistic Guttman Model

Description

This function estimates the probabilistic Guttman model which is a special case of an ordered latent trait model (Hanson, 2000; Proctor, 1970).

Usage

prob.guttman(dat, pid = NULL, guess.equal = FALSE,  slip.equal = FALSE,  itemlevel = NULL, conv1 = 0.001, glob.conv = 0.001, mmliter = 500) 
"summary"(object,...)
"anova"(object,...)
"logLik"(object,...)
"IRT.irfprob"(object,...)
"IRT.likelihood"(object,...)
"IRT.posterior"(object,...)

Arguments

dat

An $N \times I$ data frame of dichotomous item responses

pid

Optional vector of person identifiers

guess.equal

Should the same guessing parameters for all the items estimated?

slip.equal

Should the same slipping parameters for all the items estimated?

itemlevel

A vector of item levels of the Guttman scale for each item. If there are $K$ different item levels, then the Guttman scale possesses $K$ ordered trait values.

conv1

Convergence criterion for item parameters

glob.conv

Global convergence criterion for the deviance

mmliter

Maximum number of iterations

object

Object of class prob.guttman

...

Further arguments to be passed

Value

An object of class prob.guttman

References

Hanson, B. (2000). IRT parameter estimation using the EM algorithm. Technical Report.

Proctor, C. H. (1970). A probabilistic formulation and statistical analysis for Guttman scaling. Psychometrika, 35, 73-78.

Examples

#############################################################################
# EXAMPLE 1: Dataset Reading
#############################################################################
data(data.read)
dat <- data.read

#***
# Model 1: estimate probabilistic Guttman model
mod1 <- prob.guttman( dat )
summary(mod1)

#***
# Model 2: probabilistic Guttman model with equal guessing and slipping parameters
mod2 <- prob.guttman( dat , guess.equal=TRUE , slip.equal=TRUE)
summary(mod2)

#***
# Model 3: Guttman model with three a priori specified item levels
itemlevel <- rep(1,12)
itemlevel[ c(2,5,8,10,12) ] <- 2
itemlevel[ c(3,4,6) ] <- 3
mod3 <- prob.guttman( dat , itemlevel=itemlevel )
summary(mod3)

## Not run: 
# #***
# # Model3m: estimate Model 3 in mirt
# 
# library(mirt)
# # define four ordered latent classes
# Theta <- scan(nlines=1)
#    0 0 0    1 0 0   1 1 0   1 1 1    
# Theta <- matrix( Theta , nrow=4 , ncol=3,byrow=TRUE)
# 
# # define mirt model
# I <- ncol(dat)  # I = 12
# mirtmodel <- mirt::mirt.model("
#         # specify factors for each item level
#         C1 = 1,7,9,11
#         C2 = 2,5,8,10,12 
#         C3 = 3,4,6
#         ")
# # get initial parameter values
# mod.pars <- mirt::mirt(dat, model=mirtmodel ,  pars = "values")
# # redefine initial parameter values
# mod.pars[ mod.pars$name == "d" ,"value" ]  <- -1
# mod.pars[ mod.pars$name %in% paste0("a",1:3) & mod.pars$est    ,"value" ]  <- 2
# mod.pars
# # define prior for latent class analysis
# lca_prior <- function(Theta,Etable){
#   # number of latent Theta classes
#   TP <- nrow(Theta)
#   # prior in initial iteration
#   if ( is.null(Etable) ){ prior <- rep( 1/TP , TP ) }    
#   # process Etable (this is correct for datasets without missing data)
#   if ( ! is.null(Etable) ){  
#     # sum over correct and incorrect expected responses 
#     prior <- ( rowSums(Etable[ , seq(1,2*I,2)]) + rowSums(Etable[,seq(2,2*I,2)]) )/I
#                  }
#   prior <- prior / sum(prior)  
#   return(prior)
# }
# # estimate model in mirt
# mod3m <- mirt::mirt(dat, mirtmodel , pars = mod.pars , verbose=TRUE , 
#             technical = list( customTheta=Theta , customPriorFun = lca_prior) )
# # correct number of estimated parameters
# mod3m@nest <- as.integer(sum(mod.pars$est) + nrow(Theta)-1 ) 
# # extract log-likelihood and compute AIC and BIC
# mod3m@logLik
# ( AIC <- -2*mod3m@logLik+2*mod3m@nest )
# ( BIC <- -2*mod3m@logLik+log(mod3m@Data$N)*mod3m@nest )
# # compare with information criteria from prob.guttman
# mod3$ic
# # model fit in mirt
# mirt::M2(mod3m)
# # extract coefficients
# ( cmod3m <- mirt.wrapper.coef(mod3m) )
# # compare estimated distributions
# round( cbind( "sirt"  = mod3$trait$prob , "mirt" = mod3m@Prior[[1]] ) , 5 )
#   ##           sirt    mirt
#   ##   [1,] 0.13709 0.13765
#   ##   [2,] 0.30266 0.30303
#   ##   [3,] 0.15239 0.15085
#   ##   [4,] 0.40786 0.40846
# # compare estimated item parameters
# ipars <- data.frame( "guess.sirt" = mod3$item$guess , 
#                      "guess.mirt" = plogis( cmod3m$coef$d ) )
# ipars$slip.sirt <- mod3$item$slip
# ipars$slip.mirt <- 1-plogis( rowSums(cmod3m$coef[ , c("a1","a2","a3","d") ] ) )
# round( ipars , 4 )
#   ##      guess.sirt guess.mirt slip.sirt slip.mirt
#   ##   1      0.7810     0.7804    0.1383    0.1382
#   ##   2      0.4513     0.4517    0.0373    0.0368
#   ##   3      0.3203     0.3200    0.0747    0.0751
#   ##   4      0.3009     0.3007    0.3082    0.3087
#   ##   5      0.5776     0.5779    0.1800    0.1798
#   ##   6      0.3758     0.3759    0.3047    0.3051
#   ##   7      0.7262     0.7259    0.0625    0.0623
#   ##   [...]
# 
# #***
# # Model 4: Monotone item response function estimated in mirt
# 
# # define four ordered latent classes
# Theta <- scan(nlines=1)
#    0 0 0    1 0 0   1 1 0   1 1 1    
# Theta <- matrix( Theta , nrow=4 , ncol=3,byrow=TRUE)
# 
# # define mirt model
# I <- ncol(dat)  # I = 12
# mirtmodel <- mirt::mirt.model("
#         # specify factors for each item level
#         C1 = 1-12  
#         C2 = 1-12
#         C3 = 1-12
#         ")
# # get initial parameter values
# mod.pars <- mirt::mirt(dat, model=mirtmodel ,  pars = "values")
# # redefine initial parameter values
# mod.pars[ mod.pars$name == "d" ,"value" ]  <- -1
# mod.pars[ mod.pars$name %in% paste0("a",1:3) & mod.pars$est    ,"value" ]  <- .6
# # set lower bound to zero ton ensure monotonicity
# mod.pars[ mod.pars$name %in% paste0("a",1:3)     ,"lbound" ]  <- 0
# mod.pars
# # estimate model in mirt
# mod4 <- mirt::mirt(dat, mirtmodel , pars = mod.pars , verbose=TRUE , 
#             technical = list( customTheta=Theta , customPriorFun = lca_prior) )
# # correct number of estimated parameters
# mod4@nest <- as.integer(sum(mod.pars$est) + nrow(Theta)-1 ) 
# # extract coefficients
# cmod4 <- mirt.wrapper.coef(mod4)
# cmod4
# # compute item response functions
# cmod4c <- cmod4$coef[ , c("d" , "a1" , "a2" , "a3" ) ]
# probs4 <- t( apply( cmod4c , 1 , FUN = function(ll){ 
#                  plogis(cumsum(as.numeric(ll))) } ) )
# matplot( 1:4 ,  t(probs4) , type="b" , pch=1:I)
# ## End(Not run)