Nothing
set.seed(7485)
library(mstDIF)
# function to simulate 3PL data
P3pl <- function(theta, a = rep(1, length(b)), b, c = rep(0, length(b)),
generateResponses = FALSE, seed = 1){
nItem <- length(b)
nPerson <- length(theta)
A <- rep(1, nPerson) %o% a
B <- rep(1, nPerson) %o% b
C <- rep(1, nPerson) %o% c
Theta <- theta %o% rep(1, nItem)
Theta_B <- Theta - B
P <- C + (1 - C) * (1 + exp(- A * Theta_B))^(-1)
if(generateResponses){
set.seed(seed)
resp <- (P > matrix(runif(nPerson * nItem),
ncol = nItem, nrow = nPerson)) * 1
return(resp)
} else (return(P))
}
# item parameters
b <- seq(-2, 1.8, length.out = 3)
a <- runif(length(b), 1, 1.5)
c <- rep(0, length(b))
# person parameters
nPerson <- 20
theta <- rnorm(nPerson)
# generate responses
resp <- P3pl(theta, a, b, generateResponses = TRUE)
# create orders
metric <- rnorm(nPerson)
factor <- factor(sample(1:3, size = nPerson, replace = TRUE))
order <- ordered(sample(1:3, size = nPerson, replace = TRUE))
# compute scores
scores <- mstDIF:::get_scores(resp, a, b, c, theta,
slope_intercept = FALSE, return_terms = TRUE)
scores
# compute process
process <- apply(mstDIF:::scale_scores(scores$scores,
decorrelate = TRUE)[order(metric),],
2, cumsum)
process
# Compare with strucchange
library(strucchange)
gefp <- gefp(x = scale(scores$scores, scale = FALSE), fit = NULL,
scores = function(x) {x}, order.by = metric)
gefp$process[-1,]
# Is the cumsum score process the same?
max(abs(gefp$process[-1,] - process))
stopifnot(all(abs(gefp$process[-1,] - process) < 1e-8))
# bootstrap vs sctest
test_b <- bootstrap_sctest(resp = resp, a = a, b = b, nSample = 15,
item_selection = 1,
DIF_covariate = list(m1 = metric,
m2 = metric,
m3 = metric,
f = factor,
o1 = order,
o2 = order),
statistic = c("auto", "CvM", "maxLM", "auto", "auto", "WDMo"),
decorrelate = TRUE, theta = theta)
test_b$p
test_b$statistic
sctest_m1 <- sctest(x = scale(scores$scores, scale = FALSE),
scores = function(x) {x}, parm = c(1, 4), order.by = metric)
sctest_m2 <- sctest(x = scale(scores$scores, scale = FALSE),
scores = function(x) {x}, parm = c(1, 4), order.by = metric,
functional = "CvM")
sctest_m3 <- sctest(x = scale(scores$scores, scale = FALSE),
scores = function(x) {x}, parm = c(1, 4), order.by = metric,
functional = "maxLM")
sctest_f <- sctest(x = scale(scores$scores, scale = FALSE),
scores = function(x) {x}, parm = c(1, 4), order.by = factor,
functional = "LMuo")
sctest_o <- sctest(x = scale(scores$scores, scale = FALSE),
scores = function(x) {x}, parm = c(1, 4), order.by = order,
functional = "maxLMo")
sctest_o2 <- sctest(x = scale(scores$scores, scale = FALSE),
scores = function(x) {x}, parm = c(1, 4), order.by = order,
functional = "WDMo")
# Are the test statistics the same
abs(test_b$statistic / rbind(sctest_m1$statistic,
sctest_m2$statistic,
sctest_m3$statistic,
sctest_f$statistic,
sctest_o$statistic,
sctest_o2$statistic) - 1)
stopifnot(all(abs(test_b$statistic / rbind(sctest_m1$statistic,
sctest_m2$statistic,
sctest_m3$statistic,
sctest_f$statistic,
sctest_o$statistic,
sctest_o2$statistic) - 1) < 1e-8))
# permutation vs bootstrap
test_p <- permutation_sctest(resp = resp, a = a, b = b, nSample = 15,
item_selection = 1,
DIF_covariate = list(m1 = metric,
m2 = metric,
m3 = metric,
f = factor,
o1 = order,
o2 = order),
statistic = c("auto", "CvM", "maxLM", "auto", "auto", "WDMo"),
decorrelate = TRUE, theta = theta)
test_p$p
test_p$statistic
test_b$statistic
test_p$statistic
# Are the test statistics the same
stopifnot(all(test_b$statistic == test_p$statistic))
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.