#############################################################################
# EXAMPLE 1: Reading data set
#############################################################################
data(data.read)
dat <- data.read
#******
# Model 1: Rasch model with PML estimation
mod1 <- rasch.pml3( dat )
summary(mod1)
#******
# Model 2: Excluding item pairs with local dependence
# from bivariate composite likelihood
itemcluster <- rep( 1:3 , each=4)
mod2 <- rasch.pml3( dat , itemcluster = itemcluster )
summary(mod2)
## Not run:
# #*****
# # Model 3: Modelling error correlations:
# # joint residual correlations for each itemcluster
# error.corr <- diag(1,ncol(dat))
# for ( ii in 1:3){
# ind.ii <- which( itemcluster == ii )
# error.corr[ ind.ii , ind.ii ] <- ii
# }
# # estimate the model with error correlations
# mod3 <- rasch.pml3( dat , error.corr = error.corr )
# summary(mod3)
#
# #****
# # Model 4: model separate residual correlations
# I <- ncol(error.corr)
# error.corr1 <- matrix( 1:(I*I) , ncol= I )
# error.corr <- error.corr1 * ( error.corr > 0 )
# # estimate the model with error correlations
# mod4 <- rasch.pml3( dat , error.corr = error.corr )
# summary(mod4)
#
# #****
# # Model 5: assume equal item difficulties:
# # b_1 = b_7 and b_2 = b_12
# # fix item difficulty of the 6th item to .1
# est.b <- 1:I
# est.b[7] <- 1; est.b[12] <- 2 ; est.b[6] <- 0
# b.init <- rep( 0, I ) ; b.init[6] <- .1
# mod5 <- rasch.pml3( dat , est.b =est.b , b.init=b.init)
# summary(mod5)
#
# #****
# # Model 6: estimate three item slope groups
# est.a <- rep(1:3 , each=4 )
# mod6 <- rasch.pml3( dat , est.a =est.a , est.sigma=0)
# summary(mod6)
#
# #############################################################################
# # EXAMPLE 2: PISA reading
# #############################################################################
#
# data(data.pisaRead)
# dat <- data.pisaRead$data
#
# # select items
# dat <- dat[ , substring(colnames(dat),1,1)=="R" ]
#
# #******
# # Model 1: Rasch model with PML estimation
# mod1 <- rasch.pml3( as.matrix(dat) )
# ## Trait SD (Logit Link) : 1.419
#
# #******
# # Model 2: Model correlations within testlets
# error.corr <- diag(1,ncol(dat))
# testlets <- paste( data.pisaRead$item$testlet )
# itemcluster <- match( testlets , unique(testlets ) )
# for ( ii in 1:(length(unique(testlets))) ){
# ind.ii <- which( itemcluster == ii )
# error.corr[ ind.ii , ind.ii ] <- ii
# }
# # estimate the model with error correlations
# mod2 <- rasch.pml3( dat , error.corr = error.corr )
# ## Trait SD (Logit Link) : 1.384
#
# #****
# # Model 3: model separate residual correlations
# I <- ncol(error.corr)
# error.corr1 <- matrix( 1:(I*I) , ncol= I )
# error.corr <- error.corr1 * ( error.corr > 0 )
# # estimate the model with error correlations
# mod3 <- rasch.pml3( dat , error.corr = error.corr )
# ## Trait SD (Logit Link) : 1.384
#
# #############################################################################
# # EXAMPLE 3: 10 locally independent items
# #############################################################################
#
# #**********
# # simulate some data
# set.seed(554)
# N <- 500 # persons
# I <- 10 # items
# theta <- stats::rnorm(N,sd=1.3 ) # trait SD of 1.3
# b <- seq(-2 , 2 , length=I) # item difficulties
#
# # simulate data from the Rasch model
# dat <- sim.raschtype( theta = theta , b = b )
#
# # estimation with rasch.pml and probit link
# mod1 <- rasch.pml3( dat )
# summary(mod1)
#
# # estimation with rasch.mml2 function
# mod2 <- rasch.mml2( dat )
#
# # estimate item parameters for groups with five item parameters each
# est.b <- rep( 1:(I/2) , each=2 )
# mod3 <- rasch.pml3( dat , est.b=est.b )
# summary(mod3)
#
# # compare parameter estimates
# summary(mod1)
# summary(mod2)
# summary(mod3)
#
# #############################################################################
# # EXAMPLE 4: 11 items and 2 item clusters with 2 and 3 items
# #############################################################################
#
# set.seed(5698)
# I <- 11 # number of items
# n <- 5000 # number of persons
# b <- seq(-2,2, len=I) # item difficulties
# theta <- stats::rnorm( n , sd = 1 ) # person abilities
# # itemcluster
# itemcluster <- rep(0,I)
# itemcluster[c(3,5)] <- 1
# itemcluster[c(2,4,9)] <- 2
# # residual correlations
# rho <- c( .7 , .5 )
#
# # simulate data (under the logit link)
# dat <- sim.rasch.dep( theta , b , itemcluster , rho )
# colnames(dat) <- paste("I" , seq(1,ncol(dat)) , sep="")
#
# #***
# # Model 1: estimation using the Rasch model (with probit link)
# mod1 <- rasch.pml3( dat )
# #***
# # Model 2: estimation when pairs of locally dependent items are eliminated
# mod2 <- rasch.pml3( dat , itemcluster=itemcluster)
#
# #***
# # Model 3: Positive correlations within testlets
# est.corrs <- diag( 1 , I )
# est.corrs[ c(3,5) , c(3,5) ] <- 2
# est.corrs[ c(2,4,9) , c(2,4,9) ] <- 3
# mod3 <- rasch.pml3( dat , error.corr=est.corrs )
#
# #***
# # Model 4: Negative correlations between testlets
# est.corrs <- diag( 1 , I )
# est.corrs[ c(3,5) , c(2,4,9) ] <- 2
# est.corrs[ c(2,4,9) , c(3,5) ] <- 2
# mod4 <- rasch.pml3( dat , error.corr=est.corrs )
#
# #***
# # Model 5: sum constraint of zero within and between testlets
# est.corrs <- matrix( 1:(I*I) , I , I )
# cluster2 <- c(2,4,9)
# est.corrs[ setdiff( 1:I , c(cluster2)) , ] <- 0
# est.corrs[ , setdiff( 1:I , c(cluster2)) ] <- 0
# # define an error constraint matrix
# itempairs0 <- mod4$itempairs
# IP <- nrow(itempairs0)
# err.constraint <- matrix( 0 , IP , 1 )
# err.constraint[ ( itempairs0$item1 %in% cluster2 )
# & ( itempairs0$item2 %in% cluster2 ) , 1 ] <- 1
# # set sum of error covariances to 1.2
# err.constraintV <- matrix(3*.4,1,1)
#
# mod5 <- rasch.pml3( dat , error.corr=est.corrs ,
# err.constraintM=err.constraint, err.constraintV=err.constraintV)
#
# #****
# # Model 6: Constraint on sum of all correlations
# est.corrs <- matrix( 1:(I*I) , I , I )
# # define an error constraint matrix
# itempairs0 <- mod4$itempairs
# IP <- nrow(itempairs0)
# # define two side conditions
# err.constraint <- matrix( 0 , IP , 2 )
# err.constraintV <- matrix( 0 , 2 , 1)
# # sum of all correlations is zero
# err.constraint[ , 1 ] <- 1
# err.constraintV[1,1] <- 0
# # sum of items cluster c(1,2,3) is 0
# cluster2 <- c(1,2,3)
# err.constraint[ ( itempairs0$item1 %in% cluster2 )
# & ( itempairs0$item2 %in% cluster2 ) , 2 ] <- 1
# err.constraintV[2,1] <- 0
#
# mod6 <- rasch.pml3( dat , error.corr=est.corrs ,
# err.constraintM=err.constraint, err.constraintV=err.constraintV)
# summary(mod6)
#
# #############################################################################
# # EXAMPLE 5: 10 Items: Cluster 1 -> Items 1,2
# # Cluster 2 -> Items 3,4,5; Cluster 3 -> Items 7,8,9
# #############################################################################
#
# set.seed(7650)
# I <- 10 # number of items
# n <- 5000 # number of persons
# b <- seq(-2,2, len=I) # item difficulties
# bsamp <- b <- sample(b) # sample item difficulties
# theta <- stats::rnorm( n , sd = 1 ) # person abilities
# # define itemcluster
# itemcluster <- rep(0,I)
# itemcluster[ 1:2 ] <- 1
# itemcluster[ 3:5 ] <- 2
# itemcluster[ 7:9 ] <- 3
# # define residual correlations
# rho <- c( .55 , .35 , .45)
#
# # simulate data
# dat <- sim.rasch.dep( theta , b , itemcluster , rho )
# colnames(dat) <- paste("I" , seq(1,ncol(dat)) , sep="")
#
# #***
# # Model 1: residual correlation (equal within item clusters)
# # define a matrix of integers for estimating error correlations
# error.corr <- diag(1,ncol(dat))
# for ( ii in 1:3){
# ind.ii <- which( itemcluster == ii )
# error.corr[ ind.ii , ind.ii ] <- ii
# }
# # estimate the model
# mod1 <- rasch.pml3( dat , error.corr = error.corr )
#
# #***
# # Model 2: residual correlation (different within item clusters)
# # define again a matrix of integers for estimating error correlations
# error.corr <- diag(1,ncol(dat))
# for ( ii in 1:3){
# ind.ii <- which( itemcluster == ii )
# error.corr[ ind.ii , ind.ii ] <- ii
# }
# I <- ncol(error.corr)
# error.corr1 <- matrix( 1:(I*I) , ncol= I )
# error.corr <- error.corr1 * ( error.corr > 0 )
# # estimate the model
# mod2 <- rasch.pml3( dat , error.corr = error.corr )
#
# #***
# # Model 3: eliminate item pairs within itemclusters for PML estimation
# mod3 <- rasch.pml3( dat , itemcluster = itemcluster )
#
# #***
# # Model 4: Rasch model ignoring dependency
# mod4 <- rasch.pml3( dat )
#
# # compare different models
# summary(mod1)
# summary(mod2)
# summary(mod3)
# summary(mod4)
# ## End(Not run)
Run the code above in your browser using DataLab