Main simulation
The main simulation uses the PROsetta
package.
The simulation consists of 8 conditions * 20 trials = 160 tasks.
# Code: Sangdon Lim
library(PROsetta)
library(mvnfast)
library(mirt)
library(stringr)
library(parallel)
library(foreach)
library(doParallel)
<- detectCores() - 2
n_cores <- makeCluster(n_cores)
cl registerDoParallel(cl)
source("sim_get_ipar.R")
source("sim_get_data.R")
source("sim_eval_table.R")
<- loadData(
d response = "dat_dePHQ9_validation.csv",
itemmap = "imap_dePHQ9.csv",
anchor = "anchor_dePHQ9.csv",
input_dir = "data"
)
<- 2
lhs_scale <- 1
rhs_scale
<- seq(1.0, 0.6, -.05)
c1 <- 1:20
c2 <- expand.grid(c1, c2)
tasks <- nrow(tasks)
n_tasks
<- foreach(
out idx_task = 1:n_tasks,
.packages = c("mvnfast", "PROsetta", "mirt", "stringr")) %dopar% {
# Fetch task information
<- tasks[idx_task, 1]
true_corr <- tasks[idx_task, 2]
idx_trial
# Generate thetas and response data
set.seed(idx_trial)
<- list()
mu_sigma $mu <- rep(0, 2)
mu_sigma$sigma <- matrix(true_corr, 2, 2)
mu_sigmadiag(mu_sigma$sigma) <- 1
<- gen_theta(mu_sigma)
true_theta <- gen_data(true_theta, ipar_a, ipar_d)
X <- d
dX @response <- X
dX@response[[dX@person_id]] <- 990000 + 1:dim(X)[1]
dX
# Unidimensional
set.seed(1)
<- runLinking(dX, method = "FIXEDPAR")
tmp_calib <- runRSSS(
tmp_rsss min_score = 0,
dX, tmp_calib, min_theta = -4.0, max_theta = 4.0, inc = 0.05
)<- tmp_rsss[[as.character(lhs_scale)]]
table_UNI <- table_UNI[, c(1, 4:5)]
table_UNI
# Equipercentile
set.seed(1)
<- runEquateObserved(
tmp smooth = "loglinear",
dX, scale_from = lhs_scale,
scale_to = rhs_scale,
type_to = "theta",
rsss = tmp_rsss,
degrees = list(6, 1),
reps = 1000
)<- tmp$concordance
table_EQP
# Calibrated projection
if (true_corr != 1) {
set.seed(1)
<- runLinking(dX, method = "CP")
tmp_calib <- runRSSS(
tmp_rsss min_score = 0,
dX, tmp_calib, min_theta = -4.0, max_theta = 4.0, inc = 0.05
)<- tmp_rsss[[as.character(lhs_scale)]]
table_CP else {
} <- NULL
table_CP
}
# Evaluate
<- getResponse(dX, lhs_scale)
X_lhs <- eval_table(
table_out
X_lhs, true_theta,
lhs_scale, rhs_scale,
table_UNI, table_EQP, table_CP
)
# Write output
<- sprintf("out_%s_%s.csv", true_corr, idx_trial)
fn write.csv(table_out, file.path("results", fn), row.names = FALSE)
NULL
}
stopCluster(cl)