Structural Equation Modeling with Small Samples: Test Statistics

Structural equation modeling is a well-known technique for studying relationships among multivariate data. In practice, high dimensional nonnormal data with small to medium sample sizes are very common, and large sample theory, on which almost all modeling statistics are based, cannot be invoked for model evaluation with test statistics. The most natural method for nonnormal data, the asymptotically distribution free procedure, is not defined when the sample size is less than the number of nonduplicated elements in the sample covariance. Since normal theory maximum likelihood estimation remains defined for intermediate to small sample size, it may be invoked but with the probable consequence of distorted performance in model evaluation. This article studies the small sample behavior of several test statistics that are based on maximum likelihood estimator, but are designed to perform better with nonnormal data. We aim to identify statistics that work reasonably well for a range of small sample sizes and distribution conditions. Monte Carlo results indicate that Yuan and Bentler's recently proposed F-statistic performs satisfactorily.     

Simulating Univariate and Multivariate Nonnormal Distributions through the Method of Percentiles

This article derives a standard normal-based power method polynomial transformation for Monte Carlo simulation studies, approximating distributions, and fitting distributions to data based on the method of percentiles. The proposed method is used primarily when (1) conventional (or L) moment-based estimators such as skew (or L-skew) and kurtosis (or L -kurtosis) are unknown or (2) data are unavailable but percentiles are known (e.g., standardized test score reports). The proposed transformation also has the advantage that solutions to polynomial coefficients are available in simple closed form and thus obviates numerical equation solving. A procedure is also described for simulating power method distributions with specified medians, inter-decile ranges, left-right tail-weight ratios (skew function), tail-weight factors (kurtosis function), and Spearman correlations. The Monte Carlo results presented in this study indicate that the estimators based on the method of percentiles are substantially superior to their corresponding conventional product-moment estimators in terms of relative bias. It is also shown that the percentile power method can be modified for generating nonnormal distributions with specified Pearson correlations. An illustration shows the applicability of the percentile power method technique to publicly available statistics from the Idaho state educational assessment.     

Abstract: Evaluation of Test Statistics for Robust Structural Equation Modeling With Nonnormal Missing Data

Traditional structural equation modeling (SEM) techniques have trouble dealing with incomplete and/or nonnormal data that are often encountered in practice. Yuan and Zhang (2011a) developed a two-stage procedure for SEM to handle nonnormal missing data and proposed four test statistics for overall model evaluation. Although these statistics have been shown to work well with complete data, their performance for incomplete data has not been investigated in the context of robust statistics.

Focusing on a linear growth curve model, a systematic simulation study is conducted to evaluate the accuracy of the parameter estimates and the performance of five test statistics including the naive statistic derived from normal distribution based maximum likelihood (ML), the Satorra-Bentler scaled chi-square statistic (RML), the mean- and variance-adjusted chi-square statistic (AML), Yuan-Bentler residual-based test statistic (CRADF), and Yuan-Bentler residual-based F statistic (RF). Data are generated and analyzed in R using the package rsem (Yuan & Zhang, 2011b).

Based on the simulation study, we can observe the following: (a) The traditional normal distribution-based method cannot yield accurate parameter estimates for nonnormal data, whereas the robust method obtains much more accurate model parameter estimates for nonnormal data and performs almost as well as the normal distribution based method for normal distributed data. (b) With the increase of sample size, or the decrease of missing rate or the number of outliers, the parameter estimates are less biased and the empirical distributions of test statistics are closer to their nominal distributions. (c) The ML test statistic does not work well for nonnormal or missing data. (d) For nonnormal complete data, CRADF and RF work relatively better than RML and AML. (e) For missing completely at random (MCAR) missing data, in almost all the cases, RML and AML work better than CRADF and RF. (f) For nonnormal missing at random (MAR) missing data, CRADF and RF work better than AML. (g) The performance of the robust method does not seem to be influenced by the symmetry of outliers.     

Constructing a covariance matrix that yields a specified minimizer and a specified minimum discrepancy function value

A method is presented for constructing a covariance matrix Σ*₀ that is the sum of a matrix Σ(γ₀) that satisfies a specified model and a perturbation matrix,E, such that Σ*₀=Σ(γ0) +E. The perturbation matrix is chosen in such a manner that a class of discrepancy functionsF(Σ*₀, Σ(γ₀)), which includes normal theory maximum likelihood as a special case, has the prespecified parameter value γ₀ as minimizer and a prespecified minimum δ A matrix constructed in this way seems particularly valuable for Monte Carlo experiments as the covariance matrix for a population in which the model does not hold exactly. This may be a more realistic conceptualization in many instances. An example is presented in which this procedure is employed to generate a covariance matrix among nonnormal, ordered categorical variables which is then used to study the performance of a factor analysis estimator. We are grateful to Alexander Shapiro for suggesting the proof of the solution in section 2.     

