Getting started with regDIF

Using LASSO regularization to evaluate measurement bias in a 2-parameter logistic IRT model.

In this example, data in ida were generated to mimic an integrative data analysis, where data are pooled across multiple studies and the measurement model is evaluated for between-study and within-study (e.g., gender, age) measurement bias. These data include 6 item responses (binary) and 3 background characteristics – namely, age (continuous, centered), gender (categorical, effect-coded), and study (categorical, effect-coded).

DIF was generated to be on items 2 (age, gender, study), 3 (age, gender, study), 4 (age), and 5 (gender, study), for both intercepts and slopes.

library(regDIF)
head(ida)

##   item1 item2 item3 item4 item5 item6 age gender study
## 1     0     0     0     0     0     0  -2     -1    -1
## 2     0     0     0     0     0     0   0     -1    -1
## 3     0     0     0     0     0     0   3     -1    -1
## 4     0     1     1     1     1     1   1     -1    -1
## 5     0     0     0     0     0     0  -2     -1    -1
## 6     1     0     0     0     0     0   1     -1    -1

The item response data must first be separated from the predictor data (background variables) before running regDIF(). A single value of the tuning parameter, tau = 2, is then fit to the data.

item.data <- regDIF::ida[,1:6]
pred.data <- regDIF::ida[,7:9]
fit <- regDIF(item.data, pred.data, tau = 2)

summary(fit)

## Call:
## regDIF(item.data = item.data, pred.data = pred.data, tau = 2)
## 
## Optimal model (out of 1):
##      tau      bic 
##    2.000 4161.704 
## 
## Non-zero DIF effects:
##    item1.int.age item1.int.gender  item1.int.study item2.int.gender 
##          -0.0723           0.0274           0.1315          -0.0232 
##  item2.int.study    item3.int.age item3.int.gender    item5.int.age 
##          -0.0291          -0.3897           0.2731          -0.6745 
## item5.int.gender  item5.int.study    item6.int.age item6.int.gender 
##          -0.6913           0.5488          -0.1515           0.0293 
##  item6.int.study    item1.slp.age item1.slp.gender  item1.slp.study 
##           0.0590           0.0012          -0.0091          -0.4055 
##    item2.slp.age item2.slp.gender  item2.slp.study    item3.slp.age 
##           0.0321           0.0881          -0.1406          -0.0255 
## item3.slp.gender  item3.slp.study    item4.slp.age  item4.slp.study 
##          -0.2028          -0.5129           0.0631          -0.6216 
## item5.slp.gender  item5.slp.study item6.slp.gender  item6.slp.study 
##          -0.5273          -0.0769          -0.0012          -0.3516

The summary() function shows that no DIF effects remain in the model. Only the latent variable parameters and base item parameters, which were not penalized at all, are estimated to be non-zero. This is shown by using the coef method.

coef(fit)

## $tau
## [1] 2
## 
## $impact
##                [,1]
## mean.age     0.9964
## mean.gender -0.1252
## mean.study   0.9213
## var.age      0.2732
## var.gender  -0.0085
## var.study    0.4140
## 
## $base
##               [,1]
## item1.int. -0.6185
## item2.int. -1.0141
## item3.int. -1.4054
## item4.int. -1.9915
## item5.int. -1.7821
## item6.int. -2.1775
## item1.slp.  0.9661
## item2.slp.  0.9788
## item3.slp.  1.2127
## item4.slp.  1.4989
## item5.slp.  1.4708
## item6.slp.  1.3571
## 
## $dif
##                     [,1]
## item1.int.age    -0.0723
## item1.int.gender  0.0274
## item1.int.study   0.1315
## item2.int.age     0.0000
## item2.int.gender -0.0232
## item2.int.study  -0.0291
## item3.int.age    -0.3897
## item3.int.gender  0.2731
## item3.int.study   0.0000
## item4.int.age     0.0000
## item4.int.gender  0.0000
## item4.int.study   0.0000
## item5.int.age    -0.6745
## item5.int.gender -0.6913
## item5.int.study   0.5488
## item6.int.age    -0.1515
## item6.int.gender  0.0293
## item6.int.study   0.0590
## item1.slp.age     0.0012
## item1.slp.gender -0.0091
## item1.slp.study  -0.4055
## item2.slp.age     0.0321
## item2.slp.gender  0.0881
## item2.slp.study  -0.1406
## item3.slp.age    -0.0255
## item3.slp.gender -0.2028
## item3.slp.study  -0.5129
## item4.slp.age     0.0631
## item4.slp.gender  0.0000
## item4.slp.study  -0.6216
## item5.slp.age     0.0000
## item5.slp.gender -0.5273
## item5.slp.study  -0.0769
## item6.slp.age     0.0000
## item6.slp.gender -0.0012
## item6.slp.study  -0.3516

Now that the data have been properly specified in regDIF, a more thorough investigation of DIF is warranted. The regDIF() function defaults to estimating 100 values of the tuning parameter, starting with a value large enough to penalize all DIF effects to zero. However, for brevity, only 10 values of tau are specified with the num.tau argument and we reduce the tolerance for convergence using the control argument.

