---
title: "DIF detection across multiple groups"
output: 
  rmarkdown::html_document:
    toc: true
    toc_float: true
    number_sections: true
    theme: spacelab
    highlight: pygments
bibliography: references.bib
reference-section-title: "References"
vignette: >
  %\VignetteIndexEntry{DIF detection across multiple groups}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---


This document shows the usage of function multiDIF() to implement the method proposed by @Battauz2025.
It could be used to analyze:

* the effect of covariates at the within-group item level on the item parameters (e.g. item position),
* the differences in the item parameters across groups.


# Example of covariate effects: the item position

Load the package equateMultiple and the data.

```{r, message=FALSE}
library("equateMultiple")
data("timss19")
```

Estimate the item parameters separately for each booklet using the R package mirt.

```{r, message=FALSE}
library(mirt)
est.mods <- list()
for (i in 1:14)
  est.mods[[i]] <- mirt(timss19$data[[i]][, -(1:3)], 1, itemtype = "2PL", SE = TRUE)
```

Extract item parameter estimates and covariance matrices.
```{r, message=FALSE, results='hide'}
mods <- modIRT(est.mods = est.mods, names = names(timss19$data))
```

For comparison, the equating coefficients are also estimated ignoring the item position effect using the independent likelihood method
```{r}
eqINDLIK <- multiec(mods = mods, method = "indlik") 
summary(eqINDLIK)
```

Estimation of position effects on the item parameters
```{r}
dif_positions <- multiDIF(mods = mods, se = TRUE, X.DIF = timss19$positions)
summary(dif_positions)
```

Comparison of the log-likelihood without and with DIF effects
```{r}
logLik(eqINDLIK)
logLik(dif_positions)
```


# DIF effects across groups

Load the data from the equateIRT package.
```{r}
data("dataDIF", package = "equateIRT")
data1 <- dataDIF[dataDIF$group == 1, 1:20]
data2 <- dataDIF[dataDIF$group == 2, 1:20]
data3 <- dataDIF[dataDIF$group == 3, 1:20]
```

Estimate the item parameters separately for each group.
```{r}
mod1 <- mirt(data1, SE = TRUE)
mod2 <- mirt(data2, SE = TRUE)
mod3 <- mirt(data3, SE = TRUE)
mods <- list(mod1, mod2, mod3)
```

Extract item parameter estimates and covariance matrices.
```{r, message=FALSE, results='hide'}
mods_extract <- modIRT(mods)
```

Estimate the equating coefficients (without accounting for DIF)
```{r}
eqINDLIK <- multiec(mods_extract, method = "indlik")
summary(eqINDLIK)
```

For fixed values of the tuning parameter $\lambda$ estimate DIF effects and equating coefficients.
```{r}
lambda <- seq(1, 100, by = 0.5)
out <- list()
for (i in 1:length(lambda))
{
  dif_lambda <- multiDIF(mods = mods_extract, start = eqINDLIK, 
    itempar.DIF = TRUE, lambda_mcp = lambda[i])
  out[[i]] <- dif_lambda
}
```

Select the tuning parameter $\lambda$ on the basis of BIC
```{r}
bic <- lapply(out, function(x) x$BIC)
sel_bic <- which.min(bic)
dif_param <- out[[sel_bic]]
summary(dif_param)
```

Only the first item shows DIF in both the difficulty and discrimination parameters.