## File Name: linking.haberman.lq.R
## File Version: 0.202
linking.haberman.lq <- function(itempars, pow=2, eps=1e-3, a_log=TRUE,
use_nu=FALSE, est_pow=FALSE, lower_pow=.1, upper_pow=3)
{
CALL <- match.call()
s1 <- Sys.time()
#--- process item parameters
res <- linking_proc_itempars(itempars=itempars)
itempars <- res$itempars
NS <- res$NS
NI <- res$NI
items <- res$items
studies <- res$studies
wgtM <- res$wgtM
aM <- res$aM
bM <- res$bM
est_pars <- res$est_pars
weights_exist <- res$weights_exist
a.orig <- aM
b.orig <- bM
itempars1 <- itempars
items <- sirt_sup(x=itempars1[,2])
studies <- sirt_sup(x=itempars1[,1])
G <- length(studies)
I <- length(items)
ind_items <- match(itempars1[,2], items)
ind_studies <- match(itempars1[,1], studies)
IP <- nrow(itempars1)
X0 <- matrix(0, nrow=IP, ncol=G-1+I)
colnames(X0) <- c(studies[-1], items)
w <- itempars1[,5]
res_optim <- list()
eps_vec <- 10^seq(0,-10, by=-.5)
item <- data.frame(item=items)
#**** estimation of A
y <- itempars1[,3]
if (a_log){
y <- log(y)
}
X <- X0
for (gg in 2L:G){
X[ ind_studies==gg, gg-1] <- 1
}
for (ii in 1L:I){
X[ ind_items==ii, ii+G-1] <- 1
}
#- fit
res_optim$slopes <- mod0 <- lq_fit(y=y, X=X, w=w, pow=pow, eps=eps,
eps_vec=eps_vec, est_pow=est_pow, lower_pow=lower_pow,
upper_pow=upper_pow)
coef0 <- mod0$coefficients
pow_slopes <- mod0$pow
ind_groups <- 1L:(G-1)
coef0_A <- coef0[ind_groups]
a_joint <- coef0[-c(ind_groups)]
ar <- y - X %*% coef0
if (a_log){
coef0_A <- exp(c(0,coef0_A))
a_joint <- exp(a_joint)
ar <- a_joint[ind_items]*( exp(ar) - 1 )
} else {
coef0_A <- 1+c(0,coef0_A)
}
resid <- itempars1
resid$adif <- ar
At <- coef0_A
transf.personpars <- data.frame(study=studies, A_theta=coef0_A, se_A_theta=NA)
item$a <- a_joint
#**** estimation of B
y <- itempars1[,4] * coef0_A[ ind_studies ]
if (use_nu){
y <- -itempars1[,4]*itempars1[,3]
for (gg in 2L:G){
ind_gg <- which(ind_studies==gg)
X[ ind_gg, gg-1] <- itempars1[ ind_gg,3]/coef0_A[gg]
}
}
#- fit
res_optim$intercepts <- mod0 <- lq_fit(y=y, X=X, w=w, pow=pow, eps=eps,
eps_vec=eps_vec, est_pow=est_pow, lower_pow=lower_pow,
upper_pow=upper_pow)
coef0 <- mod0$coefficients
pow_intercepts <- mod0$pow
coef0_B <- -coef0[ind_groups]
b_joint <- coef0[-c(ind_groups)]
if (use_nu){
coef0_B <- -coef0_B
b_joint <- -b_joint
}
resid$b_resid <- y - X %*% coef0
Bt <- coef0_B <- c(0, coef0_B)
transf.personpars$B_theta <- coef0_B
transf.personpars$se_B_theta <- NA
rownames(transf.personpars) <- NULL
item$b <- b_joint
#- include joint item parameters
resid$a_joint <- item[ind_items, 'a']
resid$b_joint <- item[ind_items, 'b']
# transformation for item parameters
transf.itempars <- data.frame( study=studies, At=1/At, se_At=NA,
Bt=-Bt, se_Bt=NA )
converged <- res_optim$intercepts$converged & res_optim$slopes$converged
s2 <- Sys.time()
time <- list(s1=s1, s2=s2)
#-- output list
description <- 'Linking based on L_q distance according to Haberman (2009)'
res <- list( transf.personpars=transf.personpars, transf.itempars=transf.itempars,
pow=pow, eps=eps, item=item, resid=resid, description=description,
converged=converged, a_log=a_log, use_nu=use_nu, est_pow=est_pow,
pow_slopes=pow_slopes, pow_intercepts=pow_intercepts,
CALL=CALL, time=time)
class(res) <- 'linking.haberman.lq'
return(res)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.