{'fraction_non_alphanumeric': 0.03870445344129555, 'fraction_numerical': 0.015991902834008098, 'mean_word_length': 4.91782462865357, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
A Multidimensional Hierarchical Framework for Modeling Speed and Ability in Computer-based Multidimensional Tests

In psychological and educational computer-based multidimensional tests, latent speed, a rate of the amount of labor performed on the items with respect to time, may also be multidimensional. To capture the multidimensionality of latent speed, this study firstly proposed a multidimensional log-normal response time (RT) model to consider the potential multidimensional latent speed. Further, to simultaneously take into account the response accuracy (RA) and RTs in multidimensional tests, a multidimensional hierarchical modeling framework was proposed. The framework is an extension of the van der Linden (2007; doi:10.1007/s11336-006-1478-z) and allows a"plug-and-play approach"with alternative choices of multidimensional models for RA and RT. The model parameters within the framework were estimated using the Bayesian Markov chain Monte Carlo method. The 2012 Program for International Student Assessment computer-based mathematics data were analyzed first to illustrate the implications and applications of the proposed models. The results indicated that it is appropriate to simultaneously consider the multidimensionality of latent speed and latent ability for multidimensional tests. A brief simulation study was conducted to evaluate the parameter recovery of the proposed model and the consequences of ignoring the multidimensionality of latent speed.

In reality, respondents are likely to bring multiple latent abilities to bear when responding to items; meanwhile, items are likely to require various latent abilities to determine a correct response, especially in multidimensional tests (e.g., Reckase, 2009;Whitely, 1980;Tatsuoka, 1983). In addition, with the increasing demand for providing more detailed and refined feedback to test-takers, multidimensional tests have received much attention from practitioners and researchers.
In psychological and educational measurements, an appropriate notion of latent speed on test items is that of speed of labor (van der Linden, 2009). Therefore, latent speed can be defined as a rate of the amount of labor performed on the items with respect to time (van der Linden, 2011).
Actually, the definition of latent speed should be discussed in a certain dimension of latent ability, because using the required latent ability is the basis of an effective labor. Due to the multidimensional nature of latent abilities and different items may require different kinds of abilities, latent speed may also be a multidimensional concept, each dimension of which corresponds to a specific type of labor or latent ability. For example, the latent speed corresponding to the latent ability of decoding (or algebra problem solving) may be different from the latent speed corresponding to the latent ability of encoding (or geometry problem solving). For another example, when respondents, especially for non-native English speakers, take part in the GRE ® Subject Test (e.g., Mathematics), at least two abilities are needed, one for understanding the questions (e.g., the English reading ability), and one for solving the questions (e.g., the subject ability). Meanwhile, corresponding two latent speeds worked, one reflects the working speed of reading, and one reflects the working speed of problem-solving or applying subject ability.
Currently, although multidimensional models for response accuracy have been well developed (see Reckase, 2009), to our knowledge, there is a lack of multidimensional models for RTs to take account of the potential multidimensionality of latent speed. Recently, based on the hierarchical modeling, a few studies have attempted to use multidimensional models for response accuracy to capture the multidimensional structure of latent ability when multidimensional tests are involved; but only a URT model is used to capture the potentially multidimensional latent speed, as shown in Figure 1(b). For instance, Man, Jiao, Zhan, and Huang (2017) employed a compensatory multidimensional IRT model for response accuracy and the unidimensional log-normal RT model (van der Linden, 2006) for RTs. In addition, Zhan, Jiao et al. (2018) proposed a joint cognitive diagnosis modeling to simultaneously analysis RA and RTs in cognitive diagnosis. This approach was further extended to account for the paired local item dependence . However, in these studies, because of the lack of multidimensional RT (MRT) models, only the relationship among multiple latent abilities and one single latent speed can be evaluated. Logically, as aforementioned, different latent abilities may be associated with different latent speeds. Thus, assuming the latent speeds corresponding to different latent abilities to be identical, as done by Man et al. (2017) and Zhan, Jiao et al. (2018), may be too restrictive to describe complicated data and thus may lead to biased conclusions. It is desirable to release this limitation to allow each latent ability to be associated with its own latent speed. As URT models may be inappropriate in practice, it is critical to develop a new RT model that considers multidimensionality of latent speed.
To meet the demand, we firstly extend the most popular unidimensional log-normal RT model (ULRTM; van der Linden, 2006) to be multidimensional and call it the multidimensional log-normal RT model (MLRTM). Secondly, a multidimensional hierarchical modeling framework for modeling multidimensional latent speed and multidimensional latent ability was proposed, as shown in Figure 1

Overview of the Unidimensional Log-normal RT Model
Let Tni be the observed RT of person n (n = 1,..., N) to item i (i = 1,..., I). In the ULRTM, the logarithm function is used to transform the positively skewed distribution of RT to a symmetric shape and is assumed to be dominated by item i's time-intensity parameter ξi and person n's latent speed parameter τn as follows: or equivalently, ( where ξi represents the time needed to complete item i; τn represents the working speed of person n on a test and is assumed to be normally distributed with mean zero and variance 2 τ σ ; εni is the normally distributed residual error term, with mean zero and variance 2 ω  i ; ωi is the reciprocal of the standard deviation of the error term, which can be treated as a time-kurtosis parameter.
A basic assumption of the ULRTM is that logTnis are conditionally independent given the unidimensional τn, which is known as local RT independence . In other words, local RT independence is obtained when the relationship among items is fully characterized by the URT model. However, if the single latent speed is not sufficient to account for the relationship among RTs, an MRT model is needed.

The Multidimensional Log-normal RT Model
There are two types of multidimensional tests: between-item and within-item (Adams, Wilson, & Wang, 1997). In a between-item multidimensional test, each item measures a single dimension but different items may measure different dimensions so the multidimensionality occurs between items. In a within-item multidimensional test, an item may measure more than one dimension simultaneously. As the between-item multidimensionality can be seen as a special case of the within-item multidimensionality, to be general, we focus on the latter throughout this paper.
For an item that measures multiple dimensions simultaneously, the ULRTM can be extended to accommodate multiple latent speeds by replacing the scalar latent speed parameter, τn, with a vector of speed parameter, τn. The resulting MLRTM is defined as: or equivalently, where τnk represents the latent speed of person n in dimension k, and τn = (τn1, ..., τnk, ...,τnK)' is the multidimensional latent speed vector following a multivariate normal distribution: with mean vector μτ and variance and covariance matrix Στ, and μτ is set to 0-vector for identification; qik is an element of the confirmatory Q-matrix (Tatsuoka, 1983) indicating whether dimension k is required to answer item i correctly; qik = 1 if the dimension is required, and 0 otherwise. Other parameters are the same as those in the ULRTM.
In the MLRTM, it is assumed that logTnis are conditionally independent given τn. In addition, is set at a constant, m, all τ-vectors that satisfy the expression yield the same RT. This feature suggests that a low speed on one dimension can be compensated by a high speed on another dimension, so the proposed MLRTM is a compensatory MLRTM, which is in line with the compensatory assumption of latent abilities in compensatory multidimensional item response theory (MIRT) models. Logically speaking, within each dimension, the latent speed and the latent ability should be matched with each other (e.g., . Thus, when the latent abilities are compensatory, it is reasonable to assume that the corresponding latent speeds to be compensatory as well. On the other hand, if the latent abilities are non-compensatory, the corresponding latent speeds should be non-compensatory as well. However, the development of non-compensatory MLRTM is beyond the scope of this study, and can be studied in the future.

Multidimensional Hierarchical Modeling Framework
Since both response accuracy and RTs contain information about items and persons, it is advantageous to analysis them simultaneously. For example, one may be interested in the relationship between multidimensional latent ability and multidimensional latent speed of persons and the relationship between difficulty and time-intensity of items.
In the multidimensional hierarchical modeling framework, at the first level, an MIRT model can be used as the measurement model for response accuracy and an MRT model can be used as the measurement model for RTs, respectively; at the second level, all latent abilities and latent speeds are assumed to follow a multivariate normal distribution; meanwhile, the item RA parameters (e.g., item difficulty) and item RT parameters (e.g., item time-intensity) are assumed to follow a multivariate normal distribution. Given the "plug-and-play" nature of the multidimensional hierarchical modeling, various choices of MIRT models and MRT models can be adopted. For illustration purposes, in this study, the MLRTM is used as the measurement model for RTs, and according to the 2012 PISA mathematics assessment framework (OECD, 2013), the multidimensional Rasch model (MRM; Adams et al., 1997) was employed as the measurement model for response accuracy.

The Multidimensional Rasch Model for RA
Let Yni be the observed response for person n to item i. The MRM can be expressed as where logit(x) = log(x/(1-x)); P(Yni = 1) is the probability of a correct response by person n to item i; θnk is the latent ability of person n on dimension k, and θn = (θn1, ..., θnk, ..., θnK)' is assumed to follow a multivariate normal distribution as follows: , and setting μθ = 0 for identification; di is the intercept or easiness of item i; qik is an element of the confirmatory Q-matrix indicating whether dimension k is required to answer item i correctly.
Further, if only one dimension is assumed to be required by all items, the MRM reduces to the unidmensional Rasch model (URM, Rasch, 1960), where θn is the unidimensional latent ability of person n; di is the intercept or easiness of item i.
The URM can be identified by setting the mean of θn at zero.

Multidimensional Hierarchical Modeling
In Besides the measurement models, the multivariate normal distribution was used to describe the relationship among the latent abilities and latent speeds: where θn is a vector of multiple latent abilities of person n; τn is a vector of multiple latent speeds of person n; Σperson is a variance and covariance matrix of person parameters; K* is the total number of latent abilities and latent speeds. Take the structures in Figure 1 as examples, a bivariate normal distribution (K* = 2) can be employed for the UA-US; a fivefold-variate normal distribution (K* = 5) can be employed for the MA-US; and an eightfold-variate normal distribution (K* = 8) can be employed for the MA-MS.
For the item parameters, a bivariate normal distribution was used to describe the relationship between item difficulty and item time-intensity: where μd and μξ is the mean of item difficulty and item time-intensity, respectively; Σitem is a variance and covariance matrix of item parameter. The residual error variance, 2 ω  i , is assumed to be independently distributed (e.g., Zhan, Jiao et al., 2018), thus it is not included in Ψi.

Bayesian parameter estimation
Model  (Gelfand & Smith, 1990), whose code for the MRM-MLRTM is provided in online supplemental materials (runnable source code are also provided).
Under the assumption of local independence, Yni and logTni are independently distributed as The priors of the person parameters are set as In addition, the priors of item parameters are set as Furthermore, the hyper priors are specified as: where Ritem is a two-dimensional identity matrix. Finally, the posterior mean is treated as the estimate for model parameters.

Data Description
A real data analysis is conducted using the PISA 2012 computer-based mathematics data to explore whether the MRM-MLRTM fit the data better than the URM-ULRTM (e.g., van der Linden, 2007) and MRM-ULRTM (e.g., Man et al., 2017) when the test structure is multidimensional. It is also an example to illustrate the use of the proposed model. This dataset was used by Zhan, Jiao et al. (2018). There are N = 1581 respondents and I = 10 items. The logarithm of RTs have been computed before analysis, and all zero RTs were treated as missing data. Four dimensions that belong to the mathematical content knowledge were chosen in this study, namely, change and relationships (θ1), quantity (θ2), space and shape (θ3), and uncertainty and data (θ4). The Q-matrix is shown in Table 1. Due to the use of 10 items to measure four dimensions and the incompleteness and unidentifiability of the Q-matrix (Chiu, 2013;Xu & Zhang, 2016), it is expected that the estimation of model parameters may contain relatively high measurement errors. Nevertheless, such a test structure was retained because, the real data analysis, though not perfect, can still provide some information about the nature of the real data structure. In addition, Bayesian estimation works even when the model is not identified, and adherence to identifiability has been deemed to be largely superfluous (Gustafson et al., 2005;Muthén & Asparouhov, 2012). More details can be found in Zhan, Jiao et al. (2018).
Note that we also analyzed the RT data by using the MLRTM alone, the results indicated that the MLRTM fit the RT data better than the ULRTM. More details can be found in online supplemental materials.

Analysis and Model Selection
The three joint models in Figure 1, namely, the URM-ULRTM, the MRM-ULRTM, and the MRM-MLRTM were fit to the data. For each model, two Markov chains with random starting points were used, and each chain ran 10,000 iterations with the first 5,000 iterations in each chain as burn-in. Finally, the remaining 10,000 iterations were used for the model parameter inferences. The potential scale reduction factor (PSRF; Brooks & Gelman, 1998) was computed to assess the convergence of each parameter. A PSRF with values smaller than 1.1 or 1.2 indicates convergence (Brooks & Gelman, 1998;. Our studies indicated that the PSRF was smaller than 1.05 for all parameters, suggesting good convergence. Posterior predictive model checking (PPMC; Gelman, Carlin, Stern, Dunson, Vehtari, & Rubin, 2014) was used to evaluate model-data fit. A posterior predictive probability (ppp) value near 0.5 indicates that there are no systematic differences between the realized and predictive values, and thus an adequate fit of the model. In PPMC, the sum of the squared Pearson residuals for person n and item i (Yan, Mislevy, & Almond, 2003) was used as a discrepancy measure to evaluate the overall fit of the RA model, which is written as where P(Yni = 1) has the same definition as that in Equation (5). On the other hand, the sum of the standardized error function of logTni for person n and item i was employed as a discrepancy measure: Additionally, the Akaike Information Criterion (AIC; Akaike, 1974) (Gelman, Carlin, Stern, & Rubin, 2003), where D is the deviance, and D is the posterior mean of deviance (i.e., -2 log likelihood). Note that in Bayesian analysis, the AIC and BIC can be defined as (Congdon, 2003), where p is the number of estimated parameters. A smaller value of indicates a better model-data fit.

Results
The AIC, BIC, and DIC all identified the MRM-MLRTM as the best-fitting model and the URM-ULRTM as the worst-fitting model, as shown in Table 2. In the MRM-MLRTM, the ppp values of the RA model and the RT model were 0.789 and 0.616, respectively, which indicated an adequate model-data fit. Thus, it may be more appropriate to simultaneously consider the multidimensionality of the latent speeds and latent abilities. The MRM-MLRTM was used for further illustration.  To explore the relationship among multiple the latent abilities and latent speeds, Table 3 presents the estimated person variance and covariance matrix. Moderate to high positive correlations (0.70~0.92) were found among the multiple abilities, and lower correlations (0.560 .85) were found among the multiple speeds. However, low to moderate negative correlations (-0.60~-0.01) were found between the multidimensional abilities and speeds. Although these results are not consistent with common sense that more capable respondents tended to work faster, they are consistent with previous study findings van der Linden & Fox, 2015;Zhan, Jiao et al., 2018). As a low-stakes test, the PISA has more significant implications for countries or areas than individual respondents. A reasonable explanation could be that low ability respondents lacked motivation in taking the low-stakes test (Wise & Kong, 2005), which led to shorter RTs. Note that the variance of the first latent ability was quite large, indicating all respondents differ greatly in the first dimension.  dimension might be caused the unstable estimation based on a single item. According to the scatter diagrams, some respondents had low ability but high speed, found in the fourth quadrant, and they might demonstrate aberrant response behavior such as rapid-guessing (Fox & Marianti, 2017;Wang & Xu, 2015). Table 4 presents the estimates of the item parameters. The estimated mean item difficulty and mean time-intensity were -1.19 (SE = 0.54) and 4.29 (SE = 0.15), respectively. The estimates of the time-intensity and time-kurtosis were quite similar to those shown in Table S3 in online supplemental materials, which were estimated from the MLRTM alone. To further explore the relationship between the item intercept and time-intensity parameters, we present the estimated item variance and covariance matrix in Table 5. The correlation was -0.43, which implied that the items with a lower intercept (i.e., more difficult) might have higher time-intensity. This result was consistent with that in the literature in that the more difficult items often need more time to solve (Fox & Marianti, 2016;Meng et al., 2015;van der Linden, 2006;2007).

Design and Data Generation
A real data study has been provided to demonstrate the applications of the proposed joint multidimensional model. Further, to assess the parameter recovery of the proposed model, a brief simulation study was conducted. For simplicity, only the MRM-MLRTM was assessed, because other two joint hierarchical models can be seen as its special cases.
Four dimensions were measured by 30 items, which means there are four latent abilities and four corresponding latent speeds. Q matrix was presented in Figure 3. For item parameters, d and ξ parameters generated from a bivariate normal distribution with mean vector (0, 4) and covariance matrix of [1, -0.2; -0.2, 0.25], in such setting, ρdξ = -0.4. ω parameters were set to 2, which were also similar to the estimates in the real data analysis. 1,000 respondents were simulated. (Θ, Τ) generated from an eightfold-variate normal distribution with mean vector (0, 0) and covariance matrix of

Analysis
The MRM-MLRTM was fitted to each of the 30 replications. In each replication, the numbers of chains, burn-in iterations and post-burn-in iterations were the same as those set in the empirical study. It appeared that convergence was well achieved. To evaluate parameter recovery, the bias and the root mean square error (RMSE) were computed as: , where  and υ are the estimated and true values of model parameters, respectively; R is the total number of replications. The correlation between estimated and true values (Cor) was also computed. Table 6 presents the bias, RMSE, and the Cor of item and person parameters. Across the board, all model parameters are well recovered. For item parameters, the recovery of time-intensity is the best, then is time-kurtosis, and the worst is item intercept. For person parameters, the recovery of latent speeds is better than that of latent abilities. The recovery of variance and covariance was of more interest in this study; estimated bias and RMSE are given in Table 7 for Σitem and Tables 8 and 9 for Σperson, respectively. In general, all of them are well recovered. The recovery of covariances is better than that of variances, and the recovery of time-related parameters (e.g., item intensity, covariance of item difficulty and time-intensity, latent speeds, and covariance of latent ability and latent speed) is better than that of time-unrelated parameters (e.g., item intercept and latent abilities). Overall, model parameters of the MRM-MLRTM can be recovered very well via the proposed full Bayesian MCMC estimation algorithm.

Conclusion and discussion
To capture the multidimensionality of latent speed, this study proposed a multidimensional log-normal RT model and a multidimensional hierarchical modeling framework. The PISA 2012 computer-based mathematics data were analyzed to illustrate the implications and applications of the proposed models. The results indicating that it is appropriate to consider the multidimensionality of latent speed and the multidimensionality of latent ability, simultaneously, in multidimensional tests when RTs were collected. A brief simulation study was used as well to further evaluate model parameter recovery. The results indicated that model parameters could be well recovered using the Bayesian MCMC approach.
The work presented in this article is only a first attempt to deal with the multidimensionality of latent speed. Despite promising results, further exploration is encouraged. First, the proposed RT model is a multidimensional extension of the classical log-normal RT model (van der Linden, 2006), multidimensional extensions of other possible RT models (Fox & Marianti, 2016;Wang, et al., 2013;Wang & Xu, 2015) could be explored and compared in the future. Second, in the proposed MLRTM, the latent speeds are assumed to be compensatory. As non-compensatory multidimensional models for response accuracy become popular in recent decades (DeMars, 2016;Embretson, 1984;2015;Templin & Henson, 2006;Wang & Nydick, 2015;Jiao, Lissitz, & Zhan, 2017), it is important to develop corresponding non-compensatory MRT models in the future. Third, in the proposed multidimensional hierarchical modeling approach, a multivariate normal distribution was used to describe the relationships among multidimensional latent speed and multidimensional latent ability. So, the number of total dimensions is twice as many as the number of dimensions that are measured by the test. For example, in above application example, there were eight dimensions, which may pose a challenge on parameter estimation. If the multidimensional latent ability and the multidimensional latent speed can each have a second-order (or bi-factor) structure, not only the parameter estimation challenge can be largely reduced but also the structures of latent ability and latent speed can be posited and tested. Fourth, applications of the proposed RT model, like detecting aberrant responses (e.g., rapid-guessing and cheating) in multidimensional tests, need further investigation. Fifth, analyzing students' growth is an important topic in educational and psychological research (e.g., von Davier, Xu, & Carstensen, 2011;Zhan, Jiao, & Liao, 2017).
How to employ the proposed multidimensional hierarchical modeling approach into longitudinal studies is also an interesting topic (e.g., Wang, Zhang, Douglas, & Culpepper, 2018).