fit2 <- regDIF(item.data, pred.data, num.tau = 10, control = list(tol = 1e-3))

fit2

## Call:
## regDIF(item.data = item.data, pred.data = pred.data, num.tau = 10, 
##     control = list(tol = 0.001))
## 
## regDIF results:
##          tau      bic
## 1  40.926501 4133.065
## 2  28.743990 4110.741
## 3  19.256227 4104.520
## 4  12.126371 4108.451
## 5   7.017576 4118.523
## 6   3.592999 4138.299
## 7         NA       NA
## 8         NA       NA
## 9         NA       NA
## 10        NA       NA

By printing the model object, 10 rows of results appear, one for each value of the tuning parameter. The first thing to notice is that 4 rows are missing. This occurs because regDIF() automatically stops model-fitting when a small value of tau would produce a non-identified model. In focusing attention to the BIC column, it is evident that the smallest value occurs well before the model would be non-identified. This is an encouraging result. The summary() function may be used again, which, with multiple values of the tuning parameter fit to the data, produces non-zero DIF effects corresponding to the model with the minimum value of BIC.

summary(fit2)

## Call:
## regDIF(item.data = item.data, pred.data = pred.data, num.tau = 10, 
##     control = list(tol = 0.001))
## 
## Optimal model (out of 10):
##        tau        bic 
##   19.25623 4104.52000 
## 
## Non-zero DIF effects:
##    item5.int.age item5.int.gender  item5.int.study  item4.slp.study 
##          -0.0352          -0.3535           0.4341          -0.0777 
##    item5.slp.age item5.slp.gender 
##          -0.0450          -0.2071

A plot of the regularization path also shows the remaining DIF effects.

plot(fit2)

To produce other model results, the fit2 object contains lists of the impact (latent variable) parameters, base (intercept and slope) item parameters, and DIF parameters for all values of the tuning parameter. For instance, the impact parameters are printed below.

fit2$impact

##                [,1]    [,2]    [,3]    [,4]    [,5]    [,6] [,7] [,8] [,9]
## mean.age     0.7042  0.7160  0.7258  0.7663  0.8800  0.9404   NA   NA   NA
## mean.gender -0.1739 -0.1260 -0.1068 -0.0983 -0.1170 -0.1283   NA   NA   NA
## mean.study   0.9585  0.9271  0.9222  0.9262  0.9749  0.9747   NA   NA   NA
## var.age      0.4317  0.4686  0.4905  0.4734  0.4046  0.3171   NA   NA   NA
## var.gender  -0.1006 -0.0546 -0.0405 -0.0419 -0.0234 -0.0089   NA   NA   NA
## var.study   -0.1243 -0.1394 -0.1171 -0.0358  0.1544  0.3372   NA   NA   NA
##             [,10]
## mean.age       NA
## mean.gender    NA
## mean.study     NA
## var.age        NA
## var.gender     NA
## var.study      NA

EAP scores and standard deviations may also be produced.

lapply(fit2$eap, head)

## $scores
##            [,1]       [,2]       [,3]       [,4]       [,5]       [,6] [,7]
## [1,] -1.8130833 -1.8157221 -1.8217953 -1.8623577 -1.9862165 -2.0174588   NA
## [2,] -1.3342935 -1.3227238 -1.3120476 -1.2814767 -1.2176086 -1.1280142   NA
## [3,] -1.2256347 -1.2560766 -1.2536711 -1.1151152 -0.7702964 -0.4691113   NA
## [4,]  2.0793250  2.0869953  2.0768871  2.0070648  1.7744336  1.5654123   NA
## [5,] -1.8130833 -1.8157221 -1.8217953 -1.8623577 -1.9862165 -2.0174588   NA
## [6,] -0.3805602 -0.3552358 -0.3580264 -0.3680989 -0.3221837 -0.2510906   NA
##      [,8] [,9] [,10]
## [1,]   NA   NA    NA
## [2,]   NA   NA    NA
## [3,]   NA   NA    NA
## [4,]   NA   NA    NA
## [5,]   NA   NA    NA
## [6,]   NA   NA    NA
## 
## $sd
##           [,1]      [,2]      [,3]      [,4]      [,5]      [,6] [,7] [,8] [,9]
## [1,] 0.7966229 0.7758801 0.7568670 0.7374181 0.6970773 0.6676399   NA   NA   NA
## [2,] 0.9317806 0.9266071 0.9148572 0.8800996 0.7918539 0.7128887   NA   NA   NA
## [3,] 1.1450472 1.1685206 1.1736104 1.1040151 0.9101280 0.7378248   NA   NA   NA
## [4,] 0.6973359 0.6935337 0.6851802 0.6579151 0.5871376 0.5106761   NA   NA   NA
## [5,] 0.7966229 0.7758801 0.7568670 0.7374181 0.6970773 0.6676399   NA   NA   NA
## [6,] 0.8553530 0.8534927 0.8463721 0.8132854 0.7162173 0.6209409   NA   NA   NA
##      [,10]
## [1,]    NA
## [2,]    NA
## [3,]    NA
## [4,]    NA
## [5,]    NA
## [6,]    NA