Structural equation modeling with heavy tailed distributions

Data in social and behavioral sciences typically possess heavy tails. Structural equation modeling is commonly used in analyzing interrelations among variables of such data. Classical methods for structural equation modeling fit a proposed model to the sample covariance matrix, which can lead to very inefficient parameter estimates. By fitting a structural model to a robust covariance matrix for data with heavy tails, one generally gets more efficient parameter estimates. Because many robust procedures are available, we propose using the empirical efficiency of a set of invariant parameter estimates in identifying an optimal robust procedure. Within the class of elliptical distributions, analytical results show that the robust procedure leading to the most efficient parameter estimates also yields a most powerful test statistic. Examples illustrate the merit of the proposed procedure. The relevance of this procedure to data analysis in a broader context is noted. The authors thank the editor, an associate editor and four referees for their constructive comments, which led to an improved version of the paper.     

Contributions to Estimation of Polychoric Correlations

This research concerns the estimation of polychoric correlations in the context of fitting structural equation models to observed ordinal variables by multistage estimation. The first main contribution of this research is to propose and evaluate a Monte Carlo estimator for the asymptotic covariance matrix (ACM) of the polychoric correlation estimates. In multistage estimation, the ACM plays a prominent role, as overall test statistics, derived fit indices, and parameter standard errors all depend on this quantity. The ACM, however, must itself be estimated. Established approaches to estimating the ACM use a sample-based version, which can yield poor estimates with small samples. A simulation study demonstrates that the proposed Monte Carlo estimator can be more efficient than its sample-based counterpart. This leads to better calibration for established test statistics, in particular with small samples. The second main contribution of this research is a further exploration of the consequences of violating the normality assumption for the underlying response variables. We show the consequences depend on the type of nonnormality, and the number and location of thresholds. The simulation study also demonstrates that overall test statistics have little power to detect the studied forms of nonnormality, regardless of the ACM estimator.     

The Nonsingularity of Γ in Covariance Structure Analysis of Nonnormal Data

Covariance structure analysis of nonnormal data is important because in practice all data are nonnormal. When applying covariance structure analysis to nonnormal data, it is generally assumed that the asymptotic covariance matrix Γ for the nonredundant terms in the sample covariance matrix S is nonsingular. It is shown this need not be the case, which raises a question of how restrictive this assumption may be and how difficult it may be to verify it. It is shown that Γ is nonsingular whenever sampling is from a nonsingular distribution, including any distribution defined by a density function. In the discrete case necessary and sufficient conditions are given for the nonsingularity of Γ, and it is shown how to demonstrate Γ is nonsingular with high probability. Thus, the nonsingularity of Γ assumption is mild and one should feel comfortable about making it. These observations also apply to the asymptotic covariance matrix Γ that arises in structural equation modeling.     

Using the expectation maximization algorithm to estimate coefficient alpha for scales with item-level missing data

A 2-step approach for obtaining internal consistency reliability estimates with item-level missing data is outlined. In the 1st step, a covariance matrix and mean vector are obtained using the expectation maximization (EM) algorithm. In the 2nd step, reliability analyses are carried out in the usual fashion using the EM covariance matrix as input. A Monte Carlo simulation examined the impact of 6 variables (scale length, response categories, item correlations, sample size, missing data, and missing data technique) on 3 different outcomes: estimation bias, mean errors, and confidence interval coverage. The 2-step approach using EM consistently yielded the most accurate reliability estimates and produced coverage rates close to the advertised 95% rate. An easy method of implementing the procedure is outlined.     

Identifying Variables Responsible for Data not Missing at Random

When data are not missing at random (NMAR), maximum likelihood (ML) procedure will not generate consistent parameter estimates unless the missing data mechanism is correctly modeled. Understanding NMAR mechanism in a data set would allow one to better use the ML methodology. A survey or questionnaire may contain many items; certain items may be responsible for NMAR values in other items. The paper develops statistical procedures to identify the responsible items. By comparing ML estimates (MLE), statistics are developed to test whether the MLEs are changed when excluding items. The items that cause a significant change of the MLEs are responsible for the NMAR mechanism. Normal distribution is used for obtaining the MLEs; a sandwich-type covariance matrix is used to account for distribution violations. The class of nonnormal distributions within which the procedure is valid is provided. Both saturated and structural models are considered. Effect sizes are also defined and studied. The results indicate that more missing data in a sample does not necessarily imply more significant test statistics due to smaller effect sizes. Knowing the true population means and covariances or the parameter values in structural equation models may not make things easier either. The research was supported by NSF grant DMS04-37167, the James McKeen Cattell Fund.     

