Abstract
Mixtures of item response theory (IRT) models have been proposed as a technique to explore response patterns in test data related to cognitive strategies, instructional sensitivity, and differential item functioning (DIF). Estimation proves challenging due to difficulties in identification and questions of effect size needed to recover underlying structure. In particular, the impact of covariates for examinees in estimation has not been systematically explored. The goal of this study is to carry out a systematically designed simulation study to investigate the performance of mixture Rasch model (MRM) under Bayesian estimation. Insights and suggestions on model application and model estimation are discussed. The foci of this study are to use a flexible logistic regression structure to include examinees’ covariate in MRM, to study Markov chain Monte Carlo (MCMC) estimation behavior in light of effect size, and to provide an effective and applicable method for dealing with chain switching.
Keywords
Item response theory (IRT) models are used in educational assessment to compare examinees in terms of their propensities to perform well on test items and test items in terms of their dependence on this propensity (Yen & Fitzpatrick, 2006). Lack of fit can arise when relationships between proficiency and test items differ across classes of examinees, for reasons such as differential use of solution strategies, educational experiences, and cultural backgrounds. A subset of geometry items in a mathematics test may be relatively easier for boys than girls, for example, compared with algebra and probability items. This is an instance of differential item functioning or DIF (Holland & Wainer, 1993). A mixture IRT model is said to hold when different IRT models hold in different examinee classes, and examinees’ class memberships are not observed (Rost, 1997).
This article addresses a special case of IRT mixture called Mixture Rasch Model (MRM; Rost, 1990), in the case when covariates are observed for examinees and may be related to latent class membership. MRM is a synthesis of two popular statistical models, the simple Rasch model from IRT family and an unconstrained latent class model (LCM; Dayton, 1999). In particular, the author address recovery of MRM structure and parameters in the presence of an examinee covariate using Markov chain Monte Carlo (MCMC) estimation (Gilks, Richardson, & Spiegelhalter, 1996).
This study is an example of the growing interest in mixture models more generally in the statistical literature and in psychometrics more specifically. For example, Muthén (2001a, 2001b, 2004) proposed and studied a growth mixture model where the components of the mixture are growth curve models. In this model, the intercept parameter and slope parameter of the growth model are allowed to vary across latent classes of subjects. Lubke and Muthén (2005) discuss a factor mixture modelby that can be used to analyze the unobserved heterogeneity in the population. Mixture models provide flexible ways to model interesting and often theoretically important patterns in psychometric data. Mixture models are also notoriously difficult to estimate (Everitt & Hand, 1981).
Literature Review and Rationale
Early Development in MRM
The Rasch model is a central component of the proposed MRM with covariate. The Danish statistician Georg Rasch introduced this IRT model in the early 1960’s (Rasch, 1960, 1961). Three key articles in the 1990’s contributed to initial developments of MRM, namely Rost (1990), Kelderman and Macready (1990), and Mislevy and Verhelst (1990). These articles studied MRM from different perspectives, with varying formulations of the response functions and assumptions about underlying structures. Rost (1990) proposed a mixed Rasch model that combined Rasch model with LCM. The aim of integrating these two modeling techniques is to combine theoretical strengths from both approaches. Rost’s (1990) formulation is the initial building block of MRM with covariate for this current study. Kelderman and Macready used a conceptually equivalent but alternative representation, the loglinear LCM, to analyze interaction effects between grouping variables for examinees (either manifest or latent) and item parameters.
Mislevy and Verhelst (1990) proposed a mixture linear logistic test model (LLTM: Fischer, 1973) to explore the prior assumption that different strategies were adopted by different groups of examinees. In applying this particular mixed model, hypotheses about relationships between item features and item functioning in different classes were made, and item parameters and numbers of latent classes were restricted accordingly. This study presumes a general IRT framework without prespecifications of item features and hypothesized response patterns, as in studies of Rost, and of Kelderman and Macready. More recently, DeAyala, Kim, Stapleton, and Dayton (2002) supported the mixture distribution conceptualization to study DIF. They conducted simulation studies to compare various DIF detection approaches. Maij-de Meij, Kelderman, and van der Flier (2010) found that mixture IRT approach outperformed the manifest DIF approach in identifying DIF items. Their Monte Carlo study analyzed various factors in locating DIF items.
The Rationale to Include Auxiliary Examinee Variables
As Embretson (2006) pointed out, a disadvantage of using mixture IRT model was difficulty in interpreting the qualitative meaning of the latent groups identified by models. Smit, Kelderman, and Flier (1999) carried out a simulation study of MRM including collateral information from examinees. These author used MRM to identify homogeneous Rasch scalable groups in social survey data. They found that standard errors, as well as assignment of examinees to latent groups, benefit substantially by including external variables related to latent groups. In more recent research, Cohen and Bolt (2005) applied MRM to identify items that function differently across both manifest and latent groups of examinees. Cohen and Bolt included two studies. In the first study, Cohen and Bolt compared DIF items identified through traditional manifest groups approach and latent group approach. In the second study, they used a two-step procedure to identify DIF items in a college mathematics exam; they first estimated item parameters in latent class mixture model without the covariate and then used these estimated item parameters as known for modeling associations between resulting class-membership estimates and observed covariates. A limitation of their study was the use of a two-step process to detect DIF. It is an “after the fact” analysis, in contrast to directly including a covariate in the MRM as proposed in the present study.
Samuelsen (2005) incorporated covariate in MRM by creating dichotomous variables as indicators of memberships for manifest groups and estimated proportions of every manifest group in each latent class individually and separately. She found better recovery of item parameters with stronger relationship between the covariate and latent classes, especially for smaller sample size. Cho, Cohen, and Kim (2006) analyzed the effects of different kinds of prior information for detecting latent classes in MRM, also with a focus on including a covariate. The effect of a covariate on estimations was touched on only briefly. In this proposed study, the author systematically study influences of a covariate in an MRM for estimating latent classes and item parameters.
The Contribution of This Study
In this study, MRM is used to simulate data sets in which test items perform differently across latent classes, in the presence of a person covariate. The study extends work done in the studies reviewed earlier. This research builds upon several previous research studies, and differs from them in the sense of the manipulated factors in Monte Carlo study and direct inclusion of a covariate to estimate mixing proportions. The author anticipate that the inclusion of a covariate will help resolve challenges in estimating IRT mixture models.
In addition, performances of MRM with a covariate and MRM without a covariate are compared so that the results can examine advantages and possible drawbacks of incorporating collateral information.
With a careful design of this simulation study, a contribution can be made toward understanding how MRM performs under different simulated conditions, and that caveats for using this complex model can be provided for other researchers, especially those who are interested in applying MRM for DIF detection under a latent class framework. This study models the recovery of simulated structures as a function of various effect sizes on item difficulty across latent groups, and explore the effects of connections between a covariate and latent classes. Extreme simulation conditions are small DIF effect sizes, weak links between distributions of covariate, and unbalanced distributions of latent group membership. In addition, the study provides advice for dealing with the so-called label-switching problem associated with LCMs, in the particular form of mixture modeling under MCMC estimation. (Label-switching refers to the fact that for a given set of values for class memberships, proportions, and conditional probabilities, equivalent solutions can be obtained by relabeling solutions so that all parameters previously associated with Class k are now associated with Class k′.)
Method and Simulation Study Design
This section describes the model, lists manipulated factors, describes simulation conditions, and explains the decisions made in analyzing simulation results.
MRM With Covariates
Latent class membership is modeled here with a logistic model with covariates as predictors. This flexible approach makes it possible to model relations between latent groups and any number or form of manifest covariates.
The probability function for mixed Rasch model of responses of examinee j to I items, which additionally incorporates the indicator ς jg that takes the value 1 if examinee j is in group g and 0 otherwise is
where xji is the response of examinee j to item i, b ig is the difficulty of item i in class g, and latent traits are θ jg for examinee j with respect to group g (although each examinee is actually in only one group).
Suppose each class is normally distributed with mean μ
g
and variance
MRM with covariate has an extra layer of modeling for π g conditional on an examinee’s covariate. In this study, π g is modeled by a logistic regression function with the covariate as the predictor. In other words, ς jg ~ Categorical(π1, . . . , π G |yj), which will be written as ς jg ~ Categorical(πj1, . . . , π jG ).
As previously specified, this study addresses a two-class solution. Thus, in the logit model for mixing probabilities, the dependent variable π
jg
will be the probabilities for two outcome categories (i.e., g = 1 or 2 for examinee j). Specifically,
or, equivalently,
These expressions specified in Equation 4 are not identified without further restrictions. The author fix
I will focus on a covariate that is a dichotomous variable with values “0” and “1.” To be more specific, the author lay out equations for the two possible situations of an examinee j; these are (a) when his or her corresponding value in
It follows that the odds ratio between different covariate groups and latent classes is
When the probabilities of latent class membership are the same for both values of y, the covariate is independent of class membership and exp(β12) = 0. Stronger relationships are indicated by greater positive or negative values.
Primary Design of the Simulation Study
This section focuses on factors manipulated or kept constant in the simulation study. An important contribution of studies like this is to provide suggestions and cautions for researchers who apply theoretical models to analyze practical problems. To serve this purpose, the author chose levels for each manipulated factor consistent with a real-world test scenario or similar to relevant research studies. Table 1 summarizes manipulated factors and their corresponding levels.
Factors Varied in This Simulation Study.
Note: LC = latent class; DIF = differential item functioning.
Table 1 shows that there are two conditions for participants’ latent abilities. m0 standards for simulation cells in which latent abilities for both groups were simulated from Normal(0,1). m1 stands for simulation cells in which participants’ abilities were simulated from Normal(0,1) for one latent group and Normal(1,1) for the other group. There were three variations of latent class distributions: 15%:85%, 30%:70%, and 50%:50%, which are abbreviated LC15, LC30, and LC50 respectively. There were two different proportions of large DIF items.
Item parameters for one latent class were simulated from the range of −2 to +2. The sum of item parameters was adjusted to be 0. Then DIF values were added to either 6 (20%) or 12 (40%) items to create item parameters for a second latent class. Item parameters for the rest of the items for the second latent class were adjusted so that the sum of item parameters was equal to 0 as well. After summing item parameters to be zero for item parameters for a latent class, the magnitude of DIF for 40% condition was equal to 0.8 while the magnitude of DIF for 20% was equal to 1.0. The magnitude of DIF parameters before adjusting item parameters (i.e., adjust item parameters so that the sum of all item parameters is zero for a latent class), is equal to 1.3 for items with large DIF or equal to 0.3 for items with small DIF.
To express the effects of the covariate on the latent group classifications, I use the odds ratio defined above. This single index demonstrates the magnitude of relationship between the observed covariate and latent classes. Levels of exp(slope[2]) were set equal to 1, 10 and 50—that is, log-odds of 0, 2.3, and 3.9—which correspond to no, medium, and strong relationships between covariate and latent groups.
In all conditions in this simulation study, two latent classes are specified, one covariate, and 30 items and 1,000 participants.
For each cell in the simulation design, one extra replication was added to the 10 replications that were estimated using MRM. This extra replication was estimated using simple Rasch model, which is the one-class version of MRM, on response data generated from MRM specified in that cell. In other words, fit of the simple Rasch model and its single set of item difficulties are obtained based on data from the true MRM for that cell. The reason for this extra replication is that in the preliminary study, there were occasions where estimates of mixing proportions yield a very large and a very small latent group with estimated parameters far away from the generating ones. Usually, for this type of result, the estimated item difficulty parameters for the large latent classes were approximately the average of the generating item difficulties across both latent classes, while the estimated item difficulty parameters for the small classes were a random combination of two groups of generated item parameters. That is, the fitted two-class model in effect collapsed to a single class. Thus, for each cell in the design, if this type of result appears, I use this special additional replication to identify the occurrence of this situation.
In summary, the author manipulated five simulation factors and kept other four factors constant. The variable factors produce 3 × 2 × 3 × 2 × 2 = 72 conditions in total. With 11 replications for each condition in the design, there are 792 simulations runs in total.
Regarding prior distributions in the MCMC analyses, the author chose normal distribution for an average of latent abilities and gamma distribution for the precision of latent abilities. The prior for item parameters were normal distribution also. And for participants’ latent group classification, they used categorical distribution as prior.
Estimation
Fitting the MRM was accomplished using Markov chain Monte Carlo (MCMC) estimation for a full Bayesian model (Gelman, Carlin, Stern, & Rubin, 1994). This section first indicates the prior distributions used for the variables in the study. It then briefly describes MCMC estimation, highlighting the ideas of stochastic convergence and multiple chains.
A full Bayesian model requires the specification of the full joint distribution among all variables in the problem, both parameters and data. When data are obtained, inference is obtained in the form of posterior distributions from any variables which have not be observed, namely the parameters and missing values. Specifying the full model is accomplished from the inside out: A specification of the probability function for the observed data in terms of functional form and parameters, then prior distributions for the parameters, and when needed, prior distributions for parameters of the prior distributions.
In this study, the probability function for the data is the MRM, specified in Equation 2. Distributions for the other variables in the model are as follows:
Normal distributions for distributions of θ within classes, as is typical in IRT practice: i.e.,
Relatively mild priors for the parameters of the θ distributions:
For item parameters, the scale was set in each class by constraining the sum of the b
ig
s to be zero, via
As noted previously,
MCMC estimation of a full Bayesian model consists of producing a sequence, or chain, of draws for the parameters from functions that have the property that once the chain attains stationerity (i.e., it converges stochastically), a draw for the parameter is distributed in accordance with its marginal posterior distribution. The interested reader is referred to Gelman et al. (2004) for details. Certain properties of this procedure, as it is implemented in the WinBUGS computer program (Spiegelhalter, Thomas, & Best, 2000), are of particular importance to this study. The particular MCMC approach used in WinBUGS is called Gibbs sampling, where in each cycle, a draw is taken for each variable from the so-called full conditional distribution—its distribution conditional on the observed data and the value of every other unobserved variable from its preceding draw.
Under broad conditions, such a chain eventually converges stochastically to the joint posterior of the variables. In practice, determining convergence is not always easy. The chains are designed to wander around the posterior in proportion to the posterior densities but it is hard to tell simply by looking at the chains whether this has occurred. Under favorable conditions, a chain reaches a point where it appears to moving randomly but regularly through the parameter space, not “getting stuck” in some regions for long and irregular stretches. The former condition is known as “good mixing,” the latter as “bad mixing.” Bad mixing is more likely to occur in problems with models for which parameters are poorly determined, and models that accord poorly to the model being fit. To assess stochastic convergence, Gelman et al. (2004) recommend running multiple chains, starting from different sets of initial values. After running the chains from these initial values, one looks for when (if ever) the multiple chains for each variable appear to be covering the same region of the parameter space with similar densities. The Brooks-Gelman-Rubin diagnostic in WinBUGS checks whether total variance for a stretch of cycles is similar to within-chain variance, a sign of likely stochastic convergence. Sometimes, however, the chains are visually distinct.
One of the challenges in estimating mixtures is the phenomenon of label-switching mentioned earlier. It can be the case that different chains have converged to equivalent but differently-labeled solutions. If they are stochastically equivalent after relabeling, I can say the estimation process has converged. Other times nonequivalent chains show qualitatively different characteristics, in this study for example with one chain providing a stochastically converged solution for the MRM but another badly-mixed chain or a collapsed solution.
Evaluation of Estimation Outcomes
Based on findings from preliminary analyses, an initial step for each simulation run was the classification of results into different categories, for the reasons discussed in the previous section. Two kinds of analyses were carried out. One addressed percentages with which analyses in each cell recovered the generating MRM structure and design features that influenced recovery rates. The other addressed only those runs in which the correct MRM structure was recovered and in those cases the author examined the bias and root mean square error (RMSE) of parameter estimates.
Convergence and Structure Recovery
Based on experience from preliminary studies, the number of burn-in cycles for each simulation run was 4000 iterations. An additional 6000 iterations I run to achieve convergence if possible, a number that was typically more than adequate. To identify converged or non-converged MCMC chains, WinBUGS provides graphic tools to visually identify whether multiple MCMC chains converged to a single solution or multiple solutions. This function is accomplished in part by looking at history graph, which plots values for requested parameters across all MCMC iterations. I use mixing proportions as a first category of parameters to check for convergence, since the occurrence of label-switching is intuitively obvious for these parameters; for example, one chain with posterior means of 0.24 and 0.76 for π1 and π2, the other with posterior means of .76 and .24.
In this study, I consider a replication as a properly converged if the two requested MCMC chains, which start with different random values, merge into the same stable region and provide estimates with reasonable values. Some runs converged, but collapsed to essentially a single-class solution—either proportions of class membership close to 1 and 0, or distinguishable class proportions, but essentially the same item parameters for both. First, the author report and analyze proportions of replications in each cell with recovered structures. Second, they analyze parameter recovery in solutions with the correct structure.
Bias
To compare parameters from simulated data and recovered data in those solutions that recovered the generating two-class structure, the author used bias and RMSE. Bias is the average difference between estimated value and true value. RMSE is the square root of average squared differences between estimated and true values. For each cell, bias and RMSE are computed across converged replications. The author are especially interested in checking and reporting two groups of estimated parameters against the generating parameters. One is the estimated proportions of latent classes and proportions of manifest groups for covariates. The other group is items specified to be significantly different across two latent classes (i.e., the large DIF items).
Computer Programs
In this study, the author first used SAS (version 9.0; SAS Institute, 2002) to generate each dataset based on the stated MRM distributions and corresponding true parameters. Then they used WinBUGS 1.4 (Spiegelhalter, Thomas, & Best, 2000) to obtain estimated parameters for the item difficulties, the mixing proportions in logistic regression functions, and examinees’ group memberships. Finally, they extracted parameter estimates and compared them with true values to decide whether a simulation run was converged and had recovered correct structure; if so, they reported the accuracy of estimated parameters.
Summary of Results
In this study, the author generated data, estimated parameters, and summarized Bayesian MCMC estimation of simulated data. There were replications that accurately identified the underlying DIF structures while other replications were unable to recover generating DIF structures. Another phenomenon was the stability of MCMC estimation. As described in the previous section, the MCMC chains might be converged or non-converged after 4000 burn-in cycles and 6000 additional cycles. Based on these last 6000 MCMC cycles, I decided whether the two requested MCMC chains converged, and if they converged, whether one or both recovered the two-class solution.
The complexity in the simulation design offers rich information for researchers for understanding different possible results from estimation of a mixed Rasch model with a covariate, especially under known but extreme simulation conditions. On the other hand, the complexity of results added complication to identifying and categorizing each replication into different outcome groups. The author analyzed the results in two steps:
Categorized outputs into qualitative outcome groups that reflected whether a each estimation run converged and, if so, whether it recovered the generating two-class latent structures and DIF structures.
Among those simulation runs that did recover true generating structures, presented summary statistics for evaluating parameter estimations.
In both steps, the author studied the influence of the simulation factors, as independent variables for respective outcome analyses: category proportions of recovery in the first case and quantitative parameter recovery statistics in the second case.
The result classification process was multi-faceted, in which the author had to determine cut-off points at various review steps. A detailed description of the process appears in Dai (2009). Here I summarize the results of the classification process and parameter-recovery statistics from solutions that recovered generating two-class structure.
Solution to Label-Switching Runs
In this simulation study, the solution used to deal with label-switching problem was to relabel estimated parameters from each latent class in those replications where between-chain label-switching occurred. This issue was inevitable in both Bayesian estimation and maximum likelihood (ML) estimation: In ML estimation, it was manifest as multiple maxima that were identical except for labeling, and in MCMC estimation, it was manifested as between-chain or within-chain label-switching. Between-chain switching means that labeling of latent classes differs among multiple requested MCMC chains. In contrast, within-chain switching means that labeling of latent classes switches between iterations for a single MCMC chain.
In the literature, researchers have suggested several solutions to the label-switching problem. There are artificial identification constraints (e.g., Diebolt & Robert 1994), relabeling algorithms to perform a k-means type clustering of the MCMC samples (Stephens 1997; Celeux 1998), label invariant loss functions (Celeux, Hurn & Robert 2000; Hurn et al. 2003), and random permutation samplers (Fruhwirth-Schnatter, 2001). Jasra et al. (2005) provide a detailed review of available solutions.
In the preliminary analysis, the author examined the simple approach recommended by Chung, Loken, and Schafer (2004). This method requires viewing results of multiple MCMC chains in WinBUGS output, and if label switching is observed, fixing class membership of one case from each class as known in accordance with one labeling, for cases with very high posterior probabilities of belong to that class. However, this solution was not sufficient to resolve label-switching problem in the proposed MRM with a covariate model.
Recovery of the Underlying Structure
Based on five-step result reviewing process (Dai, 2009), simulation solutions were classified into categories shown in Table 2, which details the denotations and descriptions of the result categories, result types, and simulation factors. Table 3 shows corresponding decision for each result category in three specified results types:
Denotation of Result Categories, Result Types, and Simulation Factors.
Result Types With Corresponding Result Categories and Result Summary Decision.
Note: PM = poor mixing; WCLS = within chain label switching; MCMC = Markov chain Monte Carlo; C = collapsed; LSR = latent structure recovered; BCLS = between-chain label switching.
Two MCMC chains converged to a single solution: 48% of 720 simulation replications;
Two MCMC chains did not converge to a single solution and the result category for both MCMC chains is the same (includes same solution with different labeling): 37% of 720 simulation replications;
Two MCMC chains did not converge to a single solution and the result category for each MCMC chains is different: 16% of 720 simulation replications.
Table 3 shows the classification of results for each simulation cell and frequency count for each result category. Values in Table 3 are the number of replications for each simulation cell belonging to each result category per major result type. Recall that there were ten replications of generated data for each of 72 simulation cells, and for each replication, two MCMC chains were requested for each WinBUGS run. Thus there were 720 simulation replications and in total 1440 MCMC chains for the entire simulation study.
Overall, 80% of the MCMC runs recovered the underlying latent class structures and DIF structure correctly. About 20% of simulation outputs did not recover underlying latent classes and DIF structure well. For those simulation outputs where the underlying structure was not recovered, 17.29% of total simulation results belonged to the category of “collapsed solution.” This category is essentially a single-class solution. The detailed descriptions of result categories are presented in Dai (2009). In general, for simulation cells with combinations of extreme simulation conditions, within-chain label switching or “bad mixing occurred.” (When MCMC chains fluctuate nonsystematically within wide range of possible values for the parameter, it is called “poor mixing” in MCMC literature.)
As shown in Table 4, these extreme simulation conditions also lead to a low recovery rate of latent classes and DIF structure. In Table 4, 10 out of 11 simulation cells with less than 50% recovery rate are those simulation cells where there was no relation between the covariate and latent classes (i.e., odds ratio equal to 1).
List of Simulation Cells With Recovery Rate Less Than 50%.
Note: LC = latent class; dif = differential item functioning.
A regression analysis was used to model effects from simulation conditions on the percentage of successful replications for each simulation cell. In this analysis, all simulation factors were included so I can compare relative effects of the different factors. Table 5 shows how simulation factors are coded into dichotomous variables for regression analysis. The arcsine transformation of the proportion of recovered runs is used as the dependent variable to ensure the homogeneity of variance assumption in multiple regression analysis of the dichotomous dependent variable, successful structure recovery.
Denotation of Simulation Conditions in the Regression Analysis.
Table 6 presents results from the regression analysis of recovery rates. All simulation factors except distribution of covariate groups, labeled ‘C_COV’, have statistically significant effects on the proportion of recovered replications. The findings are summarized below:
More replications recovered the underlying structure when more items have large DIF effects. Keep in mind that the total number of items is equal to 30.
When latent class proportions are 30%:70% or 50%:50%, the proportion of recovered replications is significantly greater than under the 15%:85% condition. (The sample size is 1000 for all simulation replications, so the small sample as well as the lower proportion of the 15% class may both play roles. The study was not designed to disentangle these effects.)
When the log odds ratio between covariate groups and latent class proportion is equal to 50, which is the strongest connection between covariate groups and latent groups, the effect on percentage of successful replication was twice the effect from the condition when the odds ratio is equal to 10.
For those cells with different average latent ability, the percentage of successful replication is higher than those cells with equal average latent ability for both latent classes. However, this effect was smaller when the log odds ratio is equal to 50.
Regression Analysis of Recovered Replications.
The regression coefficient is statistical significant at α = .05 level.
Interactions of different pairs of simulation factors were tested. The interaction between the odds ratio of 50 and different average latent ability was found to be statistically significant. Based on the estimated coefficient in Table 6, when there is a strong relation between the covariate and latent class distribution, the recovery rate decreases when average latent abilities for two latent groups of examinees are different, compared with those for which average latent abilities for both groups are the same.
Identification of DIF Items
The author calculated correlations between estimated and generating values for all DIF parameters. Figure 1 is a histogram of calculated Pearson correlations. In both sides of Figure 1, the correlation is very high, close to −1.0 or +1.0 (equivalent solution, different labeling). There are high correlations, medium-high correlations, low correlations, and very low correlations. As Figure 1 shows, the majority of correlation between simulated and estimated DIF parameters are above 0.70.