Finally, when data include a large number of items, observations, and predictors, regDIF() can run relatively slowly. An alternative approach, which yields much faster results, is to provide an observed proxy for the latent scores. In the case of binary data, this might be sum scores. Note that using observed proxy scores is identical to performing a multivariate regression, where the item responses are regressed on the proxy scores and background variables. (The proxy scores are simultaneously regressed on the background variables as well.)

fit3 <- regDIF(item.data, pred.data, prox.data = rowSums(item.data), num.tau = 20)

summary(fit3)

## Call:
## regDIF(item.data = item.data, pred.data = pred.data, prox.data = rowSums(item.data), 
##     num.tau = 20)
## 
## Optimal model (out of 20):
##        tau        bic 
##   14.84459 3553.86250 
## 
## Non-zero DIF effects:
## item3.int.gender    item4.int.age    item5.int.age item5.int.gender 
##           0.1243           0.1735          -0.0255          -0.5186 
##  item5.int.study item2.slp.gender  item4.slp.study    item5.slp.age 
##           0.6472           0.0470          -0.1488          -0.0752 
## item5.slp.gender 
##          -0.1155

The results show more DIF on both the intercepts (Items 3 and 6) and slopes (Items 2, 3, and 5).

More modeling possibilities with regDIF.

In addition to LASSO, other penalty functions are possible with regDIF(). For instance, the elastic net penalty combines LASSO and ridge functions, which is useful when many correlated predictors are evaluated for DIF. The elastic net is controlled by a second tuning parameter, alpha, and defaults to alpha = 1, corresponding to the LASSO penalty. In contrast, alpha = 0 corresponds to the ridge penalty. When alpha is between 0 and 1, however, the elastic net is used to perform DIF selection. For brevity, observed proxy scores are used in all model fitting below.

fit_net <- regDIF(item.data, pred.data, prox.data = rowSums(item.data), num.tau = 20, alpha = .5)

summary(fit_net)

## Call:
## regDIF(item.data = item.data, pred.data = pred.data, prox.data = rowSums(item.data), 
##     num.tau = 20, alpha = 0.5)
## 
## Optimal model (out of 20):
##        tau        bic 
##   18.54054 3573.79720 
## 
## Non-zero DIF effects:
## item1.int.gender    item2.int.age    item3.int.age item3.int.gender 
##           0.0578           0.0578          -0.0651           0.1871 
##    item4.int.age    item5.int.age item5.int.gender  item5.int.study 
##           0.2223          -0.1376          -0.4994           0.5594 
## item6.int.gender item2.slp.gender  item3.slp.study    item4.slp.age 
##           0.0277           0.1513          -0.0830           0.0055 
##  item4.slp.study    item5.slp.age item5.slp.gender 
##          -0.2035          -0.1309          -0.2121

The final elastic net results yield the same DIF effects as the LASSO results, although the amount of penalization is greater for elastic net (i.e., larger tau).

Other penalty functions include the minimax concave penalty (MCP) and the group extensions of LASSO and MCP, which penalize the intercept and slope DIF effects in tandem. The group LASSO function is shown below.

fit_grp_mcp <- regDIF(item.data, pred.data, prox.data = rowSums(item.data), num.tau = 20, pen.type = "grp.mcp")

summary(fit_grp_mcp)

## Call:
## regDIF(item.data = item.data, pred.data = pred.data, prox.data = rowSums(item.data), 
##     pen.type = "grp.mcp", num.tau = 20)
## 
## Optimal model (out of 20):
##          tau          bic 
##    0.4324324 3530.2034000 
## 
## Non-zero DIF effects:
##    item4.int.age item5.int.gender  item5.int.study    item4.slp.age 
##           0.1793          -0.6400           0.8894           0.1182 
## item5.slp.gender  item5.slp.study 
##          -0.4197           0.2578

Although the MCP results appear largely the same as the LASSO results, the group MCP function included both the intercept and slope for each background variable remaining in the final model.

In summary, the regDIF R package provides a flexible implementation of using regularization to identify DIF across multiple background characteristics.

Please reach out to wbelzak@gmail.com for any questions, and remember to cite regDIF in your work. Thank you kindly!

## This package can be cited as:
## 
##   Belzak, W. C. M. (2023). The regDIF R Package: Evaluating Complex
##   Sources of Measurement Bias Using Regularized Differential Item
##   Functioning. Structural Equation Modeling: A Multidisciplinary
##   Journal, 74-984, DOI: 10.1080/10705511.2023.2170235
## 
## A BibTeX entry for LaTeX users is
## 
##   @Article{,
##     journal = {Structural Equation Modeling: A Multidisciplinary Journal},
##     entry = {manual},
##     title = {The regDIF R Package: Evaluating Complex Sources of Measurement Bias Using Regularized Differential Item Functioning},
##     author = {William C. M. Belzak},
##     organization = {Duolingo},
##     address = {Pittsburgh, PA},
##     year = {2023},
##     url = {https://doi.org/10.1080/10705511.2023.2170235},
##   }
## 
## This free open-source software implements academic research by the
## authors and co-workers. If you use it, please support the project by
## citing the appropriate journal articles.