Structural Equation Modeling Techniques for Estimating Score Dependability within Generalizability Theory-Based Univariate, Multivariate, and Bifactor Designs

1. Introduction

Although originally developed during the 1960s [1,2,3], generalizability theory (GT) continues to be used across numerous disciplines, in large part, because it can be applied to both objectively and subjectively scored measures, quantify effects of multiple sources of measurement error, and produce a wide variety of coefficients to assess the accuracy of observed scores for both norm- and criterion-referencing purposes. Introductions to doing traditional analysis of variance (ANOVA)-based GT analyses can be found in full-length books devoted exclusively to the topic [4,5,6,7,8,9]; chapters within measurement textbooks [10,11,12], research handbooks [13,14,15,16], encyclopedias [17,18,19,20,21,22,23], and edited volumes [25,26,27,28,29,30,31,32,33,34,35,36,37,38,39]; as well as articles or tutorials within professional journals [1,2,3,40,41,42,43,44,45,46,47,48,49,50,51,52,53,54,55,56,57,58,59,60,61,62,63,64]. Examples of content areas to which such analyses have been conducted over the last five years alone include medical education and training [65,66,67,68,69,70,71], radiology [72], rehabilitation [73], nursing, [74] pharmacology [75], K-12 writing skills [76,77,78,79], second language education [80,81,82], speech and hearing research [83,84], thinking skills and creativity [85], psychiatry and psychology [86,87,88,89,90], sports [91], among others.

Univariate ANOVA-based GT analyses can be run using variance component programs within comprehensive statistical packages such as SPSS, SAS, R, STATA, MATLAB, and Minitab (see, e.g., [92]) or standalone programs devoted exclusively to those purposes, which include the GENOVA suite (GENOVA [93], urGENOVA [94], and mGENOVA [95]), G_string IV (manual available at https://www.papaworx.com/HOP/Manual.pdf), EduG [9], and the gtheory package in R [62,96,97]). Most recently, applications of GT have been enhanced by conducting such analyses within structural equation model (SEM) frameworks (see, e.g., [59,60,62,64,98,99,100,101,102,103,104,105,106,107,108,109,110,111]). GT-based SEMs can be analyzed using numerous readily accessible statistical packages and provide effective methods for incorporating univariate, multivariate, and bifactor model designs, deriving confidence intervals for key parameters, adjusting for scale coarseness effects common when using binary or ordinal data, assessing scale viability, and handling missing data.

The most common indices of score accuracy reported in GT analyses include generalizability coefficients for norm-referencing purposes (e.g., rank ordering scores) and global and cut-score-specific dependability coefficients for criterion-referencing purposes (e.g., using absolute levels of scores for decision making). Initial applications of SEMs were limited to derivation of generalizability (G) coefficients within univariate designs [112,113] but later modified to allow for computation of dependability (D) coefficients [[4,98,103,114,115] and analysis of multivariate [64,108,109,110,111] and bifactor designs [64,101,102,103,104,108,109,111]. Whenever composite scores are reported in addition to subscale scores in practice, multivariate and bifactor GT designs produce more appropriate indices of score accuracy than would a univariate design for the composite because they take subscale score representation and interrelationships into account [2,4,8,64,102,111].

In a recent study by Lee and Vispoel [107]), an indicator mean-based procedure within SEMs was introduced to derive absolute error indices needed to estimate D coefficients for univariate designs that produced comparable or better results than previous methods for deriving such indices. Our primary goals in the study reported here are to replicate previous uses of SEMs within GT univariate designs using a new instrument and to expand upon such research by applying the indicator-mean method to derivation of absolute error indices for multivariate and bifactor designs based on original observed scores and ones corrected for scale coarseness effects. In doing so, the methods discussed here serve as guides to conducting complete GT analyses for common multivariate and bifactor designs with univariate analyses for individual subscales subsumed within those frameworks. We further provide computer code in R to enable readers to apply all illustrated techniques to their own data.

2. Background

2.1. GT designs

We will apply SEM procedures for analyzing persons × items (pi) single-facet and persons × items × occasions (pio) two-facet GT designs based on responses to a self-report measure (i.e., the shortened form of the Music Self-Perception Inventory (MUSPI-S) [116,117,118,119]). Each facet of measurement (items and occasions here) represents a domain to which results are generalized. Partitioning of observed score variance at the individual score level for these designs is described in Equations (1) and (2).

$p\times i\left(pi\right)\mathrm{design}:{\sigma }_{Y_{pi}}^2={\sigma }_p^2+{\sigma }_{pi,e}^2+{\sigma }_i^2.\left(1\right)$$p\times i\times o\left(pio\right)\mathrm{design}:{\mathsf{\sigma }}_{Y_{pio}}^2={\sigma }_p^2+{\sigma }_{pi}^2+{\sigma }_{po}^2+{\sigma }_{pio,e}^2+{\sigma }_i^2+{\sigma }_o^2+{\sigma }_{io}^2.\left(2\right)$ In the pi design, variance across all observed scores is partitioned into three components that represent persons $\left({\sigma }_p^2\right)$, items $\left({\sigma }_i^2\right)$, and the interaction between persons and items $\left({\sigma }_{pi,e}^2\right)$. In the pio design, observed score variance is partitioned into seven components that represent persons $\left({\sigma }_p^2\right)$, items $\left({\sigma }_i^2\right)$, occasions $\left({\sigma }_o^2\right)$and all possible interactions among persons, items, and occasions $\left({\sigma }_{pi}^2,{\sigma }_{po}^2,{\sigma }_{pio,e}^2,{\sigma }_{io}^2\right)$. The variance component for persons $\left({\sigma }_p^2\right)$in both designs represents universe score variance that parallels true score variance in classical test theory and communality in factor analysis. Interactions between persons and the measurement facets items and/or occasions represent sources of relative measurement error. The pi design includes only a single source of relative measurement error $\left({\sigma }_{pi,e}^2\right)$in contrast to the pio design that includes multiple sources of such error $\left({\sigma }_{pi}^2,{\sigma }_{po}^2,{\sigma }_{pio,e}^2\right)$. The subscript “,e” within a variance component term indicates that it also includes any remaining relative residual error not captured by other terms in the given design. Variance components involving persons are used to derive G coefficients for norm-referencing purposes, whereas those not involving persons reflect differences in mean scores for facet conditions within the design relevant for criterion-referencing purposes in which absolute score values are used for decision making. These components are combined with those for relative error to reflect absolute or total overall error when deriving D coefficients.

As already noted, GT-based analyses produce three primary indices of score accuracy that are represented in Equations (3)-(5): generalizability (G or Eρ²), global dependability (D or φ), and cut-score-specific dependability. We provide more detailed formulas to estimate these coefficients and related variance components within tables presented in later sections.

G coefficient = ${\displaystyle \frac{\mathrm{U}\mathrm{n}\mathrm{i}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{s}\mathrm{e} \mathrm{s}\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{e} \mathrm{v}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}{\mathrm{U}\mathrm{n}\mathrm{i}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{s}\mathrm{e} \mathrm{s}\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{e} \mathrm{v}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}+\mathrm{R}\mathrm{e}\mathrm{l}\mathrm{a}\mathrm{t}\mathrm{i}\mathrm{v}\mathrm{e} \mathrm{e}\mathrm{r}\mathrm{r}\mathrm{o}\mathrm{r} \mathrm{v}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}}$ (3)

