On Identification and Non-normal Simulation in Ordinal Covariance and Item Response Models

Abstract

A standard approach for handling ordinal data in covariance analysis such as structural equation modeling is to assume that the data were produced by discretizing a multivariate normal vector. Recently, concern has been raised that this approach may be less robust to violation of the normality assumption than previously reported. We propose a new perspective for studying the robustness toward distributional misspecification in ordinal models using a class of non-normal ordinal covariance models. We show how to simulate data from such models, and our simulation results indicate that standard methodology is sensitive to violation of normality. This emphasizes the importance of testing distributional assumptions in empirical studies. We include simulation results on the performance of such tests.

Appendices

Appendix A. On Dichotomous Multidimensional IRT Models

For simplicity we limit the discussion to the dichotomous IRT case. We derive a stochastic representation of the IRT model under weak assumptions, which to our knowledge is a new result, and this representation immediately shows that IRT models are of the form of Eq. (1). This representation is then applied to analyze how marginal assumptions (usually called link functions) in the IRT models transfer to the present discussion.

Assumption 1

Consider a random vector $X = (X_1, X_2, \ldots , X_d)'$ where each coordinate takes on the value 0 or 1.

• (1) There is a p-dimensional random vector f which is such that for $i \ne j$ we have that $X_i$ and $X_j$ are independent conditional on f.

• (2) We assume for $i = 1, 2, \ldots , d$ that $\pi _i(f) := P ( X_i = 1 | f) = H(\zeta _i)$ where $\zeta _i$ is a function of f, and H is a CDF with density h with respect to Lebesgue measure.

A standard assumption (Bartholomew et al. 2008) is that $f \sim N(0, I)$ and that

$$\begin{aligned} \zeta _i =\alpha _{i,0} + \sum _{j=1}^d \alpha _{i,j} f_j. \end{aligned}$$

(5)

This implies that $\zeta = (\zeta _1, \zeta _2, \ldots , \zeta _d)' \sim N(\mu , \Sigma )$ for some $\mu$ and $\Sigma$ which are functions of the $(\alpha _{i,j})$ parameters. The link function H is typically assumed to either be the normal CDF, or the logistic CDF.

The $Z = (Z_1, Z_2, \ldots , Z_d)'$ consist of IID random variables with marginal distribution H, and Z is independent from $\zeta = (\zeta _1,\zeta _2, \ldots , \zeta _d)$, where $\zeta$ is defined in Assumption 1 (2). The proof of the following result is given in the online supplementary material.

Proposition 2

A stochastic representation of X fulfilling Assumption 1 is $X = ( I \{ \xi _1 \le 0 \}, \ldots , I \{ \xi _d \le 0 \} )'$ where $\xi = Z - \zeta$.

Since $Z_1, \ldots , Z_d$ are IID and independent to $\zeta$, we have ${\text {Cov}} \, (\xi ) = {\text {Cov}} \, (Z - \zeta ) = {\text {Cov}} \, (Z) + {\text {Cov}} \, (- \zeta ) = \sigma _Z^2 I + {\text {Cov}} \, (\zeta )$ where $\sigma _Z^2 = \text {Var} \,(Z_1)$ and I is the identity matrix. This simple correspondence means that the covariance structure of Z is that of $\zeta$, except for changes in the variances. However, the choice of H influences the marginals of $\xi$, and the mathematical definition of the covariance of $\xi$ depends on both the marginals and the copula of $\xi$. Hence, H plays a major role in the interpretation of the covariance of $\xi$, since it dictates at what “scale” the covariance model is to be interpreted. When $\zeta$ is multivariate normal, $\xi$ will not be multivariate normal unless H is a normal CDF. Indeed, copulas are not preserved under marginal convolution, so that not even the copula of $\xi$ is normal when H is not a normal CDF. This means that when $\zeta$ follows a normal covariance model but when H is not the normal CDF, the resulting IRT model does not follow even a non-normal covariance model (with respect to the covariance model of $\zeta$) as defined in Definition 3, since we there insist that the marginals are standard normal.

Consider the popular choice of H given by the logistic CDF. Then, $\xi$ is not multivariate normal even when $\zeta$ is multivariate normal. Also, $\xi$ will not have a normal copula. While $\xi$ does have the same covariance matrix as $\zeta$, the covariance matrix is given at a scale where the marginals $F_{\xi _i}$ are convolutions between a logistic and a normal distribution. Since the marginal distributions are not identified when observing only copies of X, it seems difficult to interpret what the covariance matrix of $\zeta$ means. If the marginals are transformed to standard normal, one would instead of $\xi$ study the discretize equivalent variable ${\tilde{\xi }} = (\Phi ^ {-1} F_{\xi _1}(\xi _1), \ldots , \Phi ^ {-1} F_{\xi _d}(\xi _d))'$, whose covariance is neither ${\text {Cov}} \, (\xi )$ nor has a simple relation to ${\text {Cov}} \, (\zeta )$. Finally, using arguments given in the upcoming “Appendix B,” a more natural a priori class of marginals for $\xi$ is often normal, and not the convolution of a logistic and a normal.

Appendix B. An a Priori Justification for Marginal Normality of $\xi$ That May be Plausible in Certain Applications

