mirt (version 1.47)

bfactor: Full-Information Item Bi-factor and Two-Tier Analysis

Description

bfactor fits a confirmatory maximum likelihood two-tier/bifactor/testlet model to dichotomous and polytomous data under the item response theory paradigm. The IRT models are fit using a dimensional reduction EM algorithm so that regardless of the number of specific factors estimated the model only uses the number of factors in the second-tier structure plus 1. For the bifactor model the maximum number of dimensions is only 2 since the second-tier only consists of a ubiquitous unidimensional factor. See mirt for appropriate methods to be used on the objects returned from the estimation.

Usage

bfactor(
  data,
  model,
  model2 = paste0("G = 1-", ncol(data)),
  group = NULL,
  quadpts = NULL,
  invariance = "",
  ...
)

Value

function returns an object of class SingleGroupClass

(SingleGroupClass-class) or MultipleGroupClass(MultipleGroupClass-class).

Arguments

data

a matrix or data.frame that consists of numerically ordered data, organized in the form of integers, with missing data coded as NA

model

a numeric vector specifying which factor loads on which item. For example, if for a 4 item test with two specific factors, the first specific factor loads on the first two items and the second specific factor on the last two, then the vector is c(1,1,2,2). For items that should only load on the second-tier factors (have no specific component) NA values may be used as place-holders. These numbers will be translated into a format suitable for mirt.model(), combined with the definition in model2, with the letter 'S' added to the respective factor number

Alternatively, input can be specified using the mirt.model syntax with the restriction that each item must load on exactly one specific factor (or no specific factors, if it is only predicted by the general factor specified in model2)

IMPORTANT: Additional model information (e.g., keywords such as CONSTRAIN, MEAN, etc) must be specified in the model2 input instead of this specification

model2

a two-tier model specification object defined by mirt.model() or a string to be passed to mirt.model. By default the model will fit a unidimensional model in the second-tier, and therefore be equivalent to the bifactor model

group

a factor variable indicating group membership used for multiple group analyses

quadpts

number of quadrature nodes to use after accounting for the reduced number of dimensions. Scheme is the same as the one used in mirt, however it is in regards to the reduced dimensions (e.g., a bifactor model has 2 dimensions to be integrated)

invariance

see multipleGroup for details, however, the specific factor variances and means will be constrained according to the dimensional reduction algorithm

...

additional arguments to be passed to the estimation engine. See mirt for more details and examples

Author

Phil Chalmers

Details

bfactor follows the item factor analysis strategy explicated by Gibbons and Hedeker (1992), Gibbons et al. (2007), and Cai (2010). Nested models may be compared via an approximate chi-squared difference test or by a reduction in AIC or BIC (accessible via anova). See mirt for more details regarding the IRT estimation approach used in this package.

The two-tier model has a specific block diagonal covariance structure between the primary and secondary latent traits. Namely, the secondary latent traits are assumed to be orthogonal to all traits and have a fixed variance of 1, while the primary traits can be organized to vary and covary with other primary traits in the model.

$$\Sigma_{two-tier} = \left(\begin{array}{cc} G & 0 \\ 0 & diag(S) \end{array} \right)$$

The bifactor model is a special case of the two-tier model when \(G\) above is a 1x1 matrix, and therefore only 1 primary factor is being modeled. Evaluation of the numerical integrals for the two-tier model requires only \(ncol(G) + 1\) dimensions for integration since the \(S\) second order (or 'specific') factors require only 1 integration grid due to the dimension reduction technique.

Note: for multiple group two-tier analyses only the second-tier means and variances should be freed since the specific factors are not treated independently due to the dimension reduction technique.

References

Cai, L. (2010). A two-tier full-information item factor analysis model with applications. Psychometrika, 75, 581-612.

Chalmers, R. P. (2012). mirt: A Multidimensional Item Response Theory Package for the R Environment. Journal of Statistical Software, 48(6), 1-29. tools:::Rd_expr_doi("10.18637/jss.v048.i06")

Bradlow, E.T., Wainer, H., & Wang, X. (1999). A Bayesian random effects model for testlets. Psychometrika, 64, 153-168.

Gibbons, R. D., & Hedeker, D. R. (1992). Full-information Item Bi-Factor Analysis. Psychometrika, 57, 423-436.

Gibbons, R. D., Darrell, R. B., Hedeker, D., Weiss, D. J., Segawa, E., Bhaumik, D. K., Kupfer, D. J., Frank, E., Grochocinski, V. J., & Stover, A. (2007). Full-Information item bifactor analysis of graded response data. Applied Psychological Measurement, 31, 4-19.

Wainer, H., Bradlow, E.T., & Wang, X. (2007). Testlet response theory and its applications. New York, NY: Cambridge University Press.

See Also

mirt

Examples

# \donttest{

### load SAT12 and compute bifactor model with 3 specific factors
data(SAT12)
data <- key2binary(SAT12,
  key = c(1,4,5,2,3,1,2,1,3,1,2,4,2,1,5,3,4,4,1,4,3,3,4,1,3,5,1,3,1,5,4,5))
