############################################################
### EXAMPLE 1: PISA Reading
############################################################
data( data.pisaRead)
dat <- data.pisaRead$data
items <- grep("R" , colnames(dat))
# define matrix of covariates
X <- cbind( 1 , dat[ , c("female","hisei","migra" ) ] )
#***
# Model 1: Latent regression model in the Rasch model
# estimate Rasch model
mod1 <- rasch.mml2( dat[,items] )
# latent regression model
lm1 <- latent.regression.em.raschtype( data=dat[,items ],
X = X , b = mod1$item$b )
#***
# Model 2: Latent regression with generalized link function
# estimate alpha parameters for link function
mod2 <- rasch.mml2( dat[,items] , est.alpha=TRUE)
# use model estimated likelihood for latent regression model
lm2 <- latent.regression.em.raschtype( f.yi.qk=mod2$f.yi.qk ,
X = X , theta.list=mod2$theta.k)
#***
# Model 3: Latent regression model based on Rasch copula model
testlets <- paste( data.pisaRead$item$testlet)
itemclusters <- match( testlets , unique(testlets) )
# estimate Rasch copula model
mod3 <- rasch.copula2( dat[,items] , itemcluster=itemclusters )
# use model estimated likelihood for latent regression model
lm3 <- latent.regression.em.raschtype( f.yi.qk=mod3$f.yi.qk ,
X = X , theta.list=mod3$theta.k)
############################################################
### SIMULATED EXAMPLE 2
set.seed(899)
############################################################
I <- 21 # number of items
b <- seq(-2,2, len=I) # item difficulties
n <- 2000 # number of students
# simulate theta and covariates
theta <- rnorm( n )
x <- .7 * theta + rnorm( n , .5 )
y <- .2 * x+ .3*theta + rnorm( n , .4 )
dfr <- data.frame( theta , 1 , x , y )
# simulate Rasch model
dat1 <- sim.raschtype( theta = theta , b = b )
# estimate latent regression
mod <- latent.regression.em.raschtype( data = dat1 , X = dfr[,-1] , b = b )
## Regression Parameters
##
## est se.simple se t p beta fmi N.simple pseudoN.latent
## X1 -0.2554 0.0208 0.0248 -10.2853 0 0.0000 0.2972 2000 1411.322
## x 0.4113 0.0161 0.0193 21.3037 0 0.4956 0.3052 2000 1411.322
## y 0.1715 0.0179 0.0213 8.0438 0 0.1860 0.2972 2000 1411.322
##
## Residual Variance = 0.685
## Explained Variance = 0.3639
## Total Variance = 1.049
## R2 = 0.3469
# compare with linear model (based on true scores)
summary( lm( theta ~ x + y , data = dfr ) )
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.27821 0.01984 -14.02 <2e-16 ***
## x 0.40747 0.01534 26.56 <2e-16 ***
## y 0.18189 0.01704 10.67 <2e-16 ***
## ---
##
## Residual standard error: 0.789 on 1997 degrees of freedom
## Multiple R-squared: 0.3713, Adjusted R-squared: 0.3707
#***********
# define guessing parameters (lower asymptotes) and
# upper asymptotes ( 1 minus slipping parameters)
cI <- rep(.2, I) # all items get a guessing parameter of .2
cI[ c(7,9) ] <- .25 # 7th and 9th get a guessing parameter of .25
dI <- rep( .95 , I ) # upper asymptote of .95
dI[ c(7,11) ] <- 1 # 7th and 9th item have an asymptote of 1
# latent regression model
mod1 <- latent.regression.em.raschtype( data = dat1 , X = dfr[,-1] ,
b = b , c = cI , d = dI )
## Regression Parameters
##
## est se.simple se t p beta fmi N.simple pseudoN.latent
## X1 -0.7929 0.0243 0.0315 -25.1818 0 0.0000 0.4044 2000 1247.306
## x 0.5025 0.0188 0.0241 20.8273 0 0.5093 0.3936 2000 1247.306
## y 0.2149 0.0209 0.0266 8.0850 0 0.1960 0.3831 2000 1247.306
##
## Residual Variance = 0.9338
## Explained Variance = 0.5487
## Total Variance = 1.4825
## R2 = 0.3701
Run the code above in your browser using DataLab