Help from converting from Mplus code to Stan code


I have the following model in Mplus:


A BY A1-A3;
B BY B1-B3;

I am trying to write out this model in Stan, but I can’t seem to find a way to do so. Does anyone have any guidance?

I don’t know off-hand, but let me know if this is the correct translation of the MPlus code in case it helps others.

You have six observed variables

  • A1, A2, and A3 are all categorical
  • B1, B2, and B3 are all continuous

You want to estimate two continuous latent variables

  • A predicts the categorical variables
  • B predicts the continuous variables

Some of the observed variables are also correlated

  • B1 and A2 are correlated
  • B2 and A3 are correlated

That is indeed the correct description of what the code does.

One possibility would be to use blavaan for this, which uses Stan under the hood. You might be able to use mplus2lavaan() from the lavaan package to translate your model to lavaan, then pass the lavaan model to blavaan.

Or, if you wanted to do it directly in blavaan: assuming that A1-A3 are either binary or ordinal (as opposed to >2 unordered categories, which are not supported), you could do something like

m1 <- ' A =~ A1 + A2 + A3
        B =~ B1 + B2 + B3
        B1 ~~ A2
        B2 ~~ A3 '
fit <- bcfa(m1, data = mydata, ordered = paste0("A", 1:3))

Possibly with an extra line to specify that the correlation between A and B is 0 (if desired, or it is possibly required to identify the model).

The underlying Stan code is general and not tailored to your model, so not best for learning Stan. It would be worth looking at brms as well, see this thread.

1 Like

Thank you for the suggestion. I have attempted to use blavaan for this problem before, however, my parameter estimates were biased. I wanted to have a Stan code of my own to have more freedom with attempts to reparameterize the model.

I would be interested to see the evidence of parameter bias in blavaan.

For coding this model yourself, you might look at the multivariate probit example at the bottom of this page. It does the Albert-Chib thing, where you end up with z parameters that are the continuous responses underlying ordinal variables. Once you have those z parameters, you could model them jointly with the continuous data in the usual way.

That example involves binary data, and it is not trivial to extend to variables with >2 ordered categories. I could add more detail there if you need it.