specific <- c(2,3,2,3,3,2,1,2,1,1,1,3,1,3,1,2,1,1,3,3,1,1,3,1,3,3,1,3,2,3,1,2)
mod1 <- bfactor(data, specific)
summary(mod1)
itemplot(mod1, 18, drop.zeros = TRUE) #drop the zero slopes to allow plotting

# complete factor score predictions (general + specific factors)
eaps <- fscores(mod1)

# factor score predictions for general factors only (more accurate due
# to lower dimensional integration)
eaps_gen <- fscores(mod1, method = 'EAP_general')
head(cbind(eaps_gen, NA, eaps))

# similar EAP estimates, but with respect to sum-scores
fscores(mod1, method = 'EAPsum_2.0', full.scores=FALSE)

# alternative model definition via ?mirt.model syntax
specific2 <- "S1 = 7,9,10,11,13,15,17,18,21,22,24,27,31
              S2 = 1,3,6,8,16,29,32
              S3 = 2,4,5,12,14,19,20,23,25,26,28,30"
mod2 <- bfactor(data, specific2)
anova(mod1, mod2) # same

# also equivalent using item names instead (not run)
specific3 <- "S1 = Item.7, Item.9, Item.10, Item.11, Item.13, Item.15,
                Item.17, Item.18, Item.21, Item.22, Item.24, Item.27, Item.31
              S2 = Item.1, Item.3, Item.6, Item.8, Item.16, Item.29, Item.32
              S3 = Item.2, Item.4, Item.5, Item.12, Item.14, Item.19,
                Item.20, Item.23, Item.25, Item.26, Item.28, Item.30"
# mod3 <- bfactor(data, specific3)
# anova(mod1, mod2, mod3)  # all same

### Try with fixed guessing parameters added
guess <- rep(.1,32)
mod2 <- bfactor(data, specific, guess = guess)
coef(mod2)
anova(mod1, mod2)

## don't estimate specific factor for item 32
specific[32] <- NA
mod3 <- bfactor(data, specific)
anova(mod3, mod1)

# same, but with syntax (not run)
specific3 <- "S1 = 7,9,10,11,13,15,17,18,21,22,24,27,31
              S2 = 1,3,6,8,16,29
              S3 = 2,4,5,12,14,19,20,23,25,26,28,30"
# mod3b <- bfactor(data, specific3)
# anova(mod3b)


#########
# mixed itemtype example

# simulate data
a <- matrix(c(
1,0.5,NA,
1,0.5,NA,
1,0.5,NA,
1,0.5,NA,
1,0.5,NA,
1,0.5,NA,
1,0.5,NA,
1,NA,0.5,
1,NA,0.5,
1,NA,0.5,
1,NA,0.5,
1,NA,0.5,
1,NA,0.5,
1,NA,0.5),ncol=3,byrow=TRUE)

d <- matrix(c(
-1.0,NA,NA,
-1.5,NA,NA,
 1.5,NA,NA,
 0.0,NA,NA,
2.5,1.0,-1,
3.0,2.0,-0.5,
3.0,2.0,-0.5,
3.0,2.0,-0.5,
2.5,1.0,-1,
2.0,0.0,NA,
-1.0,NA,NA,
-1.5,NA,NA,
 1.5,NA,NA,
 0.0,NA,NA),ncol=3,byrow=TRUE)
items <- rep('2PL', 14)
items[5:10] <- 'graded'

sigma <- diag(3)
dataset <- simdata(a,d,5000,itemtype=items,sigma=sigma)
itemstats(dataset)

specific <- "S1 = 1-7
             S2 = 8-14"
simmod <- bfactor(dataset, specific)
coef(simmod, simplify=TRUE)


#########
# General testlet response model (Wainer, 2007)

# simulate data
set.seed(1234)
a <- matrix(0, 12, 4)
a[,1] <- rlnorm(12, .2, .3)
ind <- 1
for(i in 1:3){
   a[ind:(ind+3),i+1] <- a[ind:(ind+3),1]
   ind <- ind+4
}
print(a)
d <- rnorm(12, 0, .5)
sigma <- diag(c(1, .5, 1, .5))
dataset <- simdata(a,d,2000,itemtype=rep('2PL', 12),sigma=sigma)
itemstats(dataset)

# estimate by applying constraints and freeing the latent variances
specific <- "S1 = 1-4
             S2 = 5-8
             S3 = 9-12"
model <- "G = 1-12
          CONSTRAIN = (1, a1, a2), (2, a1, a2), (3, a1, a2), (4, a1, a2),
            (5, a1, a3), (6, a1, a3), (7, a1, a3), (8, a1, a3),
            (9, a1, a4), (10, a1, a4), (11, a1, a4), (12, a1, a4)
          COV = S1*S1, S2*S2, S3*S3"

simmod <- bfactor(dataset, specific, model)
coef(simmod, simplify=TRUE)