Bootstrap standard error and confidence intervals for the correlation corrected for range restriction: a simulation study

The standard Pearson correlation coefficient is a biased estimator of the true population correlation, rho, when the predictor and the criterion are range restricted. To correct the bias, the correlation corrected for range restriction, rc, has been recommended, and a standard formula based on asymptotic results for estimating its standard error is also available. In the present study, the bootstrap standard-error estimate is proposed as an alternative. Monte Carlo simulation studies involving both normal and nonnormal data were conducted to examine the empirical performance of the proposed procedure under different levels of rho, selection ratio, sample size, and truncation types. Results indicated that, with normal data, the bootstrap standard-error estimate is more accurate than the traditional estimate, particularly with small sample size. With nonnormal data, performance of both estimates depends critically on the distribution type. Furthermore, the bootstrap bias-corrected and accelerated interval consistently provided the most accurate coverage probability for rho.     

Evaluating Small Sample Approaches for Model Test Statistics in Structural Equation Modeling

Through Monte Carlo simulation, small sample methods for evaluating overall data-model fit in structural equation modeling were explored. Type I error behavior and power were examined using maximum likelihood (ML), Satorra-Bentler scaled and adjusted (SB; Satorra & Bentler, 1988, 1994), residual-based (Browne, 1984), and asymptotically distribution free (ADF; Browne, 1982, 1984) test statistics. To accommodate small sample sizes the ML and SB statistics were adjusted using a k-factor correction (Bartlett, 1950); the residual-based and ADF statistics were corrected using modified x² and F statistics (Yuan & Bentler, 1998, 1999). Design characteristics include model type and complexity, ratio of sample size to number of estimated parameters, and distributional form. The k-factor-corrected SB scaled test statistic was especially stable at small sample sizes with both normal and nonnormal data. Methodologists are encouraged to investigate its behavior under a wider variety of models and distributional forms.     

Bayesian analysis of two‐level nonlinear structural equation models with continuous and polytomous data

Two‐level structural equation models with mixed continuous and polytomous data and nonlinear structural equations at both the between‐groups and within‐groups levels are important but difficult to deal with. A Bayesian approach is developed for analysing this kind of model. A Markov chain Monte Carlo procedure based on the Gibbs sampler and the Metropolis‐Hasting algorithm is proposed for producing joint Bayesian estimates of the thresholds, structural parameters and latent variables at both levels. Standard errors and highest posterior density intervals are also computed. A procedure for computing Bayes factor, based on the key idea of path sampling, is established for model comparison.     

Information Matrices and Standard Errors for MLEs of Item Parameters in IRT

The paper clarifies the relationship among several information matrices for the maximum likelihood estimates (MLEs) of item parameters. It shows that the process of calculating the observed information matrix also generates a related matrix that is the middle piece of a sandwich-type covariance matrix. Monte Carlo results indicate that standard errors (SEs) based on the observed information matrix are robust to many, but not all, conditions of model/distribution misspecifications. SEs based on the sandwich-type covariance matrix perform most consistently across conditions. Results also suggest that SEs based on other matrices are either not consistent or perform not as robust as those based on the sandwich-type covariance matrix or the observed information matrix.     

Using mixture models with known class membership to address incomplete covariance structures in multiple‐group growth models

Multi‐group latent growth modelling in the structural equation modelling framework has been widely utilized for examining differences in growth trajectories across multiple manifest groups. Despite its usefulness, the traditional maximum likelihood estimation for multi‐group latent growth modelling is not feasible when one of the groups has no response at any given data collection point, or when all participants within a group have the same response at one of the time points. In other words, multi‐group latent growth modelling requires a complete covariance structure for each observed group. The primary purpose of the present study is to show how to circumvent these data problems by developing a simple but creative approach using an existing estimation procedure for growth mixture modelling. A Monte Carlo simulation study was carried out to see whether the modified estimation approach provided tangible results and to see how these results were comparable to the standard multi‐group results. The proposed approach produced results that were valid and reliable under the mentioned problematic data conditions. We also present a real data example and demonstrate that the proposed estimation approach can be used for the chi‐square difference test to check various types of measurement invariance as conducted in a standard multi‐group analysis.     

Robust Structural Equation Modeling with Missing Data and Auxiliary Variables