Histogram of correlation index.
In addition, the author calculated percentages of statistically significant DIF parameters that were generated with large DIF effect size. The logic behind this step was that a replication with a recovered DIF structure should identify most large DIF items as statistically significant. Figure 2 is a histogram of percentage of statistically significant DIF parameters out of total number of large DIF parameters. Figure 2 shows that the majority of simulation cells identified at least 60% of large DIF items.

Histogram of correlation coefficients.
Bias and RMSE: Regression Analysis
After results were categorized, those that exhibited successful MRM structure recovery were included in the calculation of average bias and RMSE for each simulation cell. There were five groups of independent variables involved in modeling process. The outcome parameters included mixing proportions, item difficulty, DIF parameters, and mean and variance of latent ability distributions. Regression analysis was used to show relations between simulation factors with average bias and RMSE as the dependent variables. Obtaining non-significant interaction terms as in the convergence terms in the convergence rate analysis, the regression analysis in this section also focuses on main effects.
Table 7 contains the results of regression analysis related to mixing proportions. The left panel shows that the bias in the estimates of mixing proportions was affected by proportions of large DIF effect items and distributions of latent groups. Though the effect was statistically significant, the magnitude of the effect of having more large DIF items was small. However, when the sizes of the latent classes were equal, the bias for the larger latent class increased significantly compared with the cases in which sizes of the latent classes were different. Regarding RMSE, as shown in right panel, equal size latent classes lead to smaller RMSE in the estimation of the mixing proportions.
Regression Analysis of Bias for the Mixing Proportion P.tot[1].
Note: P.tot[1] = mixing proportion for Latent Class 1.
The regression coefficient is statistical significant at α = .05 level.
Table 8 shows the regression analysis results for intercept and slope parameters in logistic link function. The bias of the intercept decreased when there were more large DIF effect items and when mean of the latent ability was distributed differently across latent classes. For the slope parameter, bias increased when the distribution of the latent classes was 30%:70%. When the relation between the covariate and the latent classes was moderate or strong, the bias of slopes decreased.
Regression Analysis of Bias for the Intercept and Slope in Logit Link.
The regression coefficient is statistical significant at α = .05 level.
As shown in Table 9, the RMSE of the intercept estimates was influenced by the proportion of large DIF items, the relative sizes of the covariate groups, the distributions of covariate groups among each latent class, and whether the average latent ability differed between latent classes. RMSE deceased significantly when 40% of the items were large DIF items comparing to the 20% case. When the sizes of the covariate groups were equal, the RMSE of the intercept also decreased. When there was a relationship between the covariate and latent classes, the RMSE of the intercepts were smaller. Different average latent abilities across latent classes lead to smaller RMSE for the intercept, too.
Regression Analysis of RMSE for the Intercept and Slope in Logit Link.
Note: RMSE = root mean square error.
The regression coefficient is statistical significant at α = .05 level.
For the slope parameter in the logistic regression link function, the magnitude of connection between covariate groups and latent groups, percentage of large DIF items, distributions of latent classes, and unequal average latent abilities affected the RMSE of the slope significantly. As expected, the RMSE decreased when there were more large DIF items. Comparing the latent class distribution of 15%:85% to distributions of 30%:70% or 50%:50%, the RMSE of the slope parameter was smaller. The RMSE for the slope parameter estimates increased when there was a strong relation between covariate groups and latent classes. Different average latent abilities across latent groups resulted in smaller RMSE for the slope estimates.
The last group of variables of interest were the item difficulty parameters for latent class one and latent class two. As described in previous section, the sum of item parameters was set equal to zero in each class for parameter identification. Thus, only the means of RMSE (but not bias) from item parameters for each latent class were used as dependent variables in regression analysis. Table 10 shows that average RMSEs of item parameters from latent class one (the smaller size latent class) were affected by latent class proportions, the relation between the covariate and the latent classes, and whether average latent abilities were the same across latent groups. All three factors contributed to lower RMSE for item difficulty parameters from latent class one.
Regression Analysis of RMSE for Average Item Difficulty Parameters in Each Latent Class.
Note: RMSE = root mean square error.
The regression coefficient is statistical significant at α = .05 level.
For item difficulty parameters from the larger latent class, proportions of large DIF items and latent class proportions had significant effects. The closer to equal size of latent groups, the higher were the RMSE of item parameters from latent class two. This effect was the opposite of the effect of latent class proportions on RMSE for item parameter from latent class one.
As previously described in the simulation design section, certain conditions were constant for all replications. Among them are the total sample size (N=1000) and total number of items (N=30). Under these constraints, the effects from corresponding simulation factors (latent class proportions and percent of large DIF) may also related to the number of persons in each latent class and number of large DIF items respectively. Future studies with varying sample sizes and numbers of items will help differentiate interrelations between proportions and numbers of persons and items.
Summary of Results
In result summary analysis, I used as dependent variables (1) the proportion of replications in which generating mixture structure was recovered, and (2) the bias and RMSE of parameter estimates among those cases in which the generating two-class structure was recovered. In both cases I conducted regression analyses to examine the effects of the simulation factors.
In modeling of recovery rate for each simulation cell, four simulation factors had statistically significant effects on recovery of underlying structure. The only non-significant effect was the marginal distribution of covariate group.
To be more specific, the higher the percentage of items with large DIF effect size, the higher the recovery rate of underlying DIF and two-class structure. Compared with DIF% = 20%, simulation cells with DIF% = 40% yielded higher correct recovery rates. For latent class distributions, when two latent classes were distributed as 30%:70% or 50%:50%, recovery rates from these two levels of this simulation factor were always higher than was the recovery rate from two-class structure with 15%:85% distribution. This indicates that when the sizes of the two latent classes were quite different from each other, it was harder to obtain stable and accurate estimation after controlling for other simulation factors.
As the focus of this study, the inclusion of collateral information in MRM had positive effects on correct identification of latent structure. Compared with MRM without a covariate, that is, the condition where the odds ratio between manifest groups and latent groups was one (exp[slope(2)]=1), MRM with the covariate moderately related (exp[slope(2)]=10) or strongly related (exp[slope(2)]=50) with latent class membership always achieved higher recovery rates. As expected, the stronger the connection between manifest group distributions and latent class distributions, the higher the chances of correctly identifying the latent structure.
In addition, when the examinees’ latent ability distribution was different across two latent classes, there were higher chances to obtain stable and reliable estimations of parameters in MRM with a covariate. This was consistent with the expectation that the more distinct the latent classes were, the easier the MRM model with a covariate could separate examinees into different latent groups.
Average bias ranged from −0.20 to +0.20. For RMSE, the range was from 0.15 to 0.25. In the regression analyses of bias and RMSE for mixing proportions, I found that when the percentage of items with large DIF effect size was high and when the sizes of latent classes were equal across groups, the effects of these two factors on the bias of mixing proportions were statistically significantly. When latent classes were distributed as 50%:50%, the RMSE of mixing proportions was significantly lower than that from other simulation levels and simulation factors.
The bias of intercepts in the logistic regression link function decreased when the percentage of items with large DIF was high and when latent ability distributions were different across latent groups. In contrast, for the slope parameter, the bias decreased when odds ratio between latent groups was high and increased when latent class distribution was equal to 30%:70%, compared to the other levels of these factors.
Discussion
In this study, a latent class modeling technique was combined with the Rasch model from item response theory. On one hand, by combining these two advanced statistical methods, examinees’ response data can be analyzed with respect to more complex underlying structures. On the other hand, all other things being equal, the more complicated a statistical model is, the larger the sample size is required to obtain stable estimation. The strategy was to hold sample size constant, and examine effects of these “other things” that might not be equal, but can have substantial impact on estimation.
In particular, additional information about examinees from covariates can mitigate estimation difficulties when researchers use complicated statistical models to analyze real data. The author show how a covariate can help resolve problems associated with estimating IRT mixture models. Based on the findings in this research, additional information in the form of a covariate can increase chances of recovering correct mixture structures and improve estimations of parameters within that structure.
Insights on Model Estimation
Based on our experience in this simulation study the author offer some recommendations for practitioners who are interested in applying the mixture IRT model with a covariate in simulation studies. First, it is very important to also use a single-population IRT model to estimate item parameters. The reason for doing this is to make sure that when “collapsed” run occurs using MRM, researchers can compare item parameters from MRM and item parameters from an IRT model where item parameters are constrained to be the same for all examinees. In this study, about 20% of simulation solutions were identified as “collapsed” runs using known generating parameters and item parameters estimated from simple Rasch model based on the same simulation factors.
To identify label-switching issues, the author would recommend that researchers look at detailed history plots of mixing proportions. In this study, the initial step in result reviewing process was to look at history chart of mixing proportion from WinBUGS. The reason for doing this was that most information about this parameter was available in the WinBUGS history chart. It served as a good indicator of whether label-switching occurred and, if so, what type of label-switching it was. In this study, to resolve label-switching, I relabeled estimated parameters for between-chain label-switching runs and excluded within-chain label-switching runs because there were very few of them. As noted in the Literature Review and Rationale section and discussed in Methodology and Simulation Study Design section, there were other solutions to deal with label-switching in finite mixture models. I explored some of these solutions in my preliminary analysis and chose the relabeling method for this study because this approach was much more resource-efficient in solving between-chain label-switching in the proposed model. This advice holds both for simulation studies like the current one and for applied analyses of real data.
In simulation studies, after an initial review based on mixing proportions it is next useful to calculate correlation indexes between estimated DIF parameters and generated DIF parameters. This is because comparison between the estimated and true DIF structure provides evidence about whether the underlying structure is correctly identified. Also, calculating the percentage of statistically significant DIF is meaningful in that this allows the comparison of DIF directions from estimated DIF parameters and generated ones.
Furthermore, practitioners should be cautious when estimated mixing proportions are close to zero or one. Because of the complexity of MRM, it is quite possible that the information contained in data is not sufficient to separate examinees into distinct latent classes or the sample size is not large enough to apply this model. Thus, when the true model has a very large latent class and a small one, there is a good chance that the estimation will result in ‘collapsed’ runs.
Future Research
Four simulation factors were constant in this study:
The number of latent classes was two;
The number of covariates was one;
The test length was 30;
The total number of examinees in simulated response data was 1000.
In future research studies, it would be of interest to vary these simulation factors so that researchers can study their effects on the estimation of parameters and recovery of the underlying structure.
Another direction for future research is including more than one covariate in the logistic regression function. Since this simulation is exploratory in its nature, future research can build more complicated models based on the results of this study, and having more than one covariate is an important research direction. Also in this research project, I only studied situations where the number of latent classes was two. Increasing to more than two latent classes would be another useful extension. The form of multinomial logistic multiple regression is well known, and lends itself naturally as a link function in mixture IRT models with multiple classes and multiple examinee covariates.
Recovering underlying latent structure was a main focus for this study. Another focus of studying mixture IRT models can be the classification of examinees into latent classes. In estimating MRM, it is possible that underlying latent class structure and DIF structure are well recovered, but classification of examinees into generating latent class is noisy. This can be another research direction for analysts who are interested in studying MRM. A simulation study conducted by Lu and Jiao (2009) focuses on this line of research in MRM without a covariate.
In this study, the covariate was directly incorporated in MRM and a one-step process was thus used to link the covariate with examinees’ latent class memberships. This was how this study differs from the Cohen and Bolt (2005) study which used a two-step procedure to locate the cause of DIF. Thus, another research direction would be the comparison of the performance of MRM between the one-step process and the two-step process. From the author’ experience in conducting this simulation research (and as also suggested in Samuelsen, 2008), the one-step process was more difficult to execute. The author would expect that, under certain simulation conditions, the two-step process (i.e., obtain item parameters and examinees’ latent class memberships before modeling collateral information with latent class classification) to include information from a covariate might lead to more stable but less efficient parameter estimation.
Footnotes
Acknowledgements
The author wishes to thank Dr. Robert Mislevy for his advice and contributions to this article, and the editor, Dr. Dan Bolt, and two anonymous reviewers for their insightful suggestions on an earlier draft. I am grateful to the Association of American Medical Colleges (AAMC) for their summer graduate student research program, in which the author initially studied mixture IRT models.
Declaration of Conflicting Interests
The author declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
Funding
The author received no financial support for the research, authorship, and/or publication of this article.