# Constrained testlet model (Bradlow, 1999)
model2 <- "G = 1-12
          CONSTRAIN = (1, a1, a2), (2, a1, a2), (3, a1, a2), (4, a1, a2),
            (5, a1, a3), (6, a1, a3), (7, a1, a3), (8, a1, a3),
            (9, a1, a4), (10, a1, a4), (11, a1, a4), (12, a1, a4),
            (GROUP, COV_22, COV_33, COV_44)
          COV = S1*S1, S2*S2, S3*S3"

simmod2 <- bfactor(dataset, specific, model2)
coef(simmod2, simplify=TRUE)
anova(simmod2, simmod)


#########
# Two-tier model

# simulate data
set.seed(1234)
a <- matrix(c(
  1,0,0.5,NA,NA,
  1,0,0.5,NA,NA,
  1,0,0.5,NA,NA,
  1,0,0.5,NA,NA,
  1,0,0.5,NA,NA,
  1,0,NA,0.5,NA,
  1,0,NA,0.5,NA,
  1,0,NA,0.5,NA,
  0,1,NA,0.5,NA,
  0,1,NA,0.5,NA,
  0,1,NA,0.5,NA,
  0,1,NA,NA,0.5,
  0,1,NA,NA,0.5,
  0,1,NA,NA,0.5,
  0,1,NA,NA,0.5,
  0,1,NA,NA,0.5),ncol=5,byrow=TRUE)

d <- matrix(rnorm(16))
items <- rep('2PL', 16)

sigma <- diag(5)
sigma[1,2] <- sigma[2,1] <- .4
dataset <- simdata(a,d,2000,itemtype=items,sigma=sigma)
itemstats(dataset)

specific <- "S1 = 1-5
             S2 = 6-11
             S3 = 12-16"
model <- '
    G1 = 1-8
    G2 = 9-16
    COV = G1*G2'

# quadpts dropped for faster estimation, but not as precise
simmod <- bfactor(dataset, specific, model, quadpts = 15, TOL = 1e-3)
coef(simmod, simplify=TRUE)
summary(simmod)
itemfit(simmod, QMC=TRUE)
M2(simmod, QMC=TRUE)
residuals(simmod, QMC=TRUE)

# EAP predictions for all factors (high dimensional)
eaps_all <- fscores(simmod, QMC=TRUE, quadpts=50000)
head(eaps_all)
maps <- fscores(simmod, method = 'MAP')
head(maps)

######
# Multiple-group bifactor example from Cai, Yang, Hansen (2011)

# Table 1 info from Cai, Yang, and Hansen (2011)
intercept <- c(1, .25, -.25, -1, 1, .25, -.25, -1,1,
               .25, -.25, -1,1, .25, -.25, -1)
theta0 <- c(1, 1.4,1.7,2,
            1.4,1.7,2,1,
            1.7,2,1,1.4,
            2, 1, 1.4, 1.7)
thetan <- c(.8,1.5,1.2,1,
            1,.8,1.5,1.2,
            1.2,1,.8,1.5,
            1.5,1.2,1,.8)
thetaN <- matrix(0, 16, 4)
thetaN[1:4, 1] <- thetan[1:4]
thetaN[1:4+4, 2] <- thetan[1:4+4]
thetaN[1:4+8, 3] <- thetan[1:4+8]
thetaN[1:4+12, 4] <- thetan[1:4+12]

as <- cbind(theta0, thetaN)
as

# data generation for focal group in publication does not have response
# for items 13-16. However, full-data approach presented first
N <- 1000
itemtype <- '2PL'
gmeans <- c(1, -.5, 0, .5, 0)
sigma <- diag(c(.8, 1.2, 1.5, 1, 1))

datG1 <- simdata(as, intercept, N=N, itemtype='2PL')
datG2 <- simdata(as, intercept, N=N, itemtype='2PL',
  mu = gmeans, sigma = sigma)

dat <- rbind(datG1, datG2)
group <- rep(c('G1', 'G2'), each=N)

specific <- "S1 = 1-4
             S2 = 5-8
             S3 = 9-12
             S4 = 13-16"

mod <- bfactor(dat, specific, group=group, SE=TRUE,
  invariance=c('free_means', 'free_vars', colnames(dat)))
coef(mod, simplify=TRUE)

# DIF testing for general dimension only with likelihood ratio tests (not run)
if(FALSE){
  mirtCluster()  # allocate computing cores
  DIF(mod, which.par=c('a1', 'd'), scheme='drop')
}

## same analysis, however items 13:16 do not exist in the focal group
datG2[,13:16] <- NA    # use this to match publication
dat <- rbind(datG1, datG2)
head(dat)
tail(dat)

# specify mean/cov structure explicitly
model2 <- "G = 1-16
           MEAN [G2] = G, S1, S2, S3
           COV [G2] = G*G, S1*S1, S2*S2, S3*S3"

mod2 <- bfactor(dat, specific, model2, group=group, invariance=colnames(dat))
coef(mod2, simplify=TRUE)

# }