We assume that the continuous discretized vector $\xi$ have a covariance matrix obtained from a SEM model, that is, certain equations among latent variables are to hold, and these equations have error terms that fulfill certain restrictions in terms of correlation. The covariance model $\theta \mapsto \Sigma (\theta )$ for $\xi$ is therefore motivated independently of the distributional class of $\xi$. Now, in many psychometric settings a central limit theorem argument can be used to make an a priori assumption of normality of $\xi$ plausible. Indeed, let us suppose that $\xi$ can be written as a sum of N random vectors $\xi _N^{(1)}, \ldots , \xi _N^{(N)}$. Under mild conditions, the simplest being that $\xi _{N}^{(i)} = \varepsilon _i/\sqrt{N}$ where $\varepsilon _1, \ldots , \varepsilon _N$ are IID random vectors, the multivariate distribution of $\xi$ is close to that of a multivariate normal when $N \rightarrow \infty$ by a central limit theorem. If this approximation is very good, then the normal theory ordinal model in Definition 1 is appropriate. However, the quality of the approximation need not be very good for finite N, especially when the dimensionality d is high, which is the case in many applications of ordinal covariance models: Indeed, consider the ordinal confirmatory factor analysis model underlying many standard measurement instruments in empirical psychology, containing hundreds of items. In these cases, the marginal distributions of $\xi$ may still be close to normal, since each marginal distribution is not affected by the relation between N and d, but the full distribution of $\xi$ may be far from normal. If the marginals are close to normal but the full distribution is not, then the copula of $\xi$ is not close to normal, and we have marginal normality but not joint normality. This may in certain cases make the marginal normality of $\xi$ plausible, while the full copula of $\xi$ is not normal.

Appendix C. Proofs for Sect. 2

Proof of Lemma 1

Self-discretize X, i.e., let ${\tilde{\xi }} = X$ and apply the transformation in Eq. (1). The thresholds can be chosen in such a way that the discretization transformation is the identity transformation. The discretized version of X is then equal to X, which clearly has the same distribution as X, as required by discretize equivalence. $\square$

We need the following preliminary lemma to prove Proposition 1.

Lemma 2

There exists a continuous random vector ${\tilde{\xi }}$ which is discretize equivalent to $\xi$.

Proof

We here only give a compressed version of the argument. The online supplementary material contains a detailed verification of technical details. By Lemma 1, we may without loss of generality assume that $\xi = X$. Define $x_0 = x_1 - 1$, and let ${\mathcal {Q}} = \{ x = \otimes _{l = 1}^d (x_{j_l}, x_{j_l+1}] : j_l \in \{0, 1, \ldots , K -1 \} \text { for } l = 1, 2, \ldots , d \}$ contain the hyper-rectangles contained between the points of the support $S_X^d$ of X. Let $Q_1, Q_2, \ldots , Q_N$ be the sets in ${\mathcal {Q}}$, and note that they are disjoint. We now define a density ${{\tilde{f}}}$, which smears the probability that X is in $Q_i$ uniformly over each $Q_i$, i.e., we let ${{\tilde{f}}}(x) = \sum _{i=1}^N \frac{P ( X \in Q_i )}{V_i} I \{ x \in Q_i \},$ where $V_i = \int _{{\mathbb {R}}^d} I \{ x \in Q_i \} \, \mathrm{d}x \ne 0$ for $i = 1, 2, \ldots , N$, and $I \{ A \}$ is the indicator function of A, which is one if A is true and zero otherwise. Let $\tilde{\xi }$ have ${{\tilde{f}}}$ as density. Then, ${\tilde{\xi }}$ has the same probability as $\xi$ (i.e., X) of being within the thresholds defined by the limits of the rectangles in $Q_k$ for $k = 1, 2, \ldots , N$, completing the proof. $\square$

Proof of Proposition 1

By Lemma 2, we may assume that $\xi$ is a continuous random vector. This implies that any marginal cumulative distribution function $F_i$ is continuous and increasing. Since it is illustrative, we here give a proof that assumes that $F_i$ is also strictly increasing. A proof of this special case is also given in Almeida and Mouchart (2014) (see their Eq. (12)), and our argument follows closely Section 3 in Grønneberg and Foldnes (2019). The general case, which appears to be new, is proved in the online supplementary material. Since $F_i(\xi _i)$ is uniform on [0, 1] we have that $\Phi ^ {-1} (F_i(\xi _i))$ is standard normal, where $\Phi ^ {-1}$ is the quantile function of the standard normal distribution. Since both $F_i$ and $\Phi ^ {-1}$ are strictly increasing, so is $\Phi ^ {-1} \circ F_i$. For each coordinate $X_i$ of X, we may therefore apply $\Phi ^ {-1} \circ F_i$ to each part of the inequalities defining $X_i$, and get $X_i = \sum _{j=1}^K x_j I \{ \tau _{i,j-1}< \xi _i \le \tau _{i,j} \} = \sum _{j=1}^K x_j I \{ \Phi ^ {-1} (F_i(\tau _{i,j-1}))< \Phi ^ {-1} (F_i(\xi _i)) \le \Phi ^ {-1} (F_i(\tau _{i,j})) \} = \sum _{j=1}^K x_j I \{ \tilde{\tau }_{i,j-1} < {\tilde{\xi }}_i \le {\tilde{\tau }}_{i,j} \}$

where ${\tilde{\tau }}_{i,j-1} = \Phi ^ {-1} (F_i(\tau _{i,j-1}))$, $\tilde{\xi }_i = \Phi ^ {-1} (F_i(\xi _i))$ and ${\tilde{\tau }}_{i,j} = \Phi ^ {-1} (F_i(\tau _{i,j}))$. $\square$