Clustered data WITHOUT multilevel / GEE model? I have a data-set with around 700 observations from 12 centres. Although the clustering effect as tested in a random intercept model didn't seem significant, it seems more appropriate to use a multilevel / GEE model. But it may not be easy to do so as I want to test for multiple potential mediators. (to be specific, the model in my mind is a SEM without latent variables). 
So my question is: is it theoretically advisable to fit the data in a single-level path model using ML (Huber sandwich estimator, or bootstrap SE)?
UPDATE:
Thanks for those helping out. Just a update for what I did yesterday. I tried the following models:


*

*Single-level SEM in R (lavaan package)

*Single-level SEM in MPlus 

*Two-level SEM in Mplus


(1) and (2) were fitted with ML (Huber / Bootstrap SE). (3) was fitted with Bayesian estimator as MLR did not converge. Basically all of the above gave me similar results but (1) and (2) gave me the chance to do multi-group analysis on top of the SEM and (3) seems to be more theoretically correct. I am not quite sure if the results from (1) and (2) can be trusted. 
A quick literature search told me that Huber SE will give a correct SE even if the model is specified, but the parameters estimated may be biased. Under this situation, will (1) and (2) provide biased parameter estimates?
 A: In my experience, the sensitivity of SEM to detect mediating effects is far outweighed by the cluster level variability of multicentre studies, such as samples within schools, hospitals, or counties. Furthermore, any mediator is usually measurable as a covariate and can be discussed as such. And when that fails, there are way too many assumptions involved with SEM to make me feel confident that it can be used.
With 700 observations, it may be possible to adjust for fixed effects accounting for cluster level variability (11 indicator variables). This depends, of course, upon the number of remaining adjustment variables in the model. Too many will lead to collinearity issues. This is a fully parametric approach
A: GEEs generally require a fairly large number of clusters (e.g. a minimum of around 50), so you would be better off fitting the data in a mixed model with centre as a random effect.  Any potential mediators (assuming by this you mean covariates?) can then be included simply as fixed effects.
If you wish to fit this model in R, you could use the package nlme with the following syntax:
Mod <- lme( Y ~ mediator1 + mediator2, ~1|centre)

and obtain estimates of the effects of the mediators using summary(Mod).