The paper develops a two-stage robust procedure for structural equation modeling (SEM) and an R package rsem to facilitate the use of the procedure by applied researchers. In the first stage, M-estimates of the saturated mean vector and covariance matrix of all variables are obtained. Those corresponding to the substantive variables are then fitted to the structural model in the second stage. A sandwich-type covariance matrix is used to obtain consistent standard errors (SE) of the structural parameter estimates. Rescaled, adjusted as well as corrected and F-statistics are proposed for overall model evaluation. Using R and EQS, the R package rsem combines the two stages and generates all the test statistics and consistent SEs. Following the robust analysis, multiple model fit indices and standardized solutions are provided in the corresponding output of EQS. An example with open/closed book examination data illustrates the proper use of the package. The method is further applied to the analysis of a data set from the National Longitudinal Survey of Youth 1997 cohort, and results show that the developed procedure not only gives a better endorsement of the substantive models but also yields estimates with uniformly smaller standard errors than the normal-distribution-based maximum likelihood.     

Simulating correlated multivariate nonnormal distributions: Extending the fleishman power method

A procedure for generating multivariate nonnormal distributions is proposed. Our procedure generates average values of intercorrelations much closer to population parameters than competing procedures for skewed and/or heavy tailed distributions and for small sample sizes. Also, it eliminates the necessity of conducting a factorization procedure on the population correlation matrix that underlies the random deviates, and it is simpler to code in a programming language (e.g., FORTRAN). Numerical examples demonstrating the procedures are given. Monte Carlo results indicate our procedure yields excellent agreement between population parameters and average values of intercorrelation, skew, and kurtosis.     

How to Generate Non-normal Data for Simulation of Structural Equation Models

A procedure for generating non-normal data for simulation of structural equation models is proposed. A simple transformation of univariate random variables is used for the generation of data on latent and error variables under some restrictions for the elements of the covariance matrices for these variables. Data on the observed variables is then computed from latent and error variables according to the model. It is shown that by controlling univariate skewness and kurtosis on pre-specified random latent and error variables, observed variables can be made to have a relatively wide range of univariate skewness and kurtosis characteristics according to the pre-specified model. Univariate distributions are used for the generation of data which enables a user to choose from a large number of different distributions. The use of the proposed procedure is illustrated for two different structural equation models and it is shown how PRELIS can be used to generate the data.     

Latent growth curve analysis with dichotomous items: Comparing four approaches

A Monte Carlo study was used to compare four approaches to growth curve analysis of subjects assessed repeatedly with the same set of dichotomous items: A two‐step procedure first estimating latent trait measures using MULTILOG and then using a hierarchical linear model to examine the changing trajectories with the estimated abilities as the outcome variable; a structural equation model using modified weighted least squares (WLSMV) estimation; and two approaches in the framework of multilevel item response models, including a hierarchical generalized linear model using Laplace estimation, and Bayesian analysis using Markov chain Monte Carlo (MCMC). These four methods have similar power in detecting the average linear slope across time. MCMC and Laplace estimates perform relatively better on the bias of the average linear slope and corresponding standard error, as well as the item location parameters. For the variance of the random intercept, and the covariance between the random intercept and slope, all estimates are biased in most conditions. For the random slope variance, only Laplace estimates are unbiased when there are eight time points.     

结构方程模型统计检验力分析：原理与方法

结构方程模型是心理学、管理学、社会学等学科中重要的统计工具之一。然而, 大量使用结构方程模型的研究忽视了对该方法的统计检验力进行必要的分析和报告, 在一定程度上降低了这些研究的结果的证明效力。结构方程模型的统计检验力分析方法主要有Satorra-Saris法、MacCallum法与Monte Carlo法三类。其中Satorra-Saris法适用于备择模型清晰、检验对象相对简单、检验方法基于χ²分布的情形; MacCallum法适用于基于χ²分布的模型拟合检验且备择模型不明的情形; Monte Carlo法适用于检验对象相对复杂、采用模拟或重抽样方法进行检验的情形。在实际应用中, 研究者应当首先判断检验的目的、方法以及是否有明确的备择模型, 并根据这些信息选择具体的分析方法。     

On Nonequivalence of Several Procedures of Structural Equation Modeling

The normal theory based maximum likelihood procedure is widely used in structural equation modeling. Three alternatives are: the normal theory based generalized least squares, the normal theory based iteratively reweighted least squares, and the asymptotically distribution-free procedure. When data are normally distributed and the model structure is correctly specified, the four procedures are asymptotically equivalent. However, this equivalence is often used when models are not correctly specified. This short paper clarifies conditions under which these procedures are not asymptotically equivalent. Analytical results indicate that, when a model is not correct, two factors contribute to the nonequivalence of the different procedures. One is that the estimated covariance matrices by different procedures are different, the other is that they use different scales to measure the distance between the sample covariance matrix and the estimated covariance matrix. The results are illustrated using real as well as simulated data. Implication of the results to model fit indices is also discussed using the comparative fit index as an example. The work described in this paper was supported by a grant from the Research Grants Council of Hong Kong Special Administrative Region (Project No. CUHK 4170/99M) and by NSF grant DMS04-37167.