RMSD index for non-invariance

[Erin717]

Dear Dr. Chalmers,

I am new to IRT and your package has been very helpful!

I was wondering if there is a way to calculate RMSD as a measurement of item non-invariance between groups. I read about it from PISA report and some recent literarure on using RMSD to evaluate non-invariance. However, I could not find a way to calculate this. Could you please point me the way to go about it?

Specifically, I'd like to replicate the PISA scaling procedure, where I would first estimate a fully constrained model between groups, the calculate RMSD each group and item. 

Thank you very much in advance!

Erin 

[Phil Chalmers]

Hi Erin,

I'm not quite sure what you're asking for as I'm not familiar with what is being used to compute an RMSD statistic (could be anything really, and at different levels of the analysis). Do you have any references/equations on hand that would suggest the type of invariance quantification you're after? 

Phil

--
You received this message because you are subscribed to the Google Groups "mirt-package" group.
To unsubscribe from this group and stop receiving emails from it, send an email to mirt-package...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/mirt-package/5b0bb5a7-cb25-40c3-9159-0761b4c22acen%40googlegroups.com.

[Erin717]

Hi Phil,

Thank you for your reply. As I understood it, the RMSD statistics was used since PISA 2015 to evaluate item-level non-invariance across groups. 

There are some papers mentioned it as well, please see the following list:

Lee & Von Davier, 2020, Improving measurement properties of the PISA home possessions scale through partial invariance modeling. Psychological test and assessment modeling 62:55-83 (https://www.researchgate.net/publication/339886973_Improving_measurement_properties_of_the_PISA_home_possessions_scale_through_partial_invariance_modeling)

Buchholz & Hartig, 2019, Comparing attitude across groups: An IRT based item-fit statistic for the analysis of measurement invariance

I have also attached the PISA 2015 technical report in this email. The formula is on Page 296.

Thank you very much! 

Erin

[Phil Chalmers]

Thanks. The RMSD statistic in this article suggests the fit:

- "In a first step, groups were defined whereas every country or multiple, sufficiently large samples of examinees taking the same questionnaire
language version within the country formed one group each."

- " In a second step, international item and person parameters were estimated based on all examinees across all groups"

So

RMSD = √∫ [P₀(θ) - P_M(θ)]² f(θ) dθ

is a kind of badness-of-fit measure that could suggest DIF, where P_M(θ) is the model-implied probability trace lines from a baseline constrained model (e.g., fitting a constrained multiple-group IRT model using all the data to create a pooled model across groups, though more flexible DIF-based models could be fitted later). P₀(θ) is based on the "pseudo E-table" after taking subsets of the data pertaining to the group, and I assume keeping the parameters fixed at the originally pooled MG model. f(θ) was not really specified explicitly, so I'm not sure which exact form was used by any of the authors (it could be based on the standardized version of the Etable, like the DRF statistics do, could be based on a standard normal distribution; let's assume the latter for now). Also, for models with > 2 categories an RMSD would be needed for each response curve (for dichotomous items the RMSD is actually redundant for both categories).

