Multiple groups with equal and unequal item parameters - tmatta/lsasim GitHub Wiki
The following example demonstrates how to generate item responses from multiple populations. We will generate 15 items from a Rasch model across 3 test booklets. The three booklets will then be administered to 5000 examinees from two different populations. Both means and variances will differ between the two populations. We will use TAM (Robitzsch, Kiefer, & Wu, 2017) to estimate item parameters from the pooled two-country sample. We will then estimate theta distributions through weighted likelihood estimator (Warm, 1989) using fixed item parameters.
We will provide two examples in this vignette. In the first example, item parameters are equal across two groups. In the second example, item parameters of a few items are not equal across two groups.
To begin the example, we load both TAM and lsasim as well as set the seed for reproducibility of results.
library(TAM)
library(lsasim)
set.seed(1407)
Generating multiple groups
As mentioned above, we will be generating item responses for two different populations. In an international large-scale assessment setting like PISA, this would be akin to generating item responses for two countries.
We generate 5000 (\theta)s from each of the two country's population parameters separately.
N <- 5000
theta_cnt1 <- data.frame(theta = rnorm(N, 0, 1))
theta_cnt2 <- data.frame(theta = rnorm(N, 1, 0.75))
Test booklet generation
Next, we generate the 15-item item pool using item_gen and distribute them to three item blocks, which are combined into three test booklets using the block_design and booklet_design functions, respectively. The item_gen function was described here and the block_design and booklet_design functions were described here.
I <- 15
K <- 3
Example 1: item parameters are equal across groups
We generate item parameters from a Rasch model. Items are same across two groups.
item_pool <- lsasim::item_gen(n_1pl = I, b_bounds = c(-2, 2))
print(item_pool)
## item b a c k p
## 1 1 -0.04 1 0 1 1
## 2 2 -1.19 1 0 1 1
## 3 3 0.94 1 0 1 1
## 4 4 1.49 1 0 1 1
## 5 5 -1.53 1 0 1 1
## 6 6 1.34 1 0 1 1
## 7 7 -0.78 1 0 1 1
## 8 8 0.97 1 0 1 1
## 9 9 1.35 1 0 1 1
## 10 10 -1.10 1 0 1 1
## 11 11 -1.59 1 0 1 1
## 12 12 -1.24 1 0 1 1
## 13 13 -0.27 1 0 1 1
## 14 14 0.49 1 0 1 1
## 15 15 1.15 1 0 1 1
Next, we distribute the (I = 15) items to (K = 3) blocks, and those (K) blocks to 3 test booklets.
blocks <- lsasim::block_design(n_blocks = K, item_parameters = item_pool)
print(blocks)
## $block_assignment
## b1 b2 b3
## i1 1 2 3
## i2 4 5 6
## i3 7 8 9
## i4 10 11 12
## i5 13 14 15
##
## $block_descriptives
## block length average difficulty
## b1 5 -0.140
## b2 5 -0.570
## b3 5 0.708
books <- lsasim::booklet_design(item_block_assignment = blocks$block_assignment)
print(books)
## B1 B2 B3
## i1 1 2 1
## i2 4 5 4
## i3 7 8 7
## i4 10 11 10
## i5 13 14 13
## i6 2 3 3
## i7 5 6 6
## i8 8 9 9
## i9 11 12 12
## i10 14 15 15
Finally, we administer the test booklets to the 10000 examinees across the two countries using the booklet_sample function. The booklet_sample function was described here.
book_samp <- lsasim::booklet_sample(n_subj = N*2,
book_item_design = books,
book_prob = NULL)
Generating item responses
Because we combined the two countries and item are the across two groups, generating item responses is straight forward.
cog_ex1 <- lsasim::response_gen(subject = book_samp$subject,
item = book_samp$item,
theta = rbind(theta_cnt1, theta_cnt2)$theta,
b_par = item_pool$b)
However, because we will need to know which examinees belong to which country when we estimate the population parameters, we add an additional variable, country to the cog_ex1 data frame. The first 5000 examinees are from country 1 and are thus coded country = 1 and the next 5000 examinees are from country 2 and are coded country = 2.
cog_ex1$country <- rep(1:2, each = N)
Estimation
To estimate separate theta distributions for two groups', we use the TAM function tam.mml and specify the country indicators, country for the group argument. By default, TAM fixes the mean for the first group to zero.
mod1 <- TAM::tam.mml(cog_ex1[,1:I], group = cog_ex1$country)
mod1$item gives us a data frame with item parameters. The column xsi.item denotes the estimated item difficulty.
mod1$item
## item N M xsi.item AXsi_.Cat1 B.Cat1.Dim1
## 1 i001 6545 0.6163484 -0.06702407 -0.06702407 1
## 2 i002 6746 0.8109991 -1.26437672 -1.26437672 1
## 3 i003 6709 0.4237591 0.88282556 0.88282556 1
## 4 i004 6545 0.3016043 1.51550511 1.51550511 1
## 5 i005 6746 0.8482063 -1.57249052 -1.57249052 1
## 6 i006 6709 0.3455060 1.27848527 1.27848527 1
## 7 i007 6545 0.7477464 -0.80888605 -0.80888605 1
## 8 i008 6746 0.4067596 0.94442867 0.94442867 1
## 9 i009 6709 0.3431212 1.29099706 1.29099706 1
## 10 i010 6545 0.7964859 -1.14028585 -1.14028585 1
## 11 i011 6746 0.8545805 -1.63075718 -1.63075718 1
## 12 i012 6709 0.8102549 -1.23673837 -1.23673837 1
## 13 i013 6545 0.6559206 -0.27525805 -0.27525805 1
## 14 i014 6746 0.5059294 0.46059579 0.46059579 1
## 15 i015 6709 0.3821732 1.08998169 1.08998169 1
As the results from mod1 show, the estimated mean scores for the two countries are 0 and 0.968 respectively. Recall that the mean of country 1 was fixed to 0. The estimated standard deviations for the two countries are 1.016 and 0.737. Comparing these estimates with their generating parameters: ((\theta_{1}) ~ N(0, 1); (\theta_{2}) ~ N(1, .75)) suggests that they are recovered well.
mod1_results <- rbind.data.frame(as.numeric(mod1$beta),
cbind(C1 = sqrt(mean(mod1$variance[cog_ex1$country == 1])),
C2 = sqrt(mean(mod1$variance[cog_ex1$country == 2]))))
rownames(mod1_results) <- c("Mean", "SD")
round(mod1_results, 3)
## C1 C2
## Mean 0.000 0.968
## SD 1.016 0.737
We restimate the Rasch model with fixed item difficulties using our generating item parameters, item_pool$b.
xsi.fixed <- cbind(1:nrow(mod1$item), item_pool$b)
print(xsi.fixed)
## [,1] [,2]
## [1,] 1 -0.04
## [2,] 2 -1.19
## [3,] 3 0.94
## [4,] 4 1.49
## [5,] 5 -1.53
## [6,] 6 1.34
## [7,] 7 -0.78
## [8,] 8 0.97
## [9,] 9 1.35
## [10,] 10 -1.10
## [11,] 11 -1.59
## [12,] 12 -1.24
## [13,] 13 -0.27
## [14,] 14 0.49
## [15,] 15 1.15
mod2 <- TAM::tam.mml(cog_ex1[,1:I], group = cog_ex1$country, xsi.fixed=xsi.fixed)
As the results from mod2 show, the estimated mean scores for the two countries are 0.034 and 1.003 respectively. The estimated standard deviations for the two countries are 1.017 and 0.737. Comparing these estimates with their generating parameters: ((\theta_{1}) ~ N(0, 1); (\theta_{2}) ~ N(1, .75)) suggests that they are recovered as well as mod1.
mod2_results <- rbind.data.frame(as.numeric(mod2$beta),
cbind(C1 = sqrt(mean(mod2$variance[cog_ex1$country == 1])),
C2 = sqrt(mean(mod2$variance[cog_ex1$country == 2]))))
rownames(mod2_results) <- c("Mean", "SD")
round(mod2_results, 3)
## C1 C2
## Mean 0.034 1.003
## SD 1.017 0.737
Example 2: item parameters are not equal across groups
We use the same item parameters generated from Example 1.
item_pool_cnt1 <- item_pool
print(item_pool_cnt1[1:5,])
## item b a c k p
## 1 1 -0.04 1 0 1 1
## 2 2 -1.19 1 0 1 1
## 3 3 0.94 1 0 1 1
## 4 4 1.49 1 0 1 1
## 5 5 -1.53 1 0 1 1
The first five items for country 2 are different from country 1. The first five item difficulty parameters for country 2 are -0.54, -1.69, 0.44, 0.99, and -2.03, instead of -0.04, -1.19, 0.94, 1.49, and -1.53. We generate item responses for each country separately.
item_pool_cnt2 <- item_pool
item_pool_cnt2$b[1:5] <- item_pool_cnt2$b[1:5]-.5
print(item_pool_cnt2[1:5,])
## item b a c k p
## 1 1 -0.54 1 0 1 1
## 2 2 -1.69 1 0 1 1
## 3 3 0.44 1 0 1 1
## 4 4 0.99 1 0 1 1
## 5 5 -2.03 1 0 1 1
For country 1.
blocks_cnt1 <- lsasim::block_design(n_blocks = K,
item_parameters = item_pool_cnt1,
item_block_matrix = NULL)
books_cnt1 <- lsasim::booklet_design(item_block_assignment = blocks_cnt1$block_assignment,
book_design = NULL)
book_samp_cnt1 <- lsasim::booklet_sample(n_subj = N,
book_item_design = books_cnt1,
book_prob = NULL)
cog_cnt1 <- lsasim::response_gen(subject = book_samp_cnt1$subject,
item = book_samp_cnt1$item,
theta = theta_cnt1$theta,
b_par = item_pool_cnt1$b)
And for country 2.
blocks_cnt2 <- lsasim::block_design(n_blocks = K,
item_parameters = item_pool_cnt2,
item_block_matrix = NULL)
books_cnt2 <- lsasim::booklet_design(item_block_assignment = blocks_cnt2$block_assignment,
book_design = NULL)
book_samp_cnt2 <- lsasim::booklet_sample(n_subj = N,
book_item_design = books_cnt2,
book_prob = NULL)
cog_cnt2 <- lsasim::response_gen(subject = book_samp_cnt2$subject,
item = book_samp_cnt2$item,
theta = theta_cnt2$theta,
b_par = item_pool_cnt2$b)
After item responses for each country were generated, we added an additional variable, country and combined the item responses from two countries.
cog_cnt1$country <- rep(1, each = N)
cog_cnt2$country <- rep(2, each = N)
cog_ex2 <- rbind(cog_cnt1, cog_cnt2)
Estimation
We followed the same calibration procedure described above and used tam.mml to esimtate separate theta distributions for two groups'. We specify the country indicators, country for the group argument.
mod3 <- TAM::tam.mml(cog_ex2[,1:I], group = cog_ex2$country)
As the results from mod3 show, the estimated mean scores for the two countries are 0 and 1.121 respectively. Recall that the mean of country 1 was fixed to 0. The estimated standard deviations for the two countries are 1.004 and 0.689. Comparing these estimates with their generating parameters: ((\theta_{1}) ~ N(0, 1); (\theta_{2}) ~ N(1, .75)) suggests that they are not recovered well as mod1 and mod2. When item parameters are not equal across two groups, model calibration for each group is required to avoid biased estimates.
mod3_results <- rbind.data.frame(as.numeric(mod3$beta),
cbind(C1 = sqrt(mean(mod3$variance[cog_ex2$country == 1])),
C2 = sqrt(mean(mod3$variance[cog_ex2$country == 2]))))
rownames(mod3_results) <- c("Mean", "SD")
round(mod3_results, 3)
## C1 C2
## Mean 0.000 1.121
## SD 1.004 0.689