est_mlr.Rmd

---
title: 'Estimation: MLR'
site: workflowr::wflow_site
output:
  workflowr::wflow_html:
    toc: true
editor_options:
  chunk_output_type: console
---

Maximum likelihood with robust standard errors (MLR) is a commonly used estimation method for structural equation models when observed data are continuous.
MLR is an estimation method under normal theory maximum likelihood where the observed data are assumed to follow a multivariate normal distribution.
The *robust* nature of MLR is in order to more accurately estimate standard errors. 
In this study, the standard errors are not directly of interest, so we will focus on the estimation of the fit function under MLR where the resulting $F_{ML}$ is dervied.
A discussion of the standard errors is left a paper forthcoming on parameter estimation and recover in ML-CFA. 

For ML-CFA under MLR estimation, the general idea is to find parameters ($\theta$) that maximize the likelihood function of the observed data given a distributional assumption (nromal theory typically in social science).
In ML-CFA, the model is composed of two major pieces: 1) a model for the pooled within-group covariance matrix ($\Sigma_W$), and 2) a model for the between group covariance matrix ($\Sigma_B$).

As shown in Ryu (2014), Bentler \& Liang (2003), and Liang \& Benlter (2004), the sample estimators for the two population covariance matrices for each group ($j$) are
\[S_{W_j} = {(n_j -1)}^{-1} \sum_i^{n_j} (\mathbf{y}_{ij} - \bar{\mathbf{y}}_{j}){(\mathbf{y}_{ij} - \bar{\mathbf{y}}_{j})}^{\prime}\]
\[S_{B_j} = n_j (\bar{\mathbf{y}}_{j} - \bar{\mathbf{y}}){(\bar{\mathbf{y}}_{j} - \bar{\mathbf{y}})}^{\prime}\]
where,

  * $n_j$ is the sample size of group $j$ ($N = \sum_{\forall j} n_j$);
  * $\mathbf{y}_{ij}$ is the observed vector of responses for individual $i$ in group $j$;
  * $\bar{\mathbf{y}}_{j}$ is the vector of average of observed responses in group $j$;
  * $\bar{\mathbf{y}}$ is the vector of average responses across all groups;
  * $S_{Wj}$ is the within-group covariance matrix of group $j$; and
  * $S_{Bj}$ is the between-group matrix of random effects for covariance matrix for group $j$.

The model-implied covariance matrices $\Sigma_W(\theta), \Sigma_B(\theta),$ and $\Sigma_{gj}(\theta)$ are needed for the MLR fit function.  
We need the $\Sigma_{gj}(\theta)$ model-implied covariance matrix to help identify the average difference between the observed and model-implied covariances, which is defined as
\[\Sigma_{gj}(\theta) = \Sigma_B(\theta) + n_j^{-1}\Sigma_W(\theta)\]
which can be interpreted as the group size weighted deviation of group $j$ from the average group covariance matrix.

The maximum likelihood fit function has been shown to be:

\[F_{ML}=\sum_{j=1}^{J}(n_j-1) \left\lbrace \mathrm{log}\mid \Sigma_{W}(\theta) \mid + tr\left(\Sigma_W^{-1}(\theta)S_{w_j}\right)\right\rbrace + \sum_{j=1}^{J} \left\lbrace\mathrm{log}\mid \Sigma_{gj}(\theta) \mid + tr\left(\Sigma_{gj}^{-1}(\theta)S_{Bj}\right)\right\rbrace\]

where the interested reader is refered to Muthen (1994), Bentler \& Liang (2003), Liang \& Bentler (2004), and Ryu (2014) for more information.


## Other Notes

MLR with continuous data closes matches with we have described above. 
However, when the observed data are categroical, the estimation of the covariance matrices becomes much more computationally difficult.
The additional computational demands comes from the need to perform numerical integration for each latent variable in the model.
So, when the observed data are treated as categorical, for each categorical item an additional layer of integration is required to account for the latent response variable.
Since this very quickly becomes intractable, and other schemes are needed to estimate the model, such as a Monte Carlo approximation or using a different estimation method.
In this study, we would have needed to approximate a 14 dimensionally integration... (not exactly manageable).

Given these constaints on simply the estimation of the polychoric correlations, let along the comutational constraints on the model itself and with categorical data, MLR requires the use of numerical integration for each latent variable in the model.
The assumption that ordered categorical data with at least five response options can be used to drastically spend up convergence.
Therfore, when testing MLR, we assumed the observed ordered categorical data were approximately continuous.
This was also aided by the fact that we simulated these data to be unimodal and approximately symmetric.


# References


1. Muthén, B. O. (1994). Multilevel Covariance Structure Analysis. Sociological Methods & Research, 22(3), 376–398. https://doi.org/10.1177/0049124194022003006

2. Bentler, P. M., & Liang, J. (2003). Two-level mean and covariance structures: Maximum likelihood via EM algorithm. In S. P. Reise & N. Duan (Eds.), Multilevel modeling: Methodological advances, issues, and applications (pp. 53–70). Mahwah, NJ: Lawrence Erlbaum Associates, Inc.

3. Liang, J., & Bentler, P. M. (2004). An EM algorithm for fitting two-level structural equation models. Psychometrika, 69(1), 101–122. https://doi.org/10.1007/BF02295842

4. Ryu, E. (2014). Model fit evaluation in multilevel structural equation models. Frontiers in Psychology, 5. https://doi.org/10.3389/fpsyg.2014.00081