3.1 Identification issues

Constraints must be applied to achieve model identification, and this is where different parameterizations of the auxiliary model come into play. Various sets of conditions exist that are sufficient to ensure identification, and they vary based on assumptions made about the thresholds (whether they are time-varying, time-invariant, or fixed) and the means and variances of the underlying variables (whether they are fixed on all occasions, all free but one, or freely estimated on all time points). Different software may adopt distinct identification conditions, leading to variations in parameter estimates and model fit across programs. Table 2 details the main parameterizations in the literature that we discuss.

Table 2

Alternative sets of identification constraints for the auxiliary model in presence of ordinal data

• Standard parameterization. The first row of Table 2 refers to the set of constraints commonly adopted in single-population or cross-sectional applications. It consists of fixing the expected value and standard deviation of $Y_{it}^{*}$ equal to zero and one, respectively, on each occasion. When this standard parameterization is adopted, the thresholds are freely estimated and assumed to be time-varying. These $T(C-1)$ parameters are determined as percentiles of the standard normal distribution, whereas the $T(T-1)/2$ off-diagonal elements of $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ are estimated as polychoric correlations. This is the default parametrization applied in Mplus and lavaan, known as standard delta parameterization. Jöreskog (Reference Jöreskog2001) also refers to it as standard parameterization, whereas Kamata and Bauer (Reference Kamata and Bauer2008) term it marginal parameterization since the marginal distribution of the continuous underlying variables is standardized.

One weakness of standard parameterization is that all underlying variables are standardized to have zero means and unit standard deviations. However, in the context of longitudinal data, where the response alternatives are the same across multiple time points, differences in the distribution of these variables can reflect differences in the means and/or variances of the underlying latent variables.

Researchers have proposed alternative identification constraints for categorical repeated measures, recognizing the limitations of the standard parameterization. These alternative parameterizations, which allow the means and variances of the underlying variables to vary freely, are a significant step toward a more accurate representation of the underlying data structure. Additional constraints are imposed, mainly on the thresholds of these variables, to ensure that the latent propensities are identified and comparable across time points. By imposing these constraints, analysts can establish a fixed reference point for the scale of the underlying response variables. This ensures that observed changes in the distributions reflect differences in the underlying variables, not inconsistencies in measurement.