This measure is not available in mirt (I'd never seen it before), though I've rolled a potential version below, complete with some examples. Feel free to give it a try and provide feedback, and if it works well I'll consider adding it to the package as a separate DIF function. Note that you'll need to install the latest dev version of mirt first as the extract.group() function was recently modified and would cause this code to fail with older versions of mirt. Cheers.

Phil

`````

#' RMSD statistics to quantify category-level DIF
#'
#'
#' @param dat complete dataset
#' @param group group membership vector
#' @param pooled_mod a multiple-group model (used to compute the model-implied
#'   probability in the goodness-of-fit test)
RMSD_DIF <- function(dat, group, pooled_mod){
    stopifnot(nrow(dat) == length(group))
    which.groups <- extract.mirt(pooled_mod, 'groupNames')
    ret <- vector('list', length(which.groups))
    names(ret) <- which.groups
    for(which.group in which.groups){
        smod <- extract.group(pooled_mod, which.group)
        nfact <- extract.mirt(smod, 'nfact')
        stopifnot(nfact == 1L)
        sv <- mod2values(smod)
        sv$est <- FALSE

        # make wider grid for better numerical integration
        Theta <- matrix(seq(-6,6,length.out = 201))

        mod_g <- mirt(subset(dat, group == which.group), nfact,
                      itemtype = extract.mirt(smod, 'itemtype'),
                      pars = sv, technical = list(storeEtable=TRUE, customTheta=Theta))
        Etable <- mod_g@Internals$Etable[[1]]$r1

        # standard normal dist for theta
        f_theta <- dnorm(Theta)
        f_theta <- as.vector(f_theta / sum(f_theta))
        itemloc <- extract.mirt(mod_g, 'itemloc')
        which.items <- 1L:ncol(dat)
        ret2 <- vector('list', ncol(dat))
        names(ret2) <- extract.mirt(mod_g, 'itemnames')
        for(i in seq_len(length(which.items))){
            pick <- itemloc[which.items[i]]:(itemloc[which.items[i]+1L] - 1L)
            O <- Etable[ ,pick]
            P_o <- O / rowSums(O)
            item <- extract.item(smod, which.items[i])
            P_e <- probtrace(item, Theta)
            ret2[[i]] <- sqrt(colSums((P_o - P_e)^2 * f_theta))
        }
        class(ret2) <- 'mirt_list'
        ret[[which.group]] <- ret2
    }
    ret
}

#----------------------------------------------------------------------
#----- generate some data
set.seed(12345)
a <- a2 <- matrix(abs(rnorm(15,1,.3)), ncol=1)
d <- d2 <- matrix(rnorm(15,0,.7),ncol=1)

# item 1 has DIF
d2[1] <- d[1] - .5
a2[1] <- a[1] + 1

itemtype <- rep('2PL', nrow(a))
N <- 1000
dataset1 <- simdata(a, d, N, itemtype)
dataset2 <- simdata(a2, d2, N, itemtype)
dat <- rbind(dataset1, dataset2)
group <- c(rep('D1', N), rep('D2', N))

#-----

# fully pooled model
pooled_mod <- multipleGroup(dat, 1, group=group, invariance = colnames(dat))
coef(pooled_mod, simplify=TRUE)

RMSD_DIF(dat, group=group, pooled_mod)

# more freely estimated model (item 1 has 2 parameters estimated)
MGmod <- multipleGroup(dat, 1, group=group,
                             invariance = colnames(dat)[-1])
coef(MGmod, simplify=TRUE)

# RMSD in item.1 now reduced (MG model accounts for DIF)
RMSD_DIF(dat, group=group, MGmod)



#################
# polytomous example
set.seed(12345)
a <- a2 <- matrix(rlnorm(20,.2,.3))

# for the graded model, ensure that there is enough space between the intercepts,
# otherwise closer categories will not be selected often (minimum distance of 0.3 here)
diffs <- t(apply(matrix(runif(20*4, .3, 1), 20), 1, cumsum))
diffs <- -(diffs - rowMeans(diffs))
d <- d2 <- diffs + rnorm(20)

# item 1 has slope + dif for first intercept parameter
d2[1] <- d[1] - .5
a2[1] <- a[1] + 1

itemtype <- rep('graded', nrow(a))
N <- 1000
dataset1 <- simdata(a, d, N, itemtype)
dataset2 <- simdata(a2, d2, N, itemtype)
dat <- rbind(dataset1, dataset2)
group <- c(rep('D1', N), rep('D2', N))

#-----

# fully pooled model
pooled_mod <- multipleGroup(dat, 1, group=group, invariance = colnames(dat))
coef(pooled_mod, simplify=TRUE)

RMSD_DIF(dat, group=group, pooled_mod)

# more freely estimated model (item 1 has more parameters estimated)
MGmod <- multipleGroup(dat, 1, group=group,
                       invariance = colnames(dat)[-1])
coef(MGmod, simplify=TRUE)

# RMSD in first two categories of item.1 now reduced (MG model accounts for DIF)
RMSD_DIF(dat, group=group, MGmod)

````

To view this discussion on the web visit https://groups.google.com/d/msgid/mirt-package/1a041610-5c96-4e61-bc40-dc83332f8820n%40googlegroups.com.

[QI HUANG]

Dear Phil,

I would like to follow up on this conversation. I have tried both the code you shared and the RMSD_DIF() function now built in mirt package on the same dataset and same pooled model. However, they give me different results... I am wondering if you would like kindly let me know which part you changed when composing this bulit-in function, compared to the code you shared above? As the built-in function could not accommodate some created item types (e.g., Logistic Positive Exponent model), I am trying to calculate RMSD by myself for those models but still hope to be consistent with the default procedures in mirt. So I will be extremely grateful if you would like to elaborate a bit more on this.

Thank you very much for considering this!

Sincerely,

Qi