Global D coefficient = ${\displaystyle \frac{\mathrm{U}\mathrm{n}\mathrm{i}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{s}\mathrm{e} \mathrm{s}\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{e} \mathrm{v}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}{\mathrm{U}\mathrm{n}\mathrm{i}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{s}\mathrm{e} \mathrm{s}\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{e} \mathrm{v}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}+\mathrm{A}\mathrm{b}\mathrm{s}\mathrm{o}\mathrm{l}\mathrm{u}\mathrm{t}\mathrm{e} \mathrm{e}\mathrm{r}\mathrm{r}\mathrm{o}\mathrm{r} \mathrm{v}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}}$ (4)

Cut-score-specific D coefficient =${\displaystyle \frac{\mathrm{U}\mathrm{n}\mathrm{i}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{s}\mathrm{e} \mathrm{s}\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{e} \mathrm{v}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}+{({\mu }_Y-cutscore)}^2}{\mathrm{U}\mathrm{n}\mathrm{i}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{s}\mathrm{e} \mathrm{s}\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{e} \mathrm{v}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}+{({\mu }_Y-cutscore)}^2+\mathrm{A}\mathrm{b}\mathrm{s}\mathrm{o}\mathrm{l}\mathrm{u}\mathrm{t}\mathrm{e} \mathrm{e}\mathrm{r}\mathrm{r}\mathrm{o}\mathrm{r} \mathrm{v}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}}}$ (5)

G coefficients are similar and sometimes identical to conventional alpha, split-half, parallel form, and test-retest reliability estimates in that they reflect relative differences in scores across persons (see [59] for further details) but are interpreted in relation to all possible facet conditions (items, occasion, raters, etc.) within the targeted assessment domains of interest. Equation (3) for a G coefficient would represent universe score (or person) variance divided by the sum of person variance and all error variance components involving persons $\left(\mathrm{i}.\mathrm{e}.,{\sigma }_p^2 \mathrm{a}\mathrm{n}\mathrm{d} {\sigma }_{pi,e}^2 \mathrm{i}\mathrm{n} \mathrm{t}\mathrm{h}\mathrm{e} pi \mathrm{d}\mathrm{e}\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n} \mathrm{a}\mathrm{n}\mathrm{d} {\sigma }_p^2,{\sigma }_{pi}^2,{\sigma }_{po}^2, \mathrm{a}\mathrm{n}\mathrm{d} {\sigma }_{pio,e}^2 \mathrm{i}\mathrm{n} \mathrm{t}\mathrm{h}\mathrm{e} pio \mathrm{d}\mathrm{e}\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\right)$. Equation (4) for global dependability resembles Equation (3), but with all variance components for facets and their interactions $\left(\mathrm{i}.\mathrm{e}.,{\sigma }_i^2 \mathrm{i}\mathrm{n} \mathrm{t}\mathrm{h}\mathrm{e} pi \mathrm{d}\mathrm{e}\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n} \mathrm{a}\mathrm{n}\mathrm{d} {\sigma }_i^2,{\sigma }_o^2, \mathrm{a}\mathrm{n}\mathrm{d} {\sigma }_{io}^2 \mathrm{i}\mathrm{n} \mathrm{t}\mathrm{h}\mathrm{e} pio \mathrm{d}\mathrm{e}\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\right)$combined with relative error variance components to represent absolute error in the denominator of the equation. Consequently, global D coefficients can be no larger than G coefficients and equal to them only when all facet condition means are identical.

Global D coefficients broaden the conceptualization of measurement error to include mean differences in scores and thereby reflect the contribution of the assessment procedure to overall dependability when making criterion-referenced interpretations of scores [25,28,114,115]. Finally, cut-score-specific dependability represented in Equation (5) parallels Equation (4), but with the squared difference between the grand score mean and cut score added to both the numerator and denominator of the equation. Accordingly, the value of this coefficient can change depending on the value of the cut score and represents dependability specific to that cut score. Conceptually, cut-score-specific D coefficients reflect the contribution of the assessment procedure to the decision made from the cut score over what would be expected by chance agreement [25,28,114,115]. These coefficients are especially useful for gauging accuracy in determining whether an individual’s standing truly falls above or below the targeted cut score. Like conventional reliability estimates, G and D coefficients can vary from 0 to 1 with higher values representing greater accuracy in scores for their intended purposes.

2.2. Representing univariance, multivariate, and bifactor GT designs within SEMs

In the data analyses presented here, we focus on univariate, multivariate, and bifactor pi designs based on a single occasion of administration in which items are the sole measurement facet of interest, and univariate, multivariate, and bifactor pio designs based on two occasions of administration with both items and occasions serving as measurement facets. Our illustrations represent data collected using the Instrument Playing, Reading Music, Listening, and Composing subscales from the shortened form of the Music Self-Perception Inventory (MUSPI-S; [116,117,118,119]). Each subscale consists of four items with item scores summed across all subscales to create a composite score that is used here to represent overall perceptions of music proficiency within the multivariate and bifactor designs.

2.3. Univariate GT designs

GT persons × items (pi) single-facet univariate designs. In Figure 1, we depict SEM diagrams for univariate pi and pio designs. The top diagram in Figure 1 represents a pi SEM for the Instrument Playing subscale. The SEM has a single factor for the construct of interest that is linked to all items measuring that construct. Item loadings are set equal to one, and all item uniquenesses are set equal. In total, two parameters are estimated that represent the variance component for person or universe scores (${\sigma }_p^2$) and the variance component for relative measurement error across items (${\sigma }_{pi,e}^2$).

When using the indicator-mean method [107], the remaining variance component for items (${\sigma }_i^2$) is estimated using intercepts for each item that are equivalent to their corresponding means. Specifically, the squared differences between each item mean and grand mean across items are summed and divided by the number of items minus one as shown in Table 1. Once the three variance components of interest (${\sigma }_p^2,{\sigma }_{pi,e}^2, \mathrm{a}\mathrm{n}\mathrm{d} {\sigma }_i^2$) are estimated. They can be placed in the general equations shown in Table 2 to estimate G, global D, and cut-score-specific D coefficients for the pi design.

GT persons × items × occasions (pio) two-facet univariate designs. The bottom diagram in Figure 1 represents a SEM for the MUSPI-S's Instrument Playing subscale for a pio design with four items administered on two occasions. When measuring psychological traits, the pio design is generally preferred over the pi design because it allows for separation of three key sources of measurement error and reduces confounding of trait and measurement error variance typically found in the pi design. Within the pio SEM shown in Figure 1, the factor for the construct of interest is linked to each item on each occasion. Separate factors are included for each item across all occasions and for each occasion across all items. All factor loadings are set equal to one, with item variances set equal, occasion variances set equal, and uniquenesses set equal. Collectively, these constraints result in estimation of four parameters to represent the variance of person or universe scores for the construct of interest (${\sigma }_p^2$) and three sources of relative measurement error (${\sigma }_{pi}^2,{\sigma }_{po}^2,\mathrm{a}\mathrm{n}\mathrm{d}{\sigma }_{pio,e}^2$) common when using objectively scored measures such as Likert-style questionnaires or multiple-choice tests in which all scorers would obtain the same results.

Within the psychological research literature, ${\sigma }_{pi}^2,{\sigma }_{po}^2,\mathrm{a}\mathrm{n}\mathrm{d}{\sigma }_{pio,e}^2$are often respectively referred as specific-factor error (or method effects), transient error (or state effects), and random-response error (or within-occasion “noise”; [120,121,122], also see [123]). Specific-factor error represents enduring person-specific effects on scores that are unrelated to the construct of interest such as understandings of words within items and response options. Transient error represents consistent effects on scores within a given occasion that do not generalize across occasions such as a respondent’s disposition, mindset, and physiological condition at that time as well as his or her reactions to environmental and administration conditions. Random-response error refers to fleeting motion-to-moment effects within an occasion such as momentary lapses of attention, distractions, and other effects that follow no systematic pattern. Within pi designs, universe score and transient error are confounded within the person variance component (${\sigma }_p^2$), as are specific-factor and random-response error within the relative measurement error component ${(\sigma }_{pi,e}^2).$An important advantage of the pio design is that these sources of variance can be separated to provide more appropriate estimates of score accuracy and overall measurement error. Variance component formulas for each source of measurement error are provided in Table 1.

When using the indicator-mean method to derive the remaining variance components (${\sigma }_i^2,{\sigma }_o^2,\mathrm{a}\mathrm{n}\mathrm{d}{\sigma }_{io}^2$) within the pio design, intercepts will represent means for all combinations of items and occasions. Once estimated, these means can be integrated into the formulas shown in Table 1 to obtain the remaining variance components and insert them into the formulas shown in Table 2 to estimate related indices of score dependability and proportions of measurement error (see [107] for further details).

2.3. Multivariate GT designs

GT persons × items (pi) single-facet multivariate design. Multivariate GT designs best represent indices of score accuracy when both subscale and composite scores are reported in practice. Such designs also can produce correlation coefficients corrected for all sources of measurement error included within a design to provide further insights into subscale score dimensionality, overlap, interrelationships, and validity. Embedded within the overall multivariate design are the same univariate analyses already described for each individual subscale. Variance components for the composite score, in contrast, entail formulas based on the variance components for each subscale, the covariances between each pair of subscale scores, and eventual weighting of each subscale when forming the composite (see Table 3 and [8,64,109]).

Figure 2 includes SEM diagrams for multivariate pi and pio designs that can be used to derive variance and covariance components for computing G coefficients for subscale and composite scores. The diagrams respectively represent the 4-item subscales from the MUSPI-S (Instrument Playing, Reading Music, Listening Skill, and Composing Ability) mentioned earlier administered on one or two occasions. Scores for each individual subscale are modeled and constrained in the same way as in a univariate analysis but allowed to covary/correlate with each other. Within the pi design, eight variance components (${\sigma }_p^2\mathrm{a}\mathrm{n}\mathrm{d}{\sigma }_{pi,e}^2$ for each subscale) and six covariance components (one for each possible pair of subscale scores) are estimated. The variance component for items (${\sigma }_i^2$) within each subscale is computed in the same way as described for the univariate design. ${\sigma }_p^2,{\sigma }_{pi,e}^2,\mathrm{a}\mathrm{n}\mathrm{d}{\sigma }_i^2\mathrm{i}\mathrm{n}\mathrm{d}\mathrm{i}\mathrm{c}\mathrm{e}\mathrm{s}\mathrm{f}$or the composite score can be estimated using the formulas shown in Table 3 and corresponding generalizability and dependability coefficients using the formulas shown in Table 2.

GT persons × items × occasions (pio) two-facet multivariate design. The pio univariate design for each embedded subscale within the multivariate design has additional factors for each item across occasions and for each occasion across items. Scores for each pair of subscales again are allowed to covary/correlate but to the same degree across occasions. Transient error (${\sigma }_{po}^2,$) indices also are allowed to covary/correlate when all measures are administered together within a common occasion. In all, sixteen variance components ${(\sigma }_p^2,{\sigma }_{pi}^2,{\sigma }_{po}^2,\mathrm{a}\mathrm{n}\mathrm{d}{\sigma }_{pio,e}^2\mathrm{f}\mathrm{o}\mathrm{r}\mathrm{e}\mathrm{a}\mathrm{c}\mathrm{h}\mathrm{s}\mathrm{u}\mathrm{b}\mathrm{s}\mathrm{c}\mathrm{a}\mathrm{l}\mathrm{e}),$ six covariance components for person subscale scores (one for each possible pair of subscale scores), and six covariances for transient errors (${\sigma }_{{po}_{S1}}$) are estimated. The ${\sigma }_i^2,{\sigma }_o^2,\mathrm{a}\mathrm{n}\mathrm{d}{\sigma }_{io}^2$ variance components for each subscale can be derived in the same ways described in the univariate designs, and those for the composite score using the formulas provided in Table 3. Formulas for deriving relevant generalizability and dependability coefficients for both subscale and composites scores are given in Table 2.

Correcting correlation coefficients for measurement error. In addition to providing appropriate indices of accuracy for composite and subscale scores, multivariate designs can yield correlations between all pairs of subscale scores corrected for the sources of measurement error estimated within the design. Corrected correlations can be conceptualized in relation to the formula first proposed by Spearman ([124], also see [125]), shown in Equation (6), in which the correlation coefficient between observed scores for the pair of measures of interest is divided by the square root of the product of their corresponding reliability coefficients to estimate the correlation between true scores for the targeted measures that is free of measurement error. In applications of GT, G coefficients would be substituted for conventional reliability coefficients, and universe scores for true scores.

${\widehat{\rho }}_{T_x{T}_y}=\hspace{0.33em}{\displaystyle \frac{r_{xy}}{\sqrt{r_{xx\text{'}}*r_{yy\text{'}}}}},\left(6\right)$

$\mathrm{where}{\widehat{\rho }}_{T_x{T}_y}=\mathrm{estimated}\mathrm{correlation}\mathrm{between}\mathrm{true}\mathrm{scores}\mathrm{for}\mathrm{measures}X\mathrm{and}Y,$

$r_{xy}=\mathrm{observedcorrelation}\mathrm{coefficient}\mathrm{between}\mathrm{measures}\mathrm{X}\mathrm{and}\mathrm{Y},$

$r_{xx\text{'}}=\mathrm{reliability}\mathrm{coefficient}\mathrm{for}\mathrm{measure}\mathrm{X},$

$\mathrm{and}{r}_{yy\text{'}}=\mathrm{reliability}\mathrm{coefficient}\mathrm{for}\mathrm{measure}\mathrm{Y}.$

2.4. Bifactor GT designs

Bifactor and multivariate GT designs both can be used to simultaneously partition score variance at subscale and composite levels and distinguish multiple sources of measurement error within pio designs. However, bifactor designs further allow for partitioning of universe score variance into general and group factor effects to produce indices reflecting just general factor effects, just group factor effects, or both effects combined. General factor effects reflect common explained variance shared across all indicators, whereas group factor effects reflect unique explained variance, unrelated to general factor variance, that is shared by all indicators representing a given subscale.

Bifactor models produce four key coefficients in addition ones already discussed: omega total composite, omega total subscale, omega hierarchical composite, and omega hierarchical subscale [11,126–130]. Omega total coefficients for composite and subscale scores represent proportions of variance accounted for by both general and group factor effects. They parallel overall G coefficients for the pi and pio univariate and multivariate GT designs, except that universe score variance represents the sum of general and group factor variances (see Table 2). Omega hierarchical composite score coefficients represent the proportion of variance accounted for by the general factor alone, whereas omega hierarchical subscale coefficients represent the proportion of variance accounted for by the group factor alone. We provide formulas for estimating variance components in Table 4 that can be inserted into formulas shown in Table 2 to derive G, global D, cut-score-specific D, and omega coefficients for pi and pio bifactor designs.

GT persons × items (pi) single-facet bifactor design. The pi bifactor SEM representing MUSPI-S scores is shown in the top diagram within Figure 3. The general factor is linked to all items, with independent group factors linked only to items included within each subscale. To allow for differential general factor effects across subscales, model identification constraints differ somewhat from those in the previous designs. Specifically, variances for the general factor and loadings for the group factors are set equal to one, and general factor loadings, group factor variances, and uniquenesses are estimated but set equal within but not across subscales. In all, twelve parameters are estimated (λ, ${\sigma }_{grp}^2$& ${\sigma }_{pi,e}^2$ for each subscale).

GT persons × items × occasions (pio) two-facet bifactor design. The bottom diagram in Figure 3 represents the pio bifactor design. It has the same constraints as the pi design, but with additional factors included for occasions and items. Item and occasion factor loadings are set equal to one, with item factor variances, occasion factor variances, and uniquenesses estimated. As in the multivariate pio design, separate occasion factors are included for each subscale that are allowed to covary/correlate with each other but to the same degree across occasions. Five parameters are estimated for each subscale (λ, ${\sigma }_{grp}^2,{\sigma }_{pi}^2,{\sigma }_{po}^2$, ${\&\sigma }_{pio,e}^2$) as well as six additional covariances to model possible correlated within-occasion transient error effects.

2.5. Using GT multivariate and bifactor designs to evaluate subscale viability.

An important question to consider whenever using measures that produce both composite and subscale scores in practice is the extent to which subscale scores yield useful information beyond the composite score. To address this question, Haberman ([131], also see [132,133,134,135]) devised a classical test theory-based procedure to determine whether a subscale’s true scores are better estimated using subscale or composite observed scores. Vispoel and colleagues [64,103,104,110,111] later adapted this procedure to single- and multi-facet GT multivariate and bifactor designs by replacing true score with universe score estimation.

Haberman’s method is based on comparison of indices for subscale and composite scores reflecting proportional reduction in mean-squared error (PRMSE). The PRMSE for the subscale is equivalent to its conventional reliability or GT-based generalizability coefficient, whereas the PRMSE(C) for the composite can be derived using Equation 7.

$\mathrm{PRMSE}\left(C\right)=r_{T_{S_j},T_C}^2\mathrm{*}{r}_{X_C,X_C\mathrm{\text{'}}}={\displaystyle \frac{{\widehat{\sigma }}_{T_{S_j},T_C}^2}{{\widehat{\sigma }}_{T_{S_j}}^2\mathrm{*}{\widehat{\sigma }}_{T_C}^2}}\mathrm{*}{\displaystyle \frac{{\widehat{\sigma }}_{T_C}^2}{{\widehat{\sigma }}_{X_C}^2}}={\displaystyle \frac{{\left({\widehat{\sigma }}_{T_{S_j}}^2+{\sum }_{j\ne k}{\widehat{\sigma }}_{T_{S_j},T_{S_k}}\right)}^2}{{\widehat{\sigma }}_{T_{S_j}}^2\mathrm{*}{\widehat{\sigma }}_{X_C}^2}},\left(7\right)$ where $T=\mathrm{true}\mathrm{score},X=\mathrm{observed}\mathrm{score},S=\mathrm{subscale},C=\mathrm{composite}\mathrm{score},\mathrm{and}$$r_{X_C,X_C\mathrm{\text{'}}}=\mathrm{composite}\mathrm{reliability}.$

Conceptually, a PRMSE index represents an estimate of the proportion of true or universe score variance in the present context accounted for by the targeted observed scores (subscale or composite). Once PRMSEs are obtained for a subscale and its associated composite, they can be inserted into Equation 8 to form a value-added ratio (VAR; see [136]). Subscale viability is increasingly supported as VARs deviate upwardly from 1.00.

$\mathrm{Value}-\mathrm{Added}\mathrm{Ratio}\left(\mathrm{VAR}\right)=\hspace{0.33em}{\displaystyle \frac{PRMSE\left(Subscale\right)}{PRMSE\left(Composite\right)}}\hspace{0.33em}\left(8\right)$

2.6. Comparing GT univariate, multivariate, and bifactor designs.

As previously noted, univariate and multivariate GT designs will produce the same results for individual subscales because univariate designs for subscales are embedded within the overall multivariate design. However, two additional benefits of multivariate over univariate designs described earlier are that they produce correlation coefficients between all pairs of subscale scores corrected for the sources of measurement error estimated within the design and yield more appropriate indices of generalizability and dependability for composite scores by taking subscale representation and interrelationships into account. Common findings across recent studies include stronger relationships between subscale scores when corrected for measurement error, and G and D coefficients for composite scores within multivariate designs that generally exceed those derived strictly from univariate designs in which subscale representation and interrelationships are ignored [64,110,111].

Either GT multivariate or bifactor designs can produce appropriate G and D coefficients at both composite and subscale levels as well as VARs for all subscale scores. In recent studies of personality constructs, GT multivariate and bifactor designs have produced highly comparable G coefficients, D coefficients, and subscale VARs [64,102,111] but with subscale and composite scores partitioned into general and group factor effects within the bifactor designs to provide additional insights into score dimensionality and overlap among constructs. In the vast majority of bifactor model studies, proportions of general factor exceed proportions of group factor variance at both composite and subscale levels (see, e.g., [101,102,111,128]).

2.7. Further advantages of using SEMs to perform GT analyses

Two additional benefits of conducting GT analyses using SEMs to be demonstrated here are to derive Monte Carlo-based confidence intervals for G, D, and omega coefficients and to use estimation procedures that correct for scale coarseness effects commonly encountered when analyzing dichotomous or ordinal-level data. When doing SEM analyses using the lavaan package in R [137,138], Monte Carlo-based confidence intervals [139] can be derived for nearly any parameter of interest through linkages with the semTools package [140], and dichotomous and ordinal data can be transformed to continuous latent variable metrics using diagonally weighted least squares (WLSMV in R) or other relevant estimation procedures (see e.g., [62,107]). In general, differences in G and D coefficients between observed score and continuous latent variable metrics diminish as numbers of scale points increase with the largest differences observed when items have only two response options [62,141]. Although the data analyzed here had no missing values, the anxiliary, BootMiss-class, and bsMissBoot routines within the semTools package can be linked to lavaan to handle missing data using auxiliary information, multiple imputation, bootstrapping, and related procedures (see [140] for further details).

3. This Investigation

Our main purpose within the present study is to demonstrate how indicator-mean and related procedures can be integrated into SEMs to allow for complete analyses of GT-based univariate, multivariate, and bifactor designs on both observed score and continuous latent response variable metrics. To evaluate the congruence of observed score results between SEM and ANOVA-based procedures, we compare G coefficients, D coefficients, and variance components obtained from the univariate and multivariate GT SEM designs to those obtained from the conventional package mGENOVA [95], which is often considered the gold standard when analyzing multivariate designs. Further comparisons of results involving the GT SEMs are made for composite and subscale scores across the multivariate and bifactor designs, across observed and continuous latent response variable score metrics, and across numbers of item scale points (2, 4, & 8).

In relation to the previous research studies cited here, we anticipate the following results:

• G coefficients, D coefficients, and variance components obtained from the GT-based univariate and multivariate SEMs will be highly congruent with those obtained from mGENOVA.
• Multivariate and bifactor GT SEMs will yield comparable G and D coefficients for subscale and composite scores.
• G and D coefficients for the pi designs will exceed those for the pio designs due to control of fewer sources of measurement error.
• Across all multivariate designs, correlation coefficients between scale scores will be higher after correcting for measurement error, but the difference between corrected and uncorrected coefficients will be greater in pio than in pi designs.
• General factor effects will exceed group factor effects at both subscale and composite levels within the bifactor designs.
• Similar patterns of VARs for subscales will be found across multivariate and bifactor designs.
• Composite and subscale scores will be affected by specific-factor (method), transient (state), and response-random (within-occasion noise) measurement error within the pio designs, but those effects be greater overall at the subscale than composite level due to inclusion of fewer item scores.
• Differences in G and D coefficients for two, four, and eight item scale points will be greater on observed score than on continuous latent response variable metrics.
• G and D coefficients will be greater on continuous latent response variable than on observed score metrics but to diminishing degrees with increases in numbers of item scale points.

4. Methods

4.1. Participants, Measures, and Procedure

We used the same dataset from Lee and Vispoel [107] in which 511 college students from educational psychology and statistics courses within a large Midwestern university (77.50% female, 82.00 Caucasian, mean age = 21.16) completed the full form of the adult level of the Music Self-Perception Inventory (MUSPI) [116,142,143,144] on two occasions, a week apart. However, for sake of efficiency, variety, and comparison, we analyzed responses to the same subscales (Instrument Playing, Reading Music, Listening, and Composing) using items from the shortened form of the MUSPI (MUSPI-S), all of which are included in the full form.

Each subscale within the MUSPI-S includes four positively phrased items answered along an 8-point item response metric with the following options: (1) Definitely False, (2) Mostly False, (3) Moderately False, (4) More False Than True, (5) More True Than False, (6), Moderately True, (7) Mostly True, (8) Definitely True. We computed subscale scores by adding responses to all items within the given subscale, and composite scores by summing all subscale scores to represent Overall Music Proficiency. Psychometric evidence supporting the use of the MUSPI-S scores includes alpha reliability coefficients for subscale scores no lower than 0.91, confirmatory factor analyses of responses yielding excellent fits to the data, and verification of expected relationships of MUSPI-S subscale scores with each other and with a wide variety of external criterion variables (see, e g., [116,117,118,119]). To evaluate effects of number of item scale points across analyses, we recoded original scores of 1-2, 3-4, 5-6, and 7-8, respectively, to 1, 2, 3, 4 to reduce responses to four scale points, and recoded original scores of 1-4 and 5-8, respectively, to 1 and 2 to reduce responses to two scale points.

4.2. Analyses

Initial analyses included estimation of means, standard deviations, alpha reliability coefficients, and test-retest reliability coefficients for MUSPI-S subscale and composite scores. Subsequent analyses included derivation of variance components, G coefficients, D coefficients, correlation coefficients (corrected and not corrected for measurement error), and VARs for relevant MUSPI-S scales across the pi and pio univariate, multivariate, and bifactor designs. The pi designs include data collected on the first measurement occasion only, and the pio designs include data collected on both occasions. Within the multivariate and bifactor pi and pio designs, items are nested within subscales, and occasions are crossed with subscales in the pio designs.

All SEM-based indices were estimated using procedures within the computer package R. For sake of comparison, variance components, G coefficients, and D coefficients for observed scores also were derived for univariate and multivariate designs using mGENOVA [95]. To parallel conventional ANOVA-based procedures for observed scores, SEM-based analyses were based on unweighted least squares (ULS) estimation. To convert observed score results to continuous latent response variable metrics within the SEM-based analyses, we used WLSMV estimation within the lavaan package [137,138], which is described by its authors as a diagonally weighted least squares procedure with robust standard errors and a mean and variance adjusted test statistic. We also derived 95% Monte Carlo based confidence intervals [139] using the semTools package [140] to gauge precision in estimating G, D, and omega coefficients. More detail and computer code for deriving all key indices are provided in our Supplementary Materials.

5. Results

5.1. Means, Standard Deviations, and Conventional Reliability Estimates for MUSPI-S Scores

Table 5 includes means, standard deviations, and alpha reliability estimates for all MUSPI-S scales and response metrics within each occasion, as well as test-retest reliability estimates across occasions. In relation to individual item scale metrics, means for all scales fall near or below their respective scale midpoint values of 1.5, 2.5, and 4.5 for the two-, four-, and eight-point metrics, with the Composing subscale always having the lowest mean. Given the possible range of scores on each metric (1-2, 1-4, and 1-8), standard deviations for each scale reflect a high degree of variability, respectively, ranging from 0.38-0.46, 0.93-1.13, and 1.96-2.37. These results make sense given the likely heterogeneity of music-related skills within this college student sample. Across scales and response metrics, alpha coefficients within occasions are uniformly high, ranging from 0.91 to 0.98. Test-retest coefficients range from 0.80 to 0.93 and are lower than corresponding alpha coefficients on all instances, thereby reflecting lower occasion-to-occasion than item-to-item consistency for all MUSPI-S scales represented here.

5.2. Key indices for GT pi Designs

Univariate and multivariate analyses using SEMs and mGENOVA. Results for G coefficients, global D coefficients and variance components for observed composite scores within the SEM-ULS multivariate designs and their embedded univariate designs for subscales shown in Table 6 are highly consistent with those obtained from the mGENOVA package. Between the two approaches, G coefficients are identical, global D coefficients differ by no more than 0.001, and variance components differ by no more than 0.002.

G and global D coefficients within the SEMs. Across scales and numbers of item scale points, G and global D coefficients for the SEMs are uniformly higher for WLSMV than for ULS estimation, but these differences are noticeably greater for subscale than for composite scores. This indicates that scale coarseness effects can be especially pronounced when using a small number of items. These differences are further exacerbated, but to a lesser degree here, when using more limited numbers of scale points.

In relation to widths of confidence intervals (i.e., differences between upper and lower limits) shown in Table 6, precision in estimating G and global D coefficients is weakest with two item scale points and improves with increases in numbers of item scale points for ULS but not always for WLSMV estimates. With eight-point items, widths of confidence intervals are uniformly narrower for ULS than for WLSMV estimation, thereby indicating greater relative increases in precision with increases in item scale points on the observed score metric than on the continuous latent response variable metric.

Cut-score-specific D coefficients. For purposes of comparison, cut-score-specific D coefficients for composite and subscale scores using both ULS and WLSMV estimation are depicted on Z-score metrics (M = 0, SD = 1) in Figure 4. Across scales and estimation procedures, dependability is lowest at the scale mean and progressively increases as scores deviate farther and farther away from the mean. Consistent with the results for G and global D coefficients, composite cut-score-specific D coefficients based on ULS estimates always exceed those for subscales at common standard deviation distances away from the scale mean. As expected, differences in G coefficients, global D coefficients (Table 6), and cut-score-specific D coefficients (Figure 4) for observed scores based on ULS estimation are typically greater between two and four scale points than between four and eight scale points. However, this pattern does not hold for the continuous latent variable response metric based on WLSMV estimation in which corrections for coarseness are greater for two-point than for four- or eight-point scales. However, as scores become increasingly extreme for all scales and estimation procedures, cut-score-specific coefficients across numbers of item scale points begin to coincide, thereby illustrating that numbers of scale point effects on score dependability are noticeably greater within the middle of score distributions than at the extremes.

Subscale intercorrelations. As noted earlier, an important advantage of GT multivariate analyses is to produce correlation coefficients between subscale scores corrected for measurement error. As expected, the corrected and uncorrected correlation coefficients between subscales in Table 7 reveal that the relationship between each pair for measured constructs is greater than would otherwise be inferred (ULS: uncorrected $\overline{r}$ = 0.664, mean corrected $\overline{r}$ = 0.699; WLSMV: uncorrected $\overline{r}$ = 0.748, mean corrected $\overline{r}$ = 0.764). The modest average differences between corrected and uncorrected coefficients observed here result from the generally high G coefficients for the subscales using either estimation method. Such differences are expected to increase when G coefficients are corrected for additional sources of measurement error (see Equation (6) and results of the pio designs presented in subsequent sections).

GT bifactor analyses for subscale and composite scores. GT bifactor analyses provide an additional perspective on results by generally yielding G and D coefficients comparable to those obtained from parallel multivariate analyses but further subdividing universe score variance into general and group factor effects. Comparing G and global D coefficients for the bifactor designs in Table 8 to those for the multivariate designs in Table 6 verifies that this congruence holds, with coefficients being identical with ULS estimation and varying by no more than 0.001 with WLSMV estimation. Accordingly, differences in results between ULS and WLSMV estimation and across numbers of scale points previously discussed for multivariate designs hold here, as would the differences between cut-score specific D coefficients shown in Figure 4.

Relative general and group factor effects on composite and subscale score variance can be examined using omega hierarchical composite (${\widehat{\omega }}_{H_{Composite}}$) and omega hierarchical subscale (${\widehat{\omega }}_{H_{Subscale}}$) coefficients that respectively represent proportions of explained general and group factor effects (see Table 8). Across numbers of scale points, these coefficients reveal that the strongest general and weakest group factor effects are found for the Overall Music Proficiency composite (${\widehat{\omega }}_{H_{Composite}}$: 0.869-0.938; ${\widehat{\omega }}_{H_{Subscale}}$: 0.058-0.108), Instrument Playing subscale (${\widehat{\omega }}_{H_{Composite}}$: 0.744-0.888; ${\widehat{\omega }}_{H_{Subscale}}$: 0.098-0.196), and Reading Music subscale (${\widehat{\omega }}_{H_{Composite}}$: 0.741-0.889; ${\widehat{\omega }}_{H_{Subscale}}$: 0.101-0.206), whereas the weakest general and strongest group effects are found for the Composing (${\widehat{\omega }}_{H_{Composite}}$: 0.468-0.651; ${\widehat{\omega }}_{H_{Subscale}}$: 0.327-0.443) and Listening (${\widehat{\omega }}_{H_{Composite}}$: 0.529-0.695; ${\widehat{\omega }}_{H_{Subscale}}$: 0.291-0.405) subscales.

Precision in estimating G, global D, and omega coefficients in relation to the widths of confidence intervals within the bifactor designs display patterns of G and global D coefficients parallel to those within the corresponding multivariate designs. Widths are generally narrower for composite than for subscale scores, they progressively narrow with increasing numbers of item scale points with ULS estimates, and they narrow when moving from two to either four or eight item scale points with WLSMV estimation, but not necessarily when moving between four to eight item scale points. This same relative pattern of width differences holds for omega hierarchical coefficients.

Subscale VARs within GT multivariate and bifactor designs. To use both subscale and composite scores in practice, subscale scores should provide unique information beyond that represented within the composite that subsumes the subscale scores. The VARs for MUSPI-S subscale scores in Table 9 provide a useful mechanism to verify such expectations within both GT multivariate and bifactor designs. Despite the relatively high correlations between many pairs of subscale scores shown in Table 7, the VARs for all MUSPI-S subscales examined here exceed 1.00, thereby supporting their added value beyond the composite. Consistent with results for omega hierarchical coefficients, VARs are higher for the Composing and Listening subscales than for the Instrument Playing and Reading Music subscales.

5.2. Key indices for GT pio Designs

Univariate and multivariate analyses using SEMs and mGENOVA. As previously noted, pio designs generally provide more appropriate indices of score accuracy for trait-based measures because they allow for estimation of multiple sources of measurement error (specific-factor, transient, and random-response) likely to affect scores. In Table 10, we provide G coefficients; global D coefficients; proportions of specific-factor, transient, and random-response measurement error; and variance components for observed composite scores within the ULS-SEM multivariate designs and for observed subscale scores within the embedded univariate designs. Results in the table for mGENOVA and the ULS-SEMs are essentially equivalent, with reported indices being identical in nearly all instances and differing by no than 0.002 in other instances.

G and global D coefficients within the SEMs. Consistent with results for the pi designs, G and global D coefficients within the pio designs shown in Table 10 are uniformly greater for WLSMV than for ULS estimation, again implying that score accuracy is reduced due to the dichotomous or ordinal nature of the item response metrics. G and global D coefficients based on ULS estimation are more affected by changes in numbers of item scale points than those based on WLSMV estimation, with dichotomous item scales always producing the lowest observed score coefficients. In most cases, corrections for scale coarseness using WLSMV estimates and increasing numbers of item scale points using ULS estimates reduces all sources of measurement error, and this is especially true with two and four item scale points. In line with the conventional alpha and test-retest coefficients presented earlier, transient error (i.e., occasion) effects are greater than specific-factor error (i.e., item) effects for all scales, and this holds true across both estimation procedures.

As was the case for the pi designs, confidence intervals for G and global D coefficients shown in Table 10 within the pio designs reveal different patterns in precision for ULS and WLSMV estimates as numbers of item scale points increase. Widths of the confidence intervals for these coefficients progressively narrow with increases in numbers of item scale points for ULS estimates, but this is not always the case for WLSMV estimates. With WLSMV estimates, confidence interval widths are generally wider with two item scale points and similar with four and eight item scale points. That is, precision in estimating G, global D, and omega coefficients consistently improves with increases in numbers of item scale points on the observed score metrics but not necessarily on continuous latent response variable metrics.

Cut-score-specific D coefficients. Cut-score-specific D coefficients for all scales and both estimation procedures within the pio designs are plotted in Figure 5. Overall, trends observed for the pio designs parallel those for the pi designs but with coefficients being lower on average because additional sources of measurement error are represented. As before, cut-score-specific D coefficients across scales steadily increase as scores move farther and farther away from the scale mean, are higher for composite than for subscale scores, are higher on continuous latent response variable than on observed score metrics, and vary less with changes in numbers of item scale points on continuous latent response variable than on observed score metrics.

Subscale intercorrelations. Table 11 includes correlation coefficients corrected and uncorrelated for the three sources of measurement error estimated within the pio multivariate designs. As would be expected, corrected correlation coefficients again always exceed corresponding uncorrected coefficients (ULS: uncorrected $\overline{r}$ = 0.638, corrected $\overline{r}$ = 0.721; WLSMV: uncorrected $\overline{r}$ = 0.723, corrected $\overline{r}$ = 0.776). However, due to inclusion of additional sources of measurement error, corrected coefficients on average exceed those from the pi designs (ULS: $\overline{r}$ = 0.721 versus 0.699; WLSMV: $\overline{r}$ = 0.776 versus 0.764) as do the average differences between corrected and uncorrected coefficients across those designs (ULS: mean pio difference = 0.083, mean pi difference = 0.035; WLSMV: mean pio difference = 0.053, mean pi difference = 0.016).

GT bifactor analyses for subscale and composite scores. As previously noted, GT pio bifactor designs allow for partitioning of measurement error variance into three sources (specific-factor, transient, and random-response) and universe score variance into two sources (general factor and group factor). This partitioning is reflected in G coefficients, global D coefficients, proportions of measurement error, and variance components for the bifactor pio designs that appear in Table 12. G coefficients, global D coefficients, and proportions of measurement error for bifactor and corresponding multivariate designs are identical for subscale scores and differ by no more than 0.001 for composite scores when using ULS estimation. Slightly greater differences occur between the bifactor and multivariate designs when using WLSMV estimation with the maximum difference between composite score G or global D coefficients equaling 0.008. Although not depicted here, cut-score-specific D coefficients are virtually identical to those for the corresponding multivariate design shown in Figure 5.

Due to inclusion of additional sources of measurement error, G, global D coefficients, and most omega coefficients are uniformly lower and overall measurement error uniformly higher in the pio bifactor designs than in the pi designs. Patterns of general and group factor effects for the pio designs mirror those in the pi designs, with general factor effects being stronger and group factor effects being weaker for the Overall Music Proficiency composite (${\widehat{\omega }}_{H_{Composite}}$: 0.813-0.909; ${\widehat{\omega }}_{H_{Subscale}}$: 0.052-0.091), Instrument Playing subscale (${\widehat{\omega }}_{H_{Composite}}$: 0.712-0.862; ${\widehat{\omega }}_{H_{Subscale}}$: 0.097-0.167), and Reading Music subscale (${\widehat{\omega }}_{H_{Composite}}$: 0.710-0.861; ${\widehat{\omega }}_{H_{Subscale}}$: 0.101-0.192) than for the Composing (${\widehat{\omega }}_{H_{Composite}}$: 0.419-0.600; ${\widehat{\omega }}_{H_{Subscale}}$: 0.290-0.399) and Listening (${\widehat{\omega }}_{H_{Composite}}$: 0.499-0.660; ${\widehat{\omega }}_{H_{Subscale}}$: 0.248-0.300) subscales.

Confidence interval widths for estimating G, global D, and omega coefficients in Table 12 display patterns of precision in line with those described for the pi multivariate and bifactor designs and pio multivariate designs. Widths of the intervals are narrower for composite than for subscale scores, progressively narrower as numbers of item scale points increase with ULS estimation, and wider for two item scale points but similar for four and eight points with WLSMV estimation.

Subscale VARs within GT multivariate and bifactor designs. Consistent with many of the indices already reported, VARs for the pio multivariate and bifactor designs shown in Table 13 are very similar to each other and to those reported for the pi designs in Table 9. All MUSPI-S subscale scores exceed the threshold of 1.00 to support subscale viability, with the Composing subscale always yielding higher VARs than the other subscales, and the Listening subscale yielding higher VARs than Instrument Playing and Reading Music subscales in most instances.

6. Discussion

6.1. Overview

Although introduced to the research community over 60 years ago [1], GT continues to be used widely in both research and practice to evaluate the psychometric properties for scores yielded by a broad range of assessment procedures. Recent advances in computer technology and structural equation modeling techniques have further expanded accessibility to programs for doing GT analyses and increased the scope of such analyses. In the study reported here, we sought to integrate, synthesize, and apply newly developed GT-based SEM techniques for conducting complete GT analyses of scores from univariate, multivariate, and bifactor designs with varying numbers of item scale points. Illustrations were focused on objectively scored measures with items and occasions serving as universes of generalization, but the same techniques can be readily applied to subjectively scored assessments by substituting raters for either items or occasions.

6.2. Effectiveness of the Indicator-Mean Method

Univariate and multivariate designs. A central part of the present analyses was to extend applications of the indicator-mean method for deriving absolute error variance components and related D coefficients to multivariate and bifactor designs and replicate analyses for univariate designs from Lee and Vispoel [107] using a reduced length form of the MUSPI (MUSPI-S). Across MUSPI-S scales and response metrics for observed scores, SEMs with ULS parameter estimates yielded composite and subscale score G coefficients, global D coefficients, and variance components essentially the same as those produced by the mGENOVA package, which still is considered the gold standard for performing GT multivariate analyses. However, mGENOVA is limited to analysis of one and two facet designs, and although not considered here, SEMs using the indicator-mean approach can be further extended to multivariate and bifactor designs with more than two measurement facets. When Lee and Vispoel [107] did so with three-facet univariate designs, they found that the indicator-mean approach yielded more accurate absolute error and D associated coefficients than did previous methods used in GT SEMs [98], and we would expect similar results to hold for multivariate and bifactor designs with three or more measurement facets.

Bifactor designs. Although Vispoel and colleagues [64,102–104,108,109,111] recently extended GT techniques to bifactor designs using SEMs, there currently are no other formal GT-based computer packages for analyzing such designs for purposes of comparison. Nevertheless, the results found here are congruent with those from previous studies in showing that multivariate and bifactor analyses produce essentially the same G and D coefficients within pi and pio designs [64,102,109,111]. We further demonstrated here that this congruence holds for multivariate and bifactor designs when integrating the indicator-mean method into appropriate GT-based SEMs. These overall results across studies are not surprising given that second-order hierarchical and bifactor models can be considered reparameterizations of each other (see, e.g., [126,145,146]).

6.3. Univariate GT Analyses for Shortened and Full-Length Forms of the MUSPI-S

The full length MUSPI [142,143,144] and its shortened form (MUSPI-S; [116,117,118,119]), respectively, include twelve and four items per subscale. Consequently, analyzing scores from the MUSPI and MUSPI-S using the same dataset from Lee and Vispoel [107] facilitates comparisons of effects for numbers of subscale items on conventional and GT-based indices of score accuracy. Lee and Vispoel provided alpha and test-retest reliability estimates for the same four subscales analyzed here (Instrument Playing, Reading Music, Listening, and Composing) using response metrics with two, four, and eight points as well as G and global D coefficients for pi and pio designs for the Composing subscale. In Table 14, we provide those indices for MUSPI-S observed subscale scores here and for MUSPI observed subscale scores from Lee and Vispoel [107].

Mean alpha and test-retest coefficients shown in Table 14, respectively, range from 0.933-0.970 (M =0.954) and 0.858-0.913 (M =0.890) for MUSPI-S subscales compared to 0.957-0.980 (M =0.970) and 0.912-0.936 (M =0.927) for MUSPI subscales. Across the two inventories, alpha coefficients always exceed test-retest coefficients, and the magnitude of both coefficients increases with increases in numbers of item scale points. The greatest difference between the inventories is for test-retest coefficients with two item scale points (MUSPI: 0.912 vs MUSPI-S: 0.858) and the smallest is for alpha coefficients on occasion two with eight item scale points (MUSPI: 0.980 vs MUSPI-S: 0.970). In keeping with the conventional reliability estimates, G and global D coefficients for the Composing subscale across the two inventories in both the pi and pio designs increase in magnitude with increases in numbers of item scale points. Differences in these indices range from 0.026 to 0.032 in the pi designs and from 0.030 to 0.066 in the pio designs, always favoring the full-length MUSPI. The largest differences in the pio designs are for G (0.066) and global D (0.065) coefficients with two item scale points and the smallest are for G (0.030) and global D (0.030) coefficients with eight item scale points. As a result, the best way to approximate the accuracy of original MUSPI scale scores when using its shortened form is to retain its 8-point response metric. When doing so, MUSPI-S scores come reasonably close to approximating the accuracy of MUSPI scores while using only one-third of its original items.

Another noteworthy finding within the present analyses was that G, global D, and cut-score specific D coefficients for the MUSPI-S Composing (e.g., making up your own music) and Listening (e.g., identifying characteristics of music by ear) subscales were lower than those for Instrument Playing and Reading Music subscales across designs and response metrics. This might have resulted from listening ability and composing skills being less concrete to conceptualize and/or less familiar to the respondents. Nevertheless, overall patterns of results were highly consistent across subscales for observed scores in terms of accuracy improving with increases in item scale points, but with the greater improvements occurring when moving from two to four points than when moving from four to eight points. The clear message convened by these results is to avoid using dichotomous scales when measuring constructs like those considered here. Patterns of relative effects of different sources of measurement error within the pio designs also were highly consistent across subscales. In nearly all instances, transient error was highest, followed respectively by random response and specific-factor error. These results highlight the importance of retesting when measuring the current constructs, the likely overestimation of score accuracy when relying exclusively on single-occasion data, and the value of estimating effects for multiple sources of measurement error.

6.4. Multivariate GT Analyses of MUSPI-S Scores

Composite score results. Important benefits of using multivariate GT designs include simultaneous univariate analyses for all embedded subscales as already discussed, derivation of more appropriate indices of score accuracy and measurement error for composite scores, and estimation of correlation coefficients between subscale scores corrected for all sources of measurement error estimated within a design. Over 70 years ago, Standard D 6.3 (If a test is divided into sets of items of different content, internal consistency should be determined by procedures designed for such tests.) from the Technical Recommendations for Psychological Tests and Diagnostic Techniques [147] underscored the importance of adjusting reliability estimates for composite scores to reflect subscale representation and interrelations. Accordingly, Cronbach et al. [148] subsequently developed alpha coefficients for composite scores stratified by content categories and later applied these same ideas to multi-facet multivariate GT designs to account for additional sources of measurement error [93,95].

To illustrate the effects of content stratification on score accuracy indices using the present data, we can compare G coefficients for composites from the pi multivariate designs (0.977, 0.985, 0.988; see Table 7), which are equivalent to stratified alpha coefficients, to the non-stratified alpha coefficients representing composite scores on the first measurement occasion reported in Table 5 (0.954, 0.967, 0.971). Despite non-stratified alpha coefficients and observed subscale score intercorrelations being relatively high (see Table 5 & Table 7), non-stratified alpha coefficients for composite scores are from 0.017 to 0.022 lower than corresponding stratified coefficients. As subscale score intercorrelations decrease, these differences would likely further increase [2,148]. An important advantage of multivariate GT designs is that all derived G and D coefficients for composite scores are appropriately adjusted for subscale representation and their interrelationships. Unfortunately, in research studies and practical settings, stratified alpha coefficients and score accuracy indices from GT multivariate designs for composites are still rarely reported.

Composite versus subscale scores. Due to inclusion of more items, G and D coefficients for MUSPI-S observed composite scores were uniformly higher and corresponding proportions of measurement error uniformly lower than for observed subscale scores. However, the patterns of relative effects for overall designs (pi vs pio) and numbers of item scale items on observed composite score indices mirrored those for subscale scores, with higher G and D coefficients in the pi than in the pio designs and increasingly higher values for these coefficients as numbers of item scale points increased. As with subscales, greater increases in score accuracy and reductions in measurement error occurred for observed composite scores when moving from two to four item scale points than when moving from four to eight scale points, again illustrating the disadvantages of using dichotomously scored items within these self-report measures.

Correlation coefficients. The final unique feature of multivariate designs illustrated here, and one not shared with corresponding bifactor designs, is to produce correlation coefficients between pairs of subscale scores corrected for all sources of measurement error estimated within the analyzed design. Due to measurement error being present in all designs, corrected correlation coefficients always exceeded uncorrected ones, thereby implying that the underlying constructs measured by each pair of subscale scores are more strongly related than the uncorrected coefficients would suggest. Because overall measurement error was always greater in the pio than in the pi designs, corrected correlation coefficients in the pio designs always exceeded those in the pi designs. These results illustrate the importance of including all relevant sources of measurement error affecting scores within a GT design, not only in assessing score accuracy, but also in gauging the concurrent validity of universe scores. Other consistent findings for corrected correlation coefficients across all designs were in showing that, among the constructs measured, self-perceptions of abilities to play a musical instrument and read music were most strongly related (pi design: corrected $\overline{r}$ = 0.888; pio design: corrected $\overline{r}$ = 0.889)) and self-perceptions of abilities to compose music and read music were most weakly related (pi design: corrected $\overline{r}$ = 0.646; pio design: corrected $\overline{r}$ = 0.665).

6.5. Bifactor GT Analyses of MUSPI-S Scores

Although first described by Holzinger and colleagues in the 1930s [149,150], only in recent years have uses of bifactor models truly proliferated (see, e.g., [101,126,127,128,129,130]). A bifactor model is suitable for measures that represent hierarchically structured constructs in which a broad general domain factor affects responses to all items, and additional independent group factors affect responses only to those items intended to measure narrower subdomain constructs. In the present context, the broad factor represented self-perceptions of overall music proficiency and group factors represented self-perceptions of skill in playing musical instruments, reading music, listening, and composing. In contrast to GT univariate and multivariate analyses, universe scores in GT bifactor designs represent the additive sum of independent general and group factor effects. However, unlike the typical conventional single-occasion bifactor analyses that currently dominate the research literature (see, e.g., [128]), GT-based bifactor designs can produce global and cut-score-specific D coefficients at both composite and subscale levels and distinguish multiple sources of measurement error when items are administered over two or more occasions.

G coefficients, D coefficients, and proportions of measurement error for the present GT pi and pio bifactor designs for observed scores mirrored results from their corresponding multivariate designs, with score accuracy improving more between two and four item scale points than between four and eight scale points, and transient error being the predominant source of measurement error within the pio designs. The most important additional finding from the bifactor analyses, and one consistent with most bifactor analyses reported in the research literature (see, e.g., [101,102,103,111,128]), was that general factor exceeded group factor effects at both composite and subscale levels but to a greater extent at the composite level. Among subscales, Instrument Playing and Reading Music were more affected by general factor and less affected by group factor effects than were Composing and Listening. This result along with the correlations among subscale scores discussed in the previous section further verify that perceptions of overall music proficiency were more related to perceptions of performing and reading music than to listening to or composing music.

6.6. Other Noteworthy Aspects of the GT SEM Designs

Assessing subscale added value. We chose to report value-added indices as the primary basis for evaluating subscale viability here, because they are applicable to both GT multivariate and bifactor designs. In general, viability of scores from a given subscale is undermined by its overlap with scores from other subscales that comprise the composite. Yet despite the relatively high correlations between observed subscale scores in many instances, VARs for all subscales and designs exceeded 1.00, thereby supporting subscale viability within all contexts considered here. This is likely due, in part, to the typically high G coefficients for subscale scores across designs. Based on these results, reporting of both MUSPI-S subscale and composite scores would be justified for individuals like those sampled here.

Evaluating effects of scale coarseness. Applications of common statistical procedures (ANOVA, multiple linear regression, correlational indices, etc.) typically are governed by assumptions that scores are measured on equal interval scales. However, this is unlikely to be strictly true when using Likert-style questionnaires. To evaluate effects of possible violations of this assumption more thoroughly, we analyzed MUSPI-S data using item scales with two, four, and eight points and compared the results to those obtained from SEMs using WLSMV parameter estimation in which observed item responses are transformed to continuous latent response variable metrics presumed to be equal interval in nature. Such transformations are based on estimation of item thresholds that typically alter distances between observed scale points to conform to those within continuous standard normal score distributions [62]. Results reported here revealed that G and D coefficients on observed score metrics deviated most from those on continuous latent response variable metrics when dichotomously scored items were used but increased in congruence as numbers of item scale points increased. These results again serve to discourage use of the dichotomous scales and illustrate a mechanism that can be used to evaluate effects of scale coarseness using any number of item scale points.

A finding observed here and elsewhere (see, e.g., [62,107]) was that G and D coefficients on WLSMV metrics were somewhat higher for two than for four or eight observed scale points. These differences can occur for a variety of reasons including a positive bias sometimes observed when estimating score accuracy using dichotomously scored items (see, e.g., [141]), differences in the characteristics of the observed score distributions, and the after-the-fact conversion of the original eight-point item scale metric to two and four points. However, given the added information provided by individual descriptors of responses on the eight-point scale, results obtained on that metric using WLSMV estimation might be expected to best correspond to those that would be obtained from a scale that is truly continuous in nature [151].

Confidence intervals. The semTools package [140] can be linked to the lavaan package in R to derive Monte Carlo-based confidence intervals [139] for nearly any desired parameter. For sake of illustration, we derived 95% confidence intervals for G, global D, and omega hierarchical coefficients. In most instances, widths of the intervals for these indices narrowed with increases in numbers of observed item scale points but were generally much wider for two-point than for four- and eight-point item metrics. As with many other indices already discussed, these results again emphasize drawbacks in precision when using dichotomously scored items.

7. Final Conclusions

In writing this article, we sought to provide readers with a guide to analyzing a wide variety of GT-based designs by taking advantage of procedures made possible when using SEMs to conduct such analyses. The supplementary materials associated with this article include code in R and further detail in how to analyze all illustrated designs. We hope that these resources prove useful in evaluating and understanding the quality and nature of scores used for either norm- or criterion-referencing purposes.