• Alternative parameterization 1. The first alternative freely estimates the underlying variable’s means and variances on all occasions except on the first one (Millsap & Yun-Tein, Reference Millsap and Yun-Tein2004; Muthén & Asparouhov, Reference Muthén and Asparouhov2002). That is, the underlying variable mean and variance at the first occasion, $\mu _{Y_{1}^{*}}^{alt1}$ and ${\sigma _{Y_{1}^{*}}^{2alt1}}$ , are fixed to zero and one, respectively, while $\mu _{Y_{t}^{*}}^{alt1}$ and $\sigma _{Y_{t}^{*}}^{2alt1}$ are estimated for $t> 1$ .Footnote ³ These assumptions are sufficient to identify all the thresholds on the first occasion. On subsequent occasions, two thresholds are assumed to be time-invariant, e.g., $\tau _{1t}^{alt1}=\tau _{1}^{alt1}$ and $\tau _{2t}^{alt1}=\tau _{2}^{alt1}$ for all t. A notable feature is that complete invariance of all threshold parameters is not required. As discussed before, to estimate differences in the means and variances of latent variables over time—as done with latent growth models—one must ascertain that the underlying variables are on the same scale on different occasions. This is achieved by defining a common across-time metric in terms of the standard deviation on the first occasion.

  As shown in the second row of Table 2, the number of parameters to be estimated is the same as in the standard parameterization, with $T(C-3)$ thresholds, $2(T-1)$ means and variances of the underlying variables, and $T(T-1)/2$ polychoric covariances, for a total of $T(C-1) + T(T-1)/2$ parameters. A one-to-one correspondence exists between the parameters in this alternative parameterization and those from the standard one. That is, on the first occasion,

  (2) $$ \begin{align} \begin{array}{ccc} \tau_{1}^{alt1}=\tau_{11}^{std},&\tau_{2}^{alt1}=\tau_{21}^{std},&\tau_{c1}^{alt1}=\tau_{c1}^{std}, \quad c=3, \ldots C,\\ \end{array}\end{align} $$

  whereas, on subsequent occasions,

  (3) $$ \begin{align} \begin{array}{ccc} \sigma_{Y_{t}^{*}}^{alt1}=\frac{\tau_{21}^{std}-\tau_{11}^{std}}{\tau_{2t}^{std}-\tau_{1t}^{std}},&\mu_{Y_{t}^{*}}^{alt1}=\frac{\tau_{11}^{std}\tau_{2t}^{std}-\tau_{1t}^{std}\tau_{21}^{std}}{\tau_{2t}^{std}-\tau_{1t}^{std}},&\tau_{ct}^{alt1}=\frac{\tau_{11}^{std}(\tau_{2t}^{std}-\tau_{ct}^{std})-\tau_{21}^{std}(\tau_{1t}^{std}-\tau_{ct}^{std})}{\tau_{2t}^{std}-\tau_{1t}^{std}}, \\&&\qquad \qquad\qquad\qquad \qquad\qquad c=3, \ldots C.\\ \end{array}\end{align} $$

  Some authors claim that analysts should also estimate the thresholds but constrain the same threshold to be equal over time (Muthén & Muthén, Reference Muthén and Muthén1998–2017, Example 6.5). That is, $\tau _{ct}^{alt1}=\tau _{c}^{alt1}$ for $c=1, \ldots , C-1$ , and $t=1, \ldots , T$ , and this parameterization is detailed in the third row of Table 2.Footnote ⁴

  If the assumption of threshold invariance holds, then the alternative thresholds $\tau _{c}^{alt1}$ for each category, $c=1, \ldots , C-1$ , should be equal to the corresponding threshold estimated on the first occasion under the standard parameterization ( $\tau _{c1}^{std}$ ), since the mean $\mu _{Y_{1}^{*}}^{alt1}$ and standard deviation $\sigma _{Y_{1}^{*}}^{alt1}$ on the first occasion are set equal to zero and one, respectively. Based on this relationship,

  $$ \begin{align*}\tau_{ct}^{std}=\frac{\tau_{c1}^{std}-\mu_{Y_{t}^{*}}^{alt1}}{\sigma_{Y_{t}^{*}}^{alt1}},\end{align*} $$

  or, equivalently,

  (4) $$ \begin{align}\sigma_{Y_{t}^{*}}^{alt1}\tau_{ct}^{std}=\tau_{c1}^{std}-\mu_{Y_{t}^{*}}^{alt1}.\end{align} $$

  Based on eq. (4), if thresholds are assumed to be equal over time, the ratio of the standard deviation of the underlying variable on two different occasions must be equal to the ratio of the differences between any two thresholds in the standard parameterization on these two occasions:

  $$ \begin{align*}\frac{\sigma_{Y_{t}^{*}}^{alt1}}{\sigma_{Y_{t'}^{*}}^{alt1}}=\frac{\tau_{c't'}^{std}-\tau_{ct'}^{std}}{\tau_{c't}^{std}-\tau_{ct}^{std}},\quad c < c', t \neq t'.\end{align*} $$

  In other words, if thresholds under the alternative 1 parameterization are truly invariant at times t and $t'$ , the ratio of differences between any two standardized thresholds at the two time points should be equal. If any one of the thresholds for a time point is not invariant, then the equality mentioned above should not hold for the thresholds in the standard parameterization for that point. For the alternative 1 parameterization, for any pairs ( $c',c$ ) and ( $c",c'$ ), with $c < c' < c"$ , we also obtain

  (5) $$ \begin{align}\frac{\sigma_{Y_{t}^{*}}^{alt1}}{\sigma_{Y_{1}^{*}}^{alt1}}=\sigma_{Y_{t}^{*}}^{alt1}=\frac{\tau_{c'1}^{std}-\tau_{c1}^{std}}{\tau_{c't}^{std}-\tau_{ct}^{std}}=\frac{\tau_{c"1}^{std}-\tau_{c'1}^{std}}{\tau_{c"t}^{std}-\tau_{c't}^{std}}, \quad \quad \quad t=1, \ldots, T.\end{align} $$

  Therefore, this equation can be considered as a test of the assumption regarding the lack of invariance of at least one of the thresholds for a given time point.

• Alternative parameterization 2. Jöreskog (Reference Jöreskog2001) proposed an alternative specification for longitudinal data. It is based on defining the origin and unit of measurement of $\mathbf {Y}_{i}^{*}$ in terms of thresholds (Bollen & Curran, Reference Bollen and Curran2006; Fisher & Bollen, Reference Fisher and Bollen2020; Mehta et al., Reference Mehta, Neale and Flay2004). The common practice is to fix the distance between the first two thresholds as one unit on the new scale, that is $\tau _{2t}^{alt2} - \tau _{1t}^{alt2} = 1$ , for all $t=1, \dots , T$ . Equivalently, the first threshold could be fixed to zero and the second at one on each occasion, indicating that $Y_{it}^{*}$ on this new scale must be greater than zero for an individual to be in an ordinal category greater than one. This allows us to recover, on each occasion, both the mean and variance of $Y_{it}^{*}$ as well as the other ( $C-3$ ) thresholds in this new measurement scale. Given a number of ordinal categories equal to or greater than three, these conditions are necessary and sufficient to ensure the identification of the model linking $\mathbf {Y}_{i}$ to $\mathbf {Y}_{i}^{*}$ . This set of constraints also represents an alternative but equivalent parameterization of the standard one. Indeed, it can be easily shown that

  (6) $$ \begin{align} \begin{array}{ccc} \sigma_{Y_{t}^{*}}^{alt2}=\frac{1}{\tau_{2t}^{std}-\tau_{1t}^{std}},&\mu_{Y_{t}^{*}}^{alt2}=-\frac{\tau_{1t}^{std}}{\tau_{2t}^{std}-\tau_{1t}^{std}},&\tau_{ct}^{alt2}=\frac{\tau_{ct}^{std}-\tau_{1t}^{std}}{\tau_{2t}^{std}-\tau_{1t}^{std}},\quad c=3, \ldots, C.\\ \end{array} \end{align} $$

  A one-to-one relationship also stands between the two alternative parameterizations, as illustrated in Appendix A. In applying this parameterization in conjunction with a linear growth model for categorical data, Mehta et al. (Reference Mehta, Neale and Flay2004) additionally imposed the assumption of threshold invariance over time. That is, $\tau _{ct}=\tau _{c}$ , for $c=3, \ldots , C$ , and all t. As discussed by Jöreskog (Reference Jöreskog2001) and Mehta et al. (Reference Mehta, Neale and Flay2004), this assumption implies that

  (7) $$ \begin{align} \tau_{c}^{alt2}= \mu_{Y_{t}^{*}}^{alt2}+\sigma_{Y_{t}^{*}}^{alt2}\tau_{ct}^{std}, \quad c =3, \ldots, C-1; t=1, \ldots T, \end{align} $$

  where $\tau _{ct}^{std}$ is the unconstrained threshold for the cth category at the t time point estimated under the standard parameterization. Eq. (7) sets constraints on $\mu _{Y_{t}^{*}}^{alt2}$ and $\sigma _{Y_{t}^{*}}^{alt2}$ because the right-hand side varies with t, whereas the left-hand side does not. If $C \geq 3$ , the common thresholds, $\mu _{Y_{t}^{*}}^{alt2}$ and $\sigma _{Y_{t}^{*}}^{alt2}$ can be estimated from the univariate marginal data of those variables whose thresholds are supposed to be equal. Eq. (7) can be rewritten for any two thresholds c and $c'$ for any $Y_{t}^{*}$ as $\tau _{ct}^{std}=(\tau _{c}^{alt2}-\mu _{Y_{t}^{*}}^{alt2})/\sigma _{Y_{t}^{*}}^{alt2}$ , and $\tau _{c't}^{std}=(\tau _{c'}^{alt2}-\mu _{Y_{t}^{*}}^{alt2})/\sigma _{Y_{t}^{*}}^{alt2},$ such that

  $$ \begin{align*}\sigma_{Y_{t}^{*}}^{alt2}=\frac{\tau_{c'}^{alt2}-\tau_{c}^{alt2}}{\tau_{c't}^{std}-\tau_{ct}^{std}}, \quad c < c'.\end{align*} $$

  It follows that the ratio of standard deviations on two different occasions is equal to the reciprocal of the ratio of the differences between any two standardized thresholds at those time points. For example, there are two independent equations with three thresholds at each time point. That is, for the pairs (2,1) and (3,2), we get

  (8) $$ \begin{align}\frac{\sigma_{Y_{t}^{*}}^{alt2}}{\sigma_{Y_{t'}^{*}}^{alt2}}=\frac{\tau_{2t'}^{std}-\tau_{1t'}^{std}}{\tau_{2t}^{std}-\tau_{1t}^{std}}=\frac{\tau_{3t'}^{std}-\tau_{2t'}^{std}}{\tau_{3t}^{std}-\tau_{2t}^{std}}.\end{align} $$

  If any one threshold for a time point under the alternative 2 parameterization is not invariant, eq. (8) fails to be true for that time point. Hence, this equation evaluates the lack of invariance of at least one of the thresholds at a given time point.

3.2 Binary case

Binary variables ($C=2)$ are special cases that require further comments. Indeed, with dichotomous variables, there are T available proportions, $\pi _{1t}=P(Y_{it}=1), t = 1, \ldots , T$ , to estimate the only threshold $\tau _{1t}$ , the mean and variance of the underlying variables at each occasion.

When the standard parameterization is adopted, the only available threshold is freely estimated on each occasion, whereas the means and variances of the underlying variables are always set to zero and one, respectively. Hence, no different constraints are placed with respect to the categorical case (see the first row of Table 3).

Table 3

Alternative sets of identification constraints for the auxiliary model in presence of binary data

On the other hand, additional restrictions have to be placed for the alternative parameterizations. Millsap and Yun-Tein (Reference Millsap and Yun-Tein2004) suggest freely estimating the means of the underlying variables on each occasion except on the first one, setting the variances of the underlying variables to one on each occasion and estimating the only available threshold under the assumption of time-invariance. These constraints are detailed in the second row of Table 3. Differently, Muthén and Muthén (Reference Muthén and Muthén1998–2017) suggest freely estimating the underlying variable variances on all the occasions except on the first one, where $\sigma _{Y_{i1}^{*}}^{2}$ is fixed to one. They also suggest keeping the only available threshold invariant on each occasion. In this regard, the means of the underlying variables $Y_{it}^{*}$ have to be fixed to zero for the auxiliary model to be identified (see the third row in Table 3).

In the presence of binary data, it is impossible to implement the threshold constraints proposed by Jöreskog (Reference Jöreskog2001), being only one threshold available per occasion. Even when holding thresholds invariant over time, this does not identify both means and variances on all occasions. Jöreskog (Reference Jöreskog2001) proposes to fix the thresholds equal to zero and the variances $\sigma _{Y_{t}^{*}}^{2}$ to one on all occasions, only allowing the means $\mu _{Y_{t}^{*}}$ to vary over time. This is detailed in the last row in Table 3.

3.3 Estimation

SEMs for continuous endogenous variables analyze the mean vector $\boldsymbol \mu _{\mathbf {Y}}$ and covariance matrix $\boldsymbol \Sigma _{\mathbf {Y}\mathbf {Y}}$ of the observed indicators. However, when one or more of the endogenous observed variables are categorical, the analysis shifts to the mean vector $\boldsymbol \mu _{\mathbf {Y}^{*}}$ and covariance matrix $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ corresponding to $\mathbf {Y}_{i}^{*}$ . If consistent estimators for $\boldsymbol \mu _{\mathbf {Y}^{*}}$ and $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ are available, researchers can analyze them similarly to continuous indicators. Consequently, the estimation procedure comprises two distinct steps. The first step obtains consistent estimates of the means $\boldsymbol \mu _{\mathbf {Y}^{*}}$ and covariance matrix $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ for $\mathbf {Y}_{i}^{*}$ . To perform significance testing, analysts also need the asymptotic covariance matrix of the elements in $\hat {\boldsymbol \mu }_{\mathbf {Y}^{*}}$ and $\hat {\boldsymbol \Sigma }_{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ . Once the means, variances, and covariances of $\mathbf {Y}_{i}^{*}$ are in hand, researchers can estimate the parameters of any longitudinal model they apply to $\mathbf {Y}_{i}^{*}$ .

Here, we focus on the first step needed to estimate $\boldsymbol \mu _{\mathbf {Y}^{*}}$ and the unconstrained covariance matrix $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ . As we discussed in Section 3.1, distributional assumptions and identification constraints are necessary for estimation. Assuming a bivariate standard normal distribution for each pair of variables in $\mathbf {Y}_{i}^{*}$ facilitates the estimation of thresholds and polychoric correlations/covariances for noncontinuous variables. Univariate standard normality for each underlying variable $Y_{it}^{*}$ enables threshold estimation as percentiles of the standard normal distribution:

$$ \begin{align*}\tau_{ct}^{std}=\Phi^{-1}\left(\sum_{r=1}^{c}\frac{n_{rt}}{n}\right),\end{align*} $$

where $\Phi ^{-1}$ is the quantile function of the standard normal distribution, $n_{rt}$ is the number of cases in the rth category at time t, and n is the sample size. While univariate margins are valuable for estimating the thresholds in the standard parameterization, bivariate tables are essential for estimating the correlation between $Y_{it}^{*}$ and $Y_{it'}^{*}$ . Following Olsson (Reference Olsson1979), the log-likelihood for the polychoric correlation is given by

$$ \begin{align*}\ln L=A+\sum_{l=1}^{C}\sum_{r=1}^{C}n_{lr}\ln(\pi_{lr}),\end{align*} $$

where C is the number of categories for both $Y_{it}^{*}$ and $Y_{it'}^{*}$ , $n_{lr}$ is the number of cases in the $lr$ th cell of the bivariate table, and A is an irrelevant constant that does not influence the values that maximize the likelihood. Hold the thresholds fixed at the estimates from the univariate margins,

$$ \begin{align*}\pi_{lr}=\Phi_{2}(\hat{\tau}_{lt},\hat{\tau}_{rt'})-\Phi_{2}(\hat{\tau}_{l-1,t}-\hat{\tau}_{rt'})-\Phi_{2}(\hat{\tau}_{lt},\hat{\tau}_{r-1,t'})+\Phi_{2}(\hat{\tau}_{l-1,t}-\hat{\tau}_{r-1,t'}),\end{align*} $$

where $\Phi _{2}(\cdot )$ is the bivariate normal distribution function with correlation $\rho _{Y_{it}^{*}Y_{it'}^{*}}$ . As a result of this conditional estimation procedure, $\hat {\rho }_{Y_{it}^{*}Y_{it'}^{*}}$ , for $t \neq t'$ , is a pseudo-maximum likelihood estimate.

When the means and variances of each $Y_{it}^{*}$ are of interest, it is possible to estimate these quantities rather than constraining them to 0 and 1, respectively. The alternative parameterizations of the auxiliary model allow this option. In the standard parameterization, assuming $Y_{it}^{*} \sim N(0,1)$ at each occasion allows for estimation for the mean (keeping the variance $\sigma _{Y_{t}^{*}}^{2}$ fixed to one) by shifting the distribution of $Y_{it}^{*}$ by a constant ${\mu _{Y_{it}^{*}}}^{alt}$ . This results in the transformed distribution as ${Y_{it}^{*}} \sim N({\mu _{Y_{it}^{*}}}^{alt},1)$ . If, in addition to the mean, the variance of $Y_{it}^{*}$ is estimated, the distribution of $Y_{it}^{*}$ is shifted and scaled, such that $Y_{it}^{*} \sim N({\mu _{Y_{it}^{*}}}^{alt}, {\sigma _{Y_{it}}^{*}}^{alt})$ . Hence, in the alternative parameterization of the auxiliary model, the corresponding thresholds, means, and variances of the underlying variables are based on their relationships with the thresholds of the standard parameterization. These relationships are given in Section 3.1 and detailed in eqs. (2) and (3) for the alternative parameterization proposed by Millsap and Yun-Tein (Reference Millsap and Yun-Tein2004), and in eq. (6) for the parameterization proposed by Jöreskog (Reference Jöreskog2001). If thresholds are assumed to be invariant over time, additional constraints are placed on the means and variances of the underlying variables under the alternative parameterizations, as discussed in Section 3.1.

The polychoric covariance matrix of the underlying variables under the alternative parameterization is then determined by scaling the polychoric correlation matrix $\hat {\mathbf {R}}_{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ estimated for the standard parameterization as follows:

$$ \begin{align*}\hat{\boldsymbol \Sigma}_{\mathbf{Y}^{*}\mathbf{Y}^{*}}=\mathbf{D}^{-1}\hat{\mathbf{R}}_{\mathbf{Y}^{*}\mathbf{Y}^{*}}\mathbf{D}^{-1},\end{align*} $$

where $\mathbf {D}$ is a $T \times T$ diagonal matrix with generic element $1/{\sigma _{Y_{t}^{*}}}^{alt}, t=1, \ldots , T$ .

3.4 Illustrative example: auxiliary model for NLSY97 data

Returning to our motivating NLSY97 example, we estimate the auxiliary measurement model connecting the binary and ordinal variables to their underlying continuous variables. The underlying variables $\mathbf {Y}_{i}^{*}$ represent propensities such as the inclination to use illegal drugs, experience depression, or perceive good/excellent health status. In the case of the illegal drug use variable, the threshold is the point that separates use from nonuse of illegal drugs. An individual whose propensity exceeds the threshold uses illegal drugs, whereas those who fall below it do not. Similarly, for the other two categorical variables, when the latent propensity or perception falls between thresholds $\tau _{ct}$ and $\tau _{c+1,t}$ on a given occasion, the observed ordinal response corresponds to category c.

As discussed in Section 3.1, various parameterizations of the auxiliary model are available. Both the standard auxiliary model and alternative parameterizations, where no threshold invariance constraints are imposed, result in just identified models (with zero degrees of freedom) that perfectly fit the data. When researchers fit these to the data, they obtain estimates of unknown thresholds, means of underlying variables, their variances, and polychoric correlations/covariances. For illustrative purposes, we report—in the Supplementary Material, due to space constraints—the estimated means and polychoric correlations/covariances for the auxiliary model that jointly considers all three observed variables based on the standard and alternative parameterizations. For the latter, we consider the models with and without threshold invariance.

The assumption of threshold invariance is tested by estimating the auxiliary model in which the same threshold at all time points is restricted to be equal. It is worth noting that both alternative (1 or 2) parameterizations of the auxiliary model are equivalent in their chi-square and degrees of freedom. Table 4 provides fit statistics for each categorical variable’s auxiliary model based on the alternative (1 or 2) parameterization with threshold invariance. For the binary variable, threshold invariance is mandatory for the auxiliary model to be just identified and is not testable. Under the assumption of threshold invariance, we have also estimated the multivariate auxiliary model for the three variables, which includes the estimation of thresholds, means, and variances for each variable, along with the polychoric covariances among all the underlying variables. Note that the threshold invariance hypothesis is not rejected for depressive symptoms at any conventional level and that for general health status is marginally significant. With over 5000 cases, statistical power should be high. We also checked the fit indexes and found that the CFI, TLI, and RMSEA all suggest excellent fit. The negative value of the BIC ( $= \chi ^{2} - \ln (n)df$ ) also support the invariant threshold models (Raftery, Reference Raftery1995).

Table 4

Fit statistics for the alternative(1 and 2) auxiliary models based on threshold invariance estimated for depressive symptoms, general health status, and for all the three variables jointly considered

These results suggest that the auxiliary model under alternative 1 or 2 with invariant thresholds is the most promising structure to use for these data when moving to the longitudinal model. However, these results are only with regard to the auxiliary model. In the following Sections, we consider the possible identification constraints of the different longitudinal models, and we turn to this topic next.

4.1 Identification issues

The ALT model for continuous measures requires constraints to ensure its identification, particularly with fewer than five waves of data (Bollen & Curran, Reference Bollen and Curran2004). Extending this model—and its special cases—to ordered categorical measures presents additional challenges. We can approach the identification problem by dividing the constraints needed into two parts. The first part contains constraints to identify $(\boldsymbol \mu _{\mathbf {Y}^{*}}, \boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}})$ . The second part contains conditions to identify the autoregressive latent trajectory model parameters. For the latter, we will also discuss the identification issues for the linear latent growth and autoregressive of order one models.

Once the auxiliary model is accurately identified and estimated, we treat the means and covariances between the underlying variables $\mathbf {Y}_{i}^{*}$ as known and use them to identify the longitudinal model for $\mathbf {Y}_{i}^{*}$ .

This section outlines identification conditions for the linear latent growth, the first-order autoregressive, and ALT models. We explore alternative parameterizations for these models based on the previously discussed auxiliary model specifications. Our assessment focuses on whether distinct specifications lead to equivalent models. In SEMs, equivalent models generate identical model-implied moment matrices and equally fit the data, with equal test statistics, fit indexes, and degrees of freedom (Lee & Hershberger, Reference Lee and Hershberger1990; Stelzl, Reference Stelzl1986). Equivalent models impose the same constraints on the population covariance matrix, known as $\boldsymbol \Sigma $ constraints (Steiger, Reference Steiger2002). Models with the same $\boldsymbol \Sigma $ constraints are $\boldsymbol \Sigma $ -equivalent since they cannot be empirically distinguished.

A formal definition of model equivalence has been provided by Raykov and Penev (Reference Raykov and Penev1999), who established a necessary and sufficient condition based on the existence of a parameter transformation that preserves the implied covariance matrix and covers the entire parameter space of the other model. This implies that equivalent models remain invariant under this transformation. Validating this condition involves deriving implied covariance matrices, equating corresponding elements, and solving for parameters to confirm the $\boldsymbol \Sigma $ -equivalence. If there is no mapping, the two models are not equivalent. An extension of the rule by Raykov and Penev (Reference Raykov and Penev1999) for identifying equivalent models with a mean structure has been discussed by Levy and Hancock (Reference Levy and Hancock2007) and Losardo (Reference Losardo2009) with a specific focus on the equivalence of latent curve model specifications. We follow the approach of Losardo (Reference Losardo2009), accommodating the presence of categorical repeated measures.

4.1.1 Linear growth model

We begin by considering the linear latent growth model for $\mathbf {Y}_{i}^{*}$ . Three different specifications correspond to the standard and alternative parameterizations of the auxiliary model. Despite the time-specific scale of thresholds in the standard parameterization, we consider it being the default option in Mplus 8.6 and lavaan 0.6-16.

To derive identification constraints for each alternative specification, as detailed in Table 5, we need the mean and covariance matrix implied by the linear growth model. These are given by

(13) $$ \begin{align} \boldsymbol \mu_{\mathbf{Y}^{*}}&=\mathbf{B}_{\mathbf{Y}^{*}\boldsymbol \alpha}\boldsymbol \mu_{\boldsymbol \alpha},~~~~~~~~~~~~~~~ \end{align} $$

(14) $$ \begin{align} \boldsymbol \Sigma_{\mathbf{Y}^{*}\mathbf{Y}^{*}}&=\mathbf{B}_{\mathbf{Y}^{*}\boldsymbol \alpha}\boldsymbol \Sigma_{\boldsymbol \alpha\boldsymbol \alpha}\mathbf{B}_{\mathbf{Y}^{*}\boldsymbol \alpha}' + \boldsymbol \Theta_{\boldsymbol \varepsilon}. \end{align} $$

Table 5

Identification constraints and degrees of freedom for the linear growth model for categorical repeated measures

For a deeper understanding of the relationship between these specifications, Table 5 also reports the corresponding degrees of freedom.

The WLSMV estimator only uses information coming from the first- and second-order sample statistics corresponding to the $T(C-1)$ univariate proportions and the $T(T-1)/2$ polychoric correlations. The degrees of freedom of a given specification are the difference in the number of available information, $T(C-1)+T(T-1)/2$ , and the number of parameters in the considered model. It is important to note that not all three specifications have the same degrees of freedom. Hence, the standard parameterization would not provide an equivalent specification to the alternative ones. The assumptions made for the identification of the auxiliary model impact the specification and interpretation of the linear latent growth model defined for $\mathbf {Y}_{i}^{*}$ . The two alternative parameterizations provide equivalent specifications to each other since a one-to-one relationship between their parameters exists, as detailed in Table 6.

Table 6

Parameter transformations for the linear growth models based on the parameterization proposed by Muthén and Muthén (Reference Muthén and Muthén1998–2017) (Model A) and the one illustrated by Mehta et al. (Reference Mehta, Neale and Flay2004) (Model B)

When adopting the standard parametrization (second column of Table 5), Stage 1 constraints are sufficient to identify the auxiliary model. Setting ( $\boldsymbol \mu _{\mathbf {Y}^{*}},\texttt {diag}(\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ )) to ( $\mathbf {0},\mathbf {I}$ ) identifies the thresholds $\tau _{ct}, c=1, \ldots , C-1$ , on all occasions. Stage 2 constraints involve fixing $\boldsymbol \mu _{\boldsymbol \alpha }$ to $\mathbf {0}$ to satisfy the implied mean condition in eq. (13). Appendix B shows that these conditions are sufficient to identify the linear latent growth model., followed by the fact that since the implied mean is a function of the constant matrix $\mathbf {B}_{\mathbf {Y}^{*}\boldsymbol \alpha }$ and of the vector of growth means $\boldsymbol \mu _{\boldsymbol \alpha }$ , the latter has to be fixed equal to zero for the implied mean condition in eq. (13) to be satisfied.

A key distinction between SEMs for continuous outcomes and models for categorical data lies in the identification of the error $\boldsymbol \varepsilon _{i}$ variances, $\boldsymbol \Theta _{\boldsymbol \varepsilon }$ . The variances of the errors cannot be independently identified from the variances of the underlying variables, diag( $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ ), as $\mathbf {Y}_{i}^{*}$ is a latent vector without an inherent scale. If diag( $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}})=\mathbf {I}$ , the variance–covariance matrix of the growth components $\boldsymbol \Sigma _{\boldsymbol \alpha \boldsymbol \alpha }$ is correctly identified, while the error variances $\boldsymbol \Theta _{\boldsymbol \varepsilon }$ are determined as a remainder based on eq. (14).

An alternative set of identification constraints replaces diag( $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}})=\mathbf {I}$ with $\boldsymbol \Theta _{\boldsymbol \varepsilon }=\mathbf {I}$ . It is denoted as (standard) theta parameterization in Mplus 8.6 and lavaan 0.6-16 and termed conditional parameterization by Kamata and Bauer (Reference Kamata and Bauer2008). This approach assumes standardized conditional distributions of $\mathbf {Y}_{i}^{*}$ given $\boldsymbol \eta _{i}$ , with marginal variances of $\mathbf {Y}_{i}^{*}$ obtained as the remainder based on eq. (14). In longitudinal data, this assumption could be more suitable in that the variances of $\mathbf {Y}_{i}^{*}$ are permitted to vary over time while the error variances are not, and this latter assumption might be more plausible. Of course, if there are reasons to think that these error variances differ over time, this assumption is also not desirable.

To address the conflicting conclusions drawn by Grimm and Liu (Reference Grimm and Liu2016), Lee et al. (Reference Lee, Wickrama and O’Neal2018) and Newsom and Smith (Reference Newsom and Smith2020) regarding the relationship between these standard parameterizations of the linear growth model, we provide a theoretical illustration. We limit our analyses to models with four time points to avoid complex and unproductive mathematical details. We compare the two linear growth model parameterizations (based on diag $(\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}})=\mathbf {I}$ or $\boldsymbol \Theta _{\boldsymbol \varepsilon }=\mathbf {I}$ ) in terms of $\boldsymbol \Sigma $ and $\boldsymbol \mu $ constraints. Following Steiger (Reference Steiger2002), the covariance matrix and the mean vector implied by each standard (delta and theta) model specification are first derived and reported in the Supplementary Material. Both model specifications have no $\boldsymbol \mu $ constraints and $T(T+1)/2-3$ (seven) $\boldsymbol \Sigma $ constraints corresponding to the degrees of freedom on the covariance structure of the model, as illustrated in Table 7.

Table 7

$\boldsymbol \Sigma $ constraints for the linear growth model based on the (standard) delta parameterization (diag( $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}})=\mathbf {I}$ ) and the (standard) theta parameterization ( $\boldsymbol \Theta _{\boldsymbol \varepsilon }=\mathbf {I}$ )

It is evident that these two standard specifications of the linear growth model are not empirically equivalent, as the constraints on the underlying variable variances are different, and the $\boldsymbol \Sigma $ constraints derived when all the underlying variable variances are fixed to one are defined in terms of polychoric correlations $\rho _{t,t'}$ .

When adopting the alternative parameterization proposed by Muthén and Asparouhov (Reference Muthén and Asparouhov2002) and by Millsap and Yun-Tein (Reference Millsap and Yun-Tein2004) for the auxiliary model, which fixes the mean of the underlying variable to zero only on the first occasion, the implied mean condition in eq. (13) necessitates constraining only the mean of the intercept component, $\mu _{\alpha _{0}}$ , to zero (refer to the third column of Table 5). The variance of the underlying variable on the first occasion is set to one, while the other ( $T-1$ ) variances are freely estimated. Due to the interdependence of $\sigma _{Y_{t}^{*}}^{2}$ and $\sigma _{\varepsilon _{t}}^{2}$ , an alternative specification is obtained by placing these variance restrictions on the errors rather than on the underlying variable. In Appendix B, we prove that these conditions are sufficient for the model identification.

In contrast to the standard parameterization, it can be demonstrated that placing restrictions on diag( $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ ) or $\boldsymbol \Theta _{\boldsymbol \varepsilon }$ defines two specifications of the linear growth model that are equivalent. Table 8 presents the parameter transformations for the model where all underlying variable variances but one are freely estimated (Model C) and the model where all error variances are free parameters except on the first occasion (Model D).

Table 8

Parameter transformations for the linear growth models where all underlying variable variances are free except on the first occasion (Model C) and when all error variances are free except at the first time point (Model D)

While these transformations are not identity functions, it is evident that parameters in Model C are functions of the corresponding parameters in Model D divided by the underlying variable variance on the first occasion under that specification, that is, ${\sigma _{{Y}^{*}_{1}}^{2D}}=1+ \sigma _{\alpha _{0}}^{2D}$ . Conversely, parameters in the theta specification are functions of the corresponding parameters in the Model C parameterization divided by $1- \sigma _{\alpha _{0}}^{2C}$ , where $\sigma _{\alpha _{0}}^{2C}$ represents the intercept variance in the model where all underlying variable variances are freely estimated except on the first occasion.

Finally, following the alternative 2 parameterization suggested by Jöreskog (Reference Jöreskog2001), where all means and variances of the underlying variables are freely estimated, no restrictions are needed for the identification of the growth parameters, as detailed in the fourth column on Table 5 and proven in Appendix B. For this specification, focusing on either the $\mathbf {Y}_{i}^{*}$ or error variances yields the same $\boldsymbol \Sigma $ and $\boldsymbol \mu $ constraints, as illustrated in Table 9. Hence, these two alternative linear growth model specifications are empirically equivalent, with parameters of the two specifications related via identity transformations.

Table 9

$\boldsymbol {\Sigma }$ and $\boldsymbol \mu $ constraints for the linear growth model based on the alternative 2 delta parameterization (diag $(\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}})$ freely estimated) and the alternative 2 theta parameterization $(\boldsymbol \Theta _{\boldsymbol \varepsilon }$ freely estimated)

Binary data. A special case that requires further comments arises when the measured variables are dichotomous $(C=2)$ . Four different specifications of the linear latent growth model can be derived based on the various parameterizations of the auxiliary model detailed in Table 3. Identification constraints for these specifications are provided in Table 10, along with the corresponding degrees of freedom. The degrees of freedom are computed as the difference between the number of available information, equal to $T+T(T-1)/2$ , and the number of parameters in each model.

Table 10

Identification constraints and degrees of freedom for the linear growth model for binary data

All the constraints in Table 10 are sufficient for identifying the corresponding linear latent growth model. The proof follows the same line as detailed in Appendix B in the presence of categorical data. Attention should be paid to the alternative specification suggested by Muthén and Muthén (Reference Muthén and Muthén1998–2017), whose constraints are reported in the fourth column in Table 10. Due to the dependence of the underlying variable moments $(\boldsymbol \mu _{\mathbf {Y}^{*}}, \boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ ) on the parameters of the growth model $(\boldsymbol \mu _{\boldsymbol \alpha }, \boldsymbol \Sigma _{\boldsymbol \alpha \boldsymbol \alpha })$ , based on eqs. (13) and (14), we can employ fewer Stage 1 constraints than those outlined in Table 10 by taking advantage of the model structure. In this case, the identification problem is treated simultaneously for the thresholds, $(\boldsymbol \mu _{\mathbf {Y}^{*}}, \boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ ), and for the latent growth parameters $(\boldsymbol \mu _{\boldsymbol \alpha }, \boldsymbol \Sigma _{\boldsymbol \alpha \boldsymbol \alpha }$ ). A formal proof is provided in Appendix B.

It is important to note that only the alternative parameterization proposed by Millsap and Yun-Tein (Reference Millsap and Yun-Tein2004) and the one suggested by Jöreskog (Reference Jöreskog2001) share the same number of degrees of freedom. Furthermore, they are equivalent specifications of the linear growth model, since a surjective transformation exists that expresses the parameters of one specification as a function of those of the other and vice versa. All the parameters ( $\sigma _{\alpha _{0}}^{2},\sigma _{\alpha _{0}\alpha _{1}},\sigma _{\alpha _{1}}^{2}, \mu _{\alpha _{1}}$ ) are related by identity functions except for the threshold $\tau _{1}$ in the alternative parameterization by Millsap and Yun-Tein (Reference Millsap and Yun-Tein2004) that is equal to minus the intercept mean $\mu _{\alpha _{0}}$ in the Jöreskog (Reference Jöreskog2001) parameterization and vice versa.

For each parameterization, an alternative set of sufficient conditions can be derived by replacing the constraints on diag( $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ ) with corresponding constraints on the error variances $\boldsymbol \Theta _{\boldsymbol \varepsilon }$ . The results derived for categorical data directly apply when the models are fitted to binary observations: anytime the underlying variable or error variances are all fixed to one, the two alternative parameterizations of the linear growth model are not empirically equivalent. This occurs when the standard parameterization, the alternative one proposed by Millsap and Yun-Tein (Reference Millsap and Yun-Tein2004), and that suggested by Jöreskog (Reference Jöreskog2001) are considered. On the other hand, for the parameterization suggested by Muthén and Muthén (Reference Muthén and Muthén1998–2017), the two parameterizations based on freely estimating all but the first underlying variable or error variances are empirically equivalent. The parameters of the two specifications can be expressed as one function of the other, as detailed in Table 8.

4.1.2 First-order autoregressive model

We now consider alternative parameterizations for the first-order autoregressive model. In situations where the dependent variable is influenced by or influences other dependent variables, assumptions need to be placed on the error variances due to the improper parameter constraints that come into play when focusing on the variances of underlying variables (Muthén and Muthén, Reference Muthén and Muthén1998–2017, pp. 485–486). Specifically, when autoregressive components are present, Mplus exclusively supports theta parameterizations.

The mean and covariance matrix implied by the first-order autoregressive model can be expressed as follows

(15) $$ \begin{align} \boldsymbol \mu_{\mathbf{Y}^{*}}&= (\mathbf{I}_{T}- \mathbf{B}_{\mathbf{Y}^{*}\mathbf{Y}^{*}})^{-1}\nu_{\mathbf{Y}^{*}}, ~~~~~~~~~~~~~~~~~~~~~\end{align} $$

(16) $$ \begin{align} \boldsymbol \Sigma_{\mathbf{Y}^{*}\mathbf{Y}^{*}}&=(\mathbf{I}_{T}- \mathbf{B}_{\mathbf{Y}^{*}\mathbf{Y}^{*}})^{-1}\boldsymbol \Theta_{\boldsymbol \varepsilon}(\mathbf{I}_{T}- \mathbf{B}_{\mathbf{Y}^{*}\mathbf{Y}^{*}})^{-1'}.\end{align} $$

Here, $ \mathbf {B}_{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ contains all autoregressive effects among underlying variables, such that $(\mathbf {I}-\mathbf {B}_{\mathbf {Y}^{*}\mathbf {Y}^{*}})^{-1}$ has the following specific structure

$$ \begin{align*}(\mathbf{I}-\mathbf{B}_{\mathbf{Y}^{*}\mathbf{Y}^{*}})^{-1}=\left[\begin{array}{ccccc}1&0&\ldots&0&0\\\phi_{21}&1&\ldots&0&0\\\phi_{32}\phi_{21}&\phi_{32}&\ldots&0&0\\\vdots&\vdots&\cdots&\vdots&\vdots\\\phi_{T(T-1)} \ldots \phi_{21}&\phi_{T(T-1)} \ldots \phi_{32}&\ldots&\phi_{T(T-1)}&1\end{array}\right].\end{align*} $$

Identification conditions for different autoregressive model specifications, along with corresponding degrees of freedom, are presented in Table 11. Appendix B proves the sufficiency of these constraints for model identification. Degrees of freedom are computed based on the difference between information from the first- and second-order statistics, equal to $T(C-1)+T(T-1)/2$ , and the number of parameters in each specification.

Table 11

Identification constraints and degrees of freedom for the autoregressive model for categorical repeated measures

Coherently with the analysis performed for the linear latent growth model, we still consider alternative parameterizations based on the assumption of threshold invariance. Understanding how constraints on the auxiliary model imply different or equivalent autoregressive model specifications is crucial. Alternative parameterizations endow the model with $(C-4)(T-1)+T(T-1)/2$ degrees of freedom, while the standard (theta) parameterization, which involves fixing the error means and variances on all occasions, has $(T-2)(T-1)/2$ degrees of freedom.

When the assumption of threshold invariance is relaxed, a notable distinction is observed between the linear latent growth and autoregressive models. The standard and alternative autoregressive specifications are equivalent, characterized by $(T-2)(T-1)/2$ degrees of freedom. This is unique to the autoregressive model, as the standard specification of the linear growth model for categorical data still differs from the alternative ones.

In the Supplementary Material, we illustrate the model-implied mean and covariance structures for each autoregressive model, focusing on a simplified scenario with four observed time points. Based on these implied moments, we show that the alternative 1 (Model E) and 2 (Model F) parameterizations, both based on threshold invariance, share the same $\boldsymbol \Sigma $ constraints, represented by $\sigma _{3,2}\sigma _{4,1}=\sigma _{4,2}\sigma _{3,1}$ , $\sigma _{3,2}\sigma _{2,1}=\sigma _{3,1}\sigma _{Y_{2}^{*}}^{2}$ , and $\sigma _{4,3}\sigma _{3,1}=\sigma _{4,1}\sigma _{Y_{3}^{*}}^{2}$ . Table 12 outlines the relationships between parameters of these equivalent autoregressive model specifications.

Table 12

Parameter transformations for the autoregressive models based on the parameterization proposed by Muthén and Muthén (Reference Muthén and Muthén1998–2017) (Model E) and the one illustrated by Mehta et al. (Reference Mehta, Neale and Flay2004) (Model F)

It is evident that Model F parameters are functions of the estimated error variance on the first occasion in Model E, while Model E parameters depend on the distance between the two first thresholds in Model F.

Binary data. The distinct characteristics of various parameterizations of the autoregressive model of order one become more apparent when applied to binary data. In this context, all considered specifications yield the same degrees of freedom, as presented in Table 13. It is straightforward to demonstrate that the parameters of each model can be expressed as functions of the others. To simplify the presentation, Table 14 just illustrates how parameters in the standard parameterization (Model G) relate to those in each of the alternative parameterizations.

Table 13

Identification constraints and degrees of freedom for the autoregressive model for binary data

Table 14

Parameter transformations for the autoregressive model of order one for binary data based on the standard parameterization of the auxiliary model (Model G), following the alternative parameterization proposed by Millsap and Yun-Tein (Reference Millsap and Yun-Tein2004) (Model H), that proposed by Muthén and Muthén (Reference Muthén and Muthén1998–2017) (Model I), and the one proposed by Jöreskog (Reference Jöreskog2001) (Model L)

Autoregressive latent trajectory model

The ALT model integrates components from the previously discussed models, the linear latent growth and the first-order autoregressive model. Including autoregressive effects requires constraints to be specifically placed on the error variances rather than on the underlying variable variances, as discussed earlier.

In eqs. (11) and (12), the mean vector and covariance matrix implied by the ALT model are directly influenced by the autoregressive component. This influence is evident through terms such as $(\mathbf {I}_{T}- \mathbf {B}_{\mathbf {Y}^{*}\mathbf {Y}^{*}})^{-1}\nu _{\mathbf {Y}^{*}}$ and $(\mathbf {I}_{T}- \mathbf {B}_{\mathbf {Y}^{*}\mathbf {Y}^{*}})^{-1}\boldsymbol \Theta _{\boldsymbol \varepsilon }(\mathbf {I}_{T}- \mathbf {B}_{\mathbf {Y}^{*}\mathbf {Y}^{*}})^{-1'},$ respectively.

Additionally, both moments depend on the interaction between the autoregressive and growth components. The Supplementary Material provides a detailed presentation of these moments for various ALT model specifications under the simplified scenario of a stationary autoregressive process with four observed occasions.

Table 15 presents identification conditions for different ALT parameterizations. These conditions are sufficient for identification. They can be easily proven in a manner similar to that for the linear latent growth and autoregressive model in Appendix B. The degrees of freedom for each specification are also reported. Consistent with the previous analysis, we consider alternative parameterizations based on the assumption of threshold invariance.

Table 15

Identification constraints and degrees of freedom for the autoregressive latent trajectory model for categorical repeated measures

As observed for the linear growth model, both alternative parameterizations share the same degrees of freedom, equal to $C(T-1)+\frac {T(T-7)}{2}-4$ , while the standard parameterization has $T(T-3)/2-4$ degrees of freedom. The model requires five occasions for identification without imposing additional parameter constraints beyond those outlined in Table 15.

A notable difference from both the linear latent growth and autoregressive model is that the two alternative parameterizations of the ALT model are not empirically equivalent. They exhibit distinct $\boldsymbol \mu $ and $\boldsymbol \Sigma $ constraints. To illustrate this point clearly, we consider the simplified case of a stationary first-order autoregressive model with five observed occasions. The $\boldsymbol \mu $ constraints for the different autoregressive latent trajectory specifications are illustrated in the following Table 16.

Table 16

$\boldsymbol \mu $ constraints for the different autoregressive latent trajectory model specifications

The three ALT specifications also exhibit different $\boldsymbol \Sigma $ constraints, which are detailed in the Supplementary Material due to space constraints. Attempting to establish a correspondence between the parameters of the two alternative parameterizations, as requested by Raykov and Penev (Reference Raykov and Penev1999), proves unsuccessful when equating corresponding elements in the implied covariances matrices. Consequently, the assumptions made for the identification of the auxiliary model strongly influence the specification and interpretation of the ALT model for $\mathbf {Y}_{i}^{*}$ .

Binary data. Moving on to the analysis of binary data, the ALT model inherits the discrepancies observed in the alternative specifications for categorical data. In this context, none of the proposed parameterizations of the auxiliary model leads to equivalent specifications of the ALT model. Despite the fact that the alternative specifications proposed by Millsap and Yun-Tein (Reference Millsap and Yun-Tein2004) and Jöreskog (Reference Jöreskog2001) yield the same degrees of freedom, as presented in Table 17, the parameters of each model cannot be expressed as functions of the others, as required by Raykov and Penev (Reference Raykov and Penev1999).

Table 17

Identification constraints and degrees of freedom for the autoregressive latent trajectory model for binary data

Table 17 provides sufficient identification conditions for each model, along with corresponding degrees of freedom. Five waves of data are necessary—under each parameterization—to identify the model without placing additional constraints on the model parameters. This is coherent with what was observed in the presence of continuous observations (Bollen & Curran, Reference Bollen and Curran2004). Based on these findings, researchers should carefully consider the implications of selecting specific parameterizations of the autoregressive latent trajectory model, as evidenced by the NLSY97 data analysis, where different choices significantly influenced the interpretation of the results.

4.2 Estimation

The correlation matrix (standard parameterization) or the mean vector and unconstrained covariance matrix (alternative parameterizations) estimated based on the auxiliary model specification—as described in Section 3.3—are utilized to derive point estimates for the parameters of the dynamic model specified for $\mathbf {Y}_{i}^{*}$ . Various estimators apply to estimate the structural parameters in a SEM. We specifically focus on the Diagonally Weighted Least Squares (DWLS) estimator, which is the default choice in Mplus 8.6 and lavaan 0.6-16. Once $\hat {R}_{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ or ( $\hat {\boldsymbol \mu }_{\mathbf {Y}^{*}},\hat {\boldsymbol \Sigma }_{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ ) are available, all parameters of SEM are estimated simultaneously.

To begin, we organize the estimated means $\hat {\boldsymbol \mu }_{\mathbf {Y}^{*}}$ and all the diagonal and below diagonal elements in $\hat {\boldsymbol \Sigma }_{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ (or $\hat {R}_{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ ) into a vector $\hat {\boldsymbol \rho }$ . Similarly, we place the implied moments $\boldsymbol \mu _{Y^{*}}$ and $\boldsymbol \Sigma _{\mathbf {Y}^{*}\mathbf {Y}^{*}}$ in a vector $\boldsymbol \rho (\boldsymbol \theta )$ , where $\boldsymbol \theta $ contains the SEM parameters. The DWLS estimator is then determined by minimizing

$$ \begin{align*}F_{DWLS}=[\hat{\boldsymbol \rho}-\boldsymbol \rho(\boldsymbol \theta)]'\texttt{diag}(\boldsymbol \Sigma_{\hat{\boldsymbol \rho}\hat{\boldsymbol \rho}})^{-1}[\hat{\boldsymbol \rho}-\boldsymbol \rho(\boldsymbol \theta)],\end{align*} $$

where $\boldsymbol \Sigma _{\hat {\boldsymbol \rho }\hat {\boldsymbol \rho }}$ is the asymptotic covariance matrix of $\hat {\boldsymbol \rho }$ whose derivation has been widely detailed in Muthén (Reference Muthén1984), Jöreskog (Reference Jöreskog1994), and Muthén and Satorra (Reference Muthén and Satorra1995).

The parameters $\boldsymbol \theta $ are chosen to minimize the weighted sum of squared deviations of $[\hat {\boldsymbol \rho }-\boldsymbol \rho (\boldsymbol \theta )]$ . The $F_{DWLS}$ is consistent, asymptotically unbiased, and normally distributed (Browne, Reference Browne1984). However, it lacks asymptotic efficiency. Default standard errors are no longer accurate, and goodness of fit tests no longer follow a $\chi ^{2}$ distribution. To address this, robust standard errors are obtained by considering the following sandwich-type asymptotic covariance matrix of the parameter estimates $\hat {\boldsymbol \theta }$ (Muthén et al., Reference Muthén, du Toit and Spisic1997)

$$ \begin{align*}\nonumber \boldsymbol \Sigma_{\hat{\boldsymbol \theta}\hat{\boldsymbol \theta}}&=n^{-1}\left[\left(\frac{\partial \boldsymbol \rho(\boldsymbol \theta)}{\partial \boldsymbol \theta}\right)'\texttt{diag}(\boldsymbol \Sigma_{\hat{\boldsymbol \rho}\hat{\boldsymbol \rho}})^{-1}\left(\frac{\partial \boldsymbol \rho(\boldsymbol \theta)}{\partial \boldsymbol \theta}\right)\right]^{-1}\left(\frac{\partial \boldsymbol \rho(\boldsymbol \theta)}{\partial \boldsymbol \theta}\right)'\texttt{diag}(\boldsymbol \Sigma_{\hat{\boldsymbol \rho}\hat{\boldsymbol \rho}})^{-1}\boldsymbol \Sigma_{\hat{\boldsymbol \rho}\hat{\boldsymbol \rho}}^{-1}\\\nonumber &\quad\times\texttt{diag}(\boldsymbol \Sigma_{\hat{\boldsymbol \rho}\hat{\boldsymbol \rho}})^{-1}\left(\frac{\partial \boldsymbol \rho(\boldsymbol \theta)}{\partial \boldsymbol \theta}\right)\left[\left(\frac{\partial \boldsymbol \rho(\boldsymbol \theta)}{\partial \boldsymbol \theta}\right)'\texttt{diag}(\boldsymbol \Sigma_{\hat{\boldsymbol \rho}\hat{\boldsymbol \rho}})^{-1}\left(\frac{\partial \boldsymbol \rho(\boldsymbol \theta)}{\partial \boldsymbol \theta}\right)\right]^{-1}. \end{align*} $$

The square root of the main diagonal at the estimated parameters represents the robust standard errors of the parameter estimates. Mean and variance-adjusted chi-square statistics have been proposed to approximate the shape of the test statistics to the reference chi-square distribution with the associated degrees of freedom. The WLSMV estimator, which is the default estimator in Mplus 8.6 and lavaan 0.6-16 for models with endogenous categorical variables, relies on the Satterthwaite (Reference Satterthwaite1941) type correction. We refer the reader to Satorra and Bentler (Reference Satorra, Bentler, von Eye and Clogg1994) and Muthén et al. (Reference Muthén, du Toit and Spisic1997) for its detailed description.