Estimating the Model in the Paper

2022-02-21

Intro

The objective of this vignette is to provide clarity as to the estimation procedure used to obtain the results in the paper. This vignette breaks down each section in the example displayed in cIRT() function.

Data

Within this vignette, we used two different data sets to generate a Choice Item Response Theory Model with the routine located in cIRT.

library(cIRT)

The first data set, trial_matrix, contains whether or not the subject correctly identified the spatial rotation. The second dataset, choice_matrix, provides information regarding the choice decision subjects were asked to make.

data(trial_matrix)
data(choice_matrix)

Constructing the Model Matrix

Here we construct a thurstone design matrix by obtaining the IDs of hard and easy questions presented for the subject to make a decision upon.

# Create the Thurstone Design Matrices
hard_items = choice_matrix$hard_q_id
easy_items = choice_matrix$easy_q_id

D_easy = model.matrix( ~ -1 + factor(easy_items))
D_hard = -1 * model.matrix( ~ -1 + factor(hard_items))[, -c(5, 10, 15)]

Within this setting, we setup the effect-codes for different constraints.

# Defining effect-coded contrasts
high_contrasts = rbind(-1, diag(4))
rownames(high_contrasts) = 12:16
low_contrasts = rbind(-1, diag(2))
rownames(low_contrasts) = 4:6

# Creating high & low factors
high = factor(choice_matrix[, 'high_value'])
low = factor(choice_matrix[, 'low_value'])
contrasts(high) = high_contrasts
contrasts(low) = low_contrasts

fixed_effects = model.matrix( ~ high + low)
fixed_effects_base = fixed_effects[, 1]
fixed_effects_int = model.matrix( ~ high * low)

Modeling the Data

Generate the cIRT model using a Thurstone Design Matrix generated above.

# Model with Thurstone D matrix
system.time({
  out_model_thurstone = cIRT(
    choice_matrix[, 'subject_id'],
    cbind(fixed_effects[, -1], D_easy, D_hard),
    c(1:ncol(fixed_effects)),
    as.matrix(fixed_effects),
    as.matrix(trial_matrix),
    choice_matrix[, 'choose_hard_q'],
    20000,
    25000
  )
})
##    user  system elapsed 
## 327.137   6.405 359.415

We recommend saving the model object as a .rda file even though the total computational time is less than 2.5 minutes.

## Save model output to an rda file.
# save(out_model_thurstone, file='choiceMCMCoutput.rda')

## Load model output back into R.
# load(file='choiceMCMCoutput.rda')

Parameter Estimates

Next up, we obtain the parameter estimates of the model by averaging over the different estimates obtained via the Gibbs sampling technique employed.

vlabels_thurstone = colnames(cbind(fixed_effects[, -1], D_easy, D_hard))

G_thurstone = t(apply(
  out_model_thurstone$gs0,
  2,
  FUN = quantile,
  probs = c(.5, .025, .975)
))
rownames(G_thurstone) = vlabels_thurstone

B_thurstone = t(apply(
  out_model_thurstone$beta,
  2,
  FUN = quantile,
  probs = c(.5, 0.025, .975)
))
rownames(B_thurstone) = colnames(fixed_effects)

S_thurstone = solve(
  apply(out_model_thurstone$Sigma_zeta_inv, c(1, 2), FUN = mean)
)

inv_sd = diag(1 / sqrt(diag(solve(
  apply(out_model_thurstone$Sigma_zeta_inv, c(1, 2), FUN = mean)
))))

corrmat = inv_sd %*% S_thurstone %*% inv_sd
as = apply(out_model_thurstone$as, 2, FUN = mean)
bs = apply(out_model_thurstone$bs, 2, FUN = mean)

Thus, we have the following results:

# gs0
G_thurstone
##                               50%         2.5%        97.5%
## high1                -0.011986505 -0.118096320  0.095107122
## high2                -0.079234376 -0.183565722  0.025500545
## high3                 0.129954364  0.023954544  0.234103673
## high4                 0.154603793  0.044170320  0.265609784
## low1                  0.003644706 -0.075253794  0.082219112
## low2                  0.086774878  0.006642201  0.165868729
## factor(easy_items)1  -0.255902418 -0.532232477  0.029533337
## factor(easy_items)2  -0.360954584 -0.642187540 -0.073575511
## factor(easy_items)3  -0.242949397 -0.522774877  0.042362258
## factor(easy_items)4  -0.398500907 -0.674870956 -0.124139855
## factor(easy_items)5  -0.197921825 -0.484754367  0.088724567
## factor(easy_items)6  -0.443660385 -0.727736648 -0.166942035
## factor(easy_items)7  -0.295547931 -0.580973028 -0.008978882
## factor(easy_items)8   0.056440781 -0.215232580  0.329255011
## factor(easy_items)9  -0.370512198 -0.652698716 -0.091148667
## factor(easy_items)10 -0.034019803 -0.314521261  0.242682479
## factor(easy_items)11 -0.297678473 -0.575926391 -0.020115218
## factor(easy_items)12  0.035772787 -0.237049486  0.309729138
## factor(easy_items)13 -0.262907572 -0.550424105  0.021391375
## factor(easy_items)14 -0.512496412 -0.787338382 -0.234069012
## factor(easy_items)15 -0.092065468 -0.370120950  0.182016774
## factor(hard_items)16 -0.518383831 -0.784117644 -0.252111125
## factor(hard_items)17 -0.569133477 -0.835826309 -0.304754958
## factor(hard_items)18 -1.391503528 -1.681528441 -1.106541498
## factor(hard_items)19 -0.206118019 -0.472978241  0.058422458
## factor(hard_items)21 -0.716516073 -0.991132072 -0.444580149
## factor(hard_items)22 -0.249177487 -0.519388628  0.015502068
## factor(hard_items)23 -0.221818958 -0.487172904  0.042926480
## factor(hard_items)24 -1.205749013 -1.494873751 -0.924907859
## factor(hard_items)26 -0.155672145 -0.423201723  0.110397390
## factor(hard_items)27 -0.256375185 -0.520746891  0.003902338
## factor(hard_items)28  0.368892433  0.097393703  0.637294233
## factor(hard_items)29  0.752571195  0.474303203  1.035419614
# betas
B_thurstone
##                      50%        2.5%       97.5%
## (Intercept)  0.695148869  0.56470595  0.83213822
## high1       -0.107902888 -0.23864627  0.02227350
## high2       -0.008217892 -0.13901126  0.12997900
## high3        0.036147793 -0.09418339  0.17644129
## high4        0.102310293 -0.03642197  0.24900377
## low1         0.023120289 -0.07095377  0.12498188
## low2        -0.116479876 -0.21337626 -0.02293902
# Sigma Thurstone
S_thurstone
##              [,1]         [,2]         [,3]         [,4]         [,5]
## [1,]  0.384416241 -0.007365872  0.010874414  0.003774921 -0.015137776
## [2,] -0.007365872  0.077353168 -0.011380124 -0.002565672 -0.024912286
## [3,]  0.010874414 -0.011380124  0.063993915 -0.007866015 -0.006398075
## [4,]  0.003774921 -0.002565672 -0.007866015  0.062647588 -0.008309805
## [5,] -0.015137776 -0.024912286 -0.006398075 -0.008309805  0.084759017
## [6,] -0.012319126  0.005282207 -0.002257572 -0.005987035  0.002229131
## [7,] -0.016812867 -0.006756094  0.001544687  0.003783982  0.006225211
##              [,6]         [,7]
## [1,] -0.012319126 -0.016812867
## [2,]  0.005282207 -0.006756094
## [3,] -0.002257572  0.001544687
## [4,] -0.005987035  0.003783982
## [5,]  0.002229131  0.006225211
## [6,]  0.064425714 -0.013224306
## [7,] -0.013224306  0.062298035
## Item parameters ----

# a
as
##  [1] 0.6294940 0.5004532 1.1337972 0.3215831 1.0083180 1.0777463 0.8164757
##  [8] 0.3753405 0.8190618 0.7419126 0.6163115 1.2075325 0.8378053 0.9734849
## [15] 0.7943470 1.0149385 0.7102592 0.5866380 0.9448062 1.0784960 1.0144621
## [22] 0.6586769 0.6333581 0.5813622 0.6158422 0.6159609 0.6653505 0.7934624
## [29] 0.5041641 0.6428622
# b
bs
##  [1] -1.13721344 -1.15805252 -1.53520016 -0.88584178 -1.01008772 -0.93959630
##  [7] -0.84002304 -0.43503494 -0.47496244 -0.74840458 -0.69300641 -0.59557173
## [13] -0.66524799 -0.48926592 -0.30385062 -0.59018710 -0.49888004 -0.26667694
## [19] -0.43041537 -0.52830544 -0.32355055 -0.30627683 -0.40235259 -0.22419274
## [25] -0.06060915  0.06529378 -0.31580061  0.13138751  0.48997732  0.93040382