Correted pcm weights

This commit is contained in:
2025-05-26 16:00:29 +02:00
parent 1a259b6b89
commit d071500c44
2 changed files with 94 additions and 76 deletions

168
R/pcm.R
View File

@ -1,5 +1,5 @@
## File Name: pcm.R ## File Name: pcm.R
## File version: 1.0 ## File version: 1.1
#' Compute Partial Credit Model (PCM) for polytomous and dichotomous items #' Compute Partial Credit Model (PCM) for polytomous and dichotomous items
#' #'
@ -10,7 +10,7 @@
#' @param grp string containing the name of the column where an optional group membership variable is stored in df #' @param grp string containing the name of the column where an optional group membership variable is stored in df
#' @param dif.items vector containing the list of indexes in "items" corresponding to dif items #' @param dif.items vector containing the list of indexes in "items" corresponding to dif items
#' @param type.dif vector containing DIF form for each item specified in dif.items. 1 is homogeneous DIF, 0 is heterogeneous DIF #' @param type.dif vector containing DIF form for each item specified in dif.items. 1 is homogeneous DIF, 0 is heterogeneous DIF
#' @param weights string containing the name of the column where optional weights are stored in df #' @param weights string containing the name of the column where an optional variable containing weights is stored in df
#' @param verbose set to TRUE to print a detailed output, FALSE otherwise #' @param verbose set to TRUE to print a detailed output, FALSE otherwise
#' @param fit string determining the optimization algorithm. Values "ucminf" or "nlminb" ar recommended #' @param fit string determining the optimization algorithm. Values "ucminf" or "nlminb" ar recommended
#' @param method.theta string determining the estimation method for individual latent variable values. Either "eap", "mle" or "wle" #' @param method.theta string determining the estimation method for individual latent variable values. Either "eap", "mle" or "wle"
@ -19,6 +19,7 @@
#' @import PP #' @import PP
#' @export #' @export
pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights=NULL,verbose=T,fit="ucminf",method.theta="eap") { pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights=NULL,verbose=T,fit="ucminf",method.theta="eap") {
##### Detecting errors ##### Detecting errors
@ -64,39 +65,42 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
cat("######################################### FITTING MODEL #########################################\n") cat("######################################### FITTING MODEL #########################################\n")
cat("#################################################################################################\n") cat("#################################################################################################\n")
} }
grp <- NULL grp <- NULL
# prepare data # prepare data
if (is.null(weights)) { if (is.null(weights)) {
df <- df[,c('id',items)] df <- df[,c('id',items)]
} else {
df <- df[,c('id',items,weights)]
}
print(df)
colnames(df)[2:(length(colnames(df)))] <- paste0("item",seq(1,length(colnames(df))-1)) colnames(df)[2:(length(colnames(df)))] <- paste0("item",seq(1,length(colnames(df))-1))
df.long <- reshape(df,v.names=c("item"),direction="long",varying=c(items)) } else {
if (is.null(weights)) { df <- df[,c('id',items,weights)]
colnames(df.long) <- c("id","item","resp") colnames(df)[2:(length(colnames(df)-1))] <- paste0("item",seq(1,length(colnames(df))-1))
} else { }
colnames(df.long) <- c("id","item","resp","weights") df.long <- reshape(df,v.names=c("item"),direction="long",varying=c(items))
} if (is.null(weights)) {
colnames(df.long) <- c("id","item","resp")
nbitems <- length(2:(length(colnames(df)))) nbitems <- length(2:(length(colnames(df))))
maxmod <- max(df[,2:(length(colnames(df)))]) maxmod <- max(df[,2:(length(colnames(df)))])
df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-1),ordered = F) df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-1),ordered = F)
df.long$resp <- factor(df.long$resp,0:maxmod,ordered=T) } else {
df.long$id <- factor(df.long$id) colnames(df.long) <- c("id","weights","item","resp")
# fit pcm nbitems <- length(2:(length(colnames(df))-1))
if (is.null(weights)) { maxmod <- max(df[,2:(length(colnames(df))-1)])
mod <- olmm(resp ~ 0 + ce(item) + re(0|id),data=df.long,family = adjacent(link = "logit")) df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-2),ordered = F)
} else { }
mod <- olmm(resp ~ 0 + ce(item) + re(0|id),data=df.long,family = adjacent(link = "logit"),weights = df.long$weights) df.long$resp <- factor(df.long$resp,0:maxmod,ordered=T)
} df.long$id <- factor(df.long$id)
comod <- coef(mod) # fit pcm
# output results if (is.null(weights)) {
restab <- t(sapply(1:nbitems,function(x) comod[seq(x,length(comod)-1,nbitems)])) mod <- olmm(resp ~ 0 + ce(item) + re(0|id),data=df.long,family = adjacent(link = "logit"))
rownames(restab) <- paste0("item",1:nbitems) } else {
colnames(restab) <- paste0("delta_",1:maxmod) mod <- olmm(resp ~ 0 + ce(item) + re(0|id),data=df.long,family = adjacent(link = "logit"),weights = df.long$weights)
restab.dif <- NULL }
beta <- NULL comod <- coef(mod)
# output results
restab <- t(sapply(1:nbitems,function(x) comod[seq(x,length(comod)-1,nbitems)]))
rownames(restab) <- paste0("item",1:nbitems)
colnames(restab) <- paste0("delta_",1:maxmod)
restab.dif <- NULL
beta <- NULL
} }
# If group # If group
else { else {
@ -114,19 +118,23 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
# prepare data # prepare data
if (is.null(weights)) { if (is.null(weights)) {
df <- df[,c('id',items,"grp")] df <- df[,c('id',items,"grp")]
colnames(df)[2:(length(colnames(df))-1)] <- paste0("item",seq(1,length(colnames(df))-2))
} else { } else {
df <- df[,c('id',items,"grp",weights)] df <- df[,c('id',items,"grp",weights)]
colnames(df)[2:(length(colnames(df))-2)] <- paste0("item",seq(1,length(colnames(df))-3))
} }
colnames(df)[2:(length(colnames(df))-1)] <- paste0("item",seq(1,length(colnames(df))-2))
df.long <- reshape(df,v.names=c("item"),direction="long",varying=c(items)) df.long <- reshape(df,v.names=c("item"),direction="long",varying=c(items))
if (is.null(weights)) { if (is.null(weights)) {
colnames(df.long) <- c("id","grp","item","resp") colnames(df.long) <- c("id","grp","item","resp")
nbitems <- length(2:(length(colnames(df))-1))
maxmod <- max(df[,2:(length(colnames(df))-1)])
df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-2),ordered = F)
} else { } else {
colnames(df.long) <- c("id","grp","item","resp","weights") colnames(df.long) <- c("id","grp","weights","item","resp")
nbitems <- length(2:(length(colnames(df))-2))
maxmod <- max(df[,2:(length(colnames(df))-2)])
df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-3),ordered = F)
} }
nbitems <- length(2:(length(colnames(df))-1))
maxmod <- max(df[,2:(length(colnames(df))-1)])
df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-2),ordered = F)
df.long$resp <- factor(df.long$resp,0:maxmod,ordered=T) df.long$resp <- factor(df.long$resp,0:maxmod,ordered=T)
df.long$id <- factor(df.long$id) df.long$id <- factor(df.long$id)
@ -148,7 +156,11 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
comod <- coef(mod) comod <- coef(mod)
# output results # output results
nbcoef <- nbitems+length(difvar.nonunif) nbcoef <- nbitems+length(difvar.nonunif)
restab <- t(sapply(1:nbcoef,function(x) comod[seq(x,length(comod)-2-length(difvar.unif),nbitems+length(difvar.nonunif))])) if (is.null(weights)) {
restab <- t(sapply(1:nbcoef,function(x) comod[seq(x,length(comod)-2-length(difvar.unif),nbitems+length(difvar.nonunif))]))
} else {
restab <- t(sapply(1:nbcoef,function(x) comod[seq(x,length(comod)-2-length(difvar.unif),nbitems+length(difvar.nonunif))]))
}
difcoef.unif <- NULL difcoef.unif <- NULL
if (length(difvar.unif)>0) { if (length(difvar.unif)>0) {
difcoef.unif <- comod[(length(comod)-length(difvar.unif)):(length(comod)-1)] difcoef.unif <- comod[(length(comod)-length(difvar.unif)):(length(comod)-1)]
@ -204,45 +216,54 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
cat("######################################### FITTING MODEL #########################################\n") cat("######################################### FITTING MODEL #########################################\n")
cat("#################################################################################################\n") cat("#################################################################################################\n")
} }
# prepare data # prepare data
if (is.null(weights)) { if (is.null(weights)) {
df <- df[,c('id',items,"grp")] df <- df[,c('id',items,"grp")]
} else {
df <- df[,c('id',items,"grp",weights)]
}
colnames(df)[2:(length(colnames(df))-1)] <- paste0("item",seq(1,length(colnames(df))-2)) colnames(df)[2:(length(colnames(df))-1)] <- paste0("item",seq(1,length(colnames(df))-2))
df.long <- reshape(df,v.names=c("item"),direction="long",varying=c(items))
if (is.null(weights)) { } else {
colnames(df.long) <- c("id","grp","item","resp") df <- df[,c('id',items,"grp",weights)]
} else { colnames(df)[2:(length(colnames(df))-2)] <- paste0("item",seq(1,length(colnames(df))-3))
colnames(df.long) <- c("id","grp","item","resp","weights") }
} df.long <- reshape(df,v.names=c("item"),direction="long",varying=c(items))
if (is.null(weights)) {
colnames(df.long) <- c("id","grp","item","resp")
nbitems <- length(2:(length(colnames(df))-1)) nbitems <- length(2:(length(colnames(df))-1))
maxmod <- max(df[,2:(length(colnames(df))-1)]) maxmod <- max(df[,2:(length(colnames(df))-1)])
df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-2),ordered = F) df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-2),ordered = F)
df.long$resp <- factor(df.long$resp,0:maxmod,ordered=T) } else {
df.long$id <- factor(df.long$id) colnames(df.long) <- c("id","grp","weights","item","resp")
# fit pcm nbitems <- length(2:(length(colnames(df))-2))
if (is.null(weights)) { maxmod <- max(df[,2:(length(colnames(df))-2)])
mod <- olmm(resp ~ 0 + ge(grp) + ce(item) + re(0|id),data=df.long,family = adjacent(link = "logit"),control=olmm_control(fit=fit)) df.long$item <- factor(df.long$item,levels=seq(1,length(colnames(df))-3),ordered = F)
} else { }
mod <- olmm(resp ~ 0 + ge(grp) + ce(item) + re(0|id),data=df.long,family = adjacent(link = "logit"),control=olmm_control(fit=fit),weights=df.long$weights) df.long$resp <- factor(df.long$resp,0:maxmod,ordered=T)
} df.long$id <- factor(df.long$id)
comod <- coef(mod) # fit pcm
# output results if (is.null(weights)) {
mod <- olmm(resp ~ 0 + ge(grp) + ce(item) + re(0|id),data=df.long,family = adjacent(link = "logit"),control=olmm_control(fit=fit))
} else {
mod <- olmm(resp ~ 0 + ge(grp) + ce(item) + re(0|id),data=df.long,family = adjacent(link = "logit"),control=olmm_control(fit=fit),weights=df.long$weights)
}
comod <- coef(mod)
# output results
if (is.null(weights)) {
restab <- t(sapply(1:nbitems,function(x) comod[seq(x,length(comod)-2,nbitems)])) restab <- t(sapply(1:nbitems,function(x) comod[seq(x,length(comod)-2,nbitems)]))
rownames(restab) <- items_o } else {
colnames(restab) <- paste0("delta_",1:maxmod) restab <- t(sapply(1:nbitems,function(x) comod[seq(x,length(comod)-2,nbitems)]))
restab.dif <- NULL }
beta <- comod[length(comod)-1] rownames(restab) <- items_o
se.beta <- (confint(mod)["grp",2]-beta)/1.96 colnames(restab) <- paste0("delta_",1:maxmod)
beta.ci <- confint(mod)["grp",] restab.dif <- NULL
beta.p <- 2*pnorm(-abs(beta/se.beta)) beta <- comod[length(comod)-1]
beta <- as.numeric(beta) se.beta <- (confint(mod)["grp",2]-beta)/1.96
se.beta <- as.numeric(se.beta) beta.ci <- confint(mod)["grp",]
beta.p <- as.numeric(beta.p) beta.p <- 2*pnorm(-abs(beta/se.beta))
beta <- -1*beta beta <- as.numeric(beta)
beta.ci <- -1*c(beta.ci[2],beta.ci[1]) se.beta <- as.numeric(se.beta)
beta.p <- as.numeric(beta.p)
beta <- -1*beta
beta.ci <- -1*c(beta.ci[2],beta.ci[1])
} }
} }
@ -275,9 +296,6 @@ pcm <- function(df=NULL,items=NULL,grp=NULL,dif.items=NULL,type.dif=NULL,weights
dif.param=restab.dif, dif.param=restab.dif,
theta=theta, theta=theta,
residuals=resid residuals=resid
) )
return(out) return(out)
} }

View File

@ -27,7 +27,7 @@ pcm(
\item{type.dif}{vector containing DIF form for each item specified in dif.items. 1 is homogeneous DIF, 0 is heterogeneous DIF} \item{type.dif}{vector containing DIF form for each item specified in dif.items. 1 is homogeneous DIF, 0 is heterogeneous DIF}
\item{weights}{string containing the name of the column where optional weights are stored in df} \item{weights}{string containing the name of the column where an optional variable containing weights is stored in df}
\item{verbose}{set to TRUE to print a detailed output, FALSE otherwise} \item{verbose}{set to TRUE to print a detailed output, FALSE otherwise}