Modeling Multivariate Multilevel Model in STAN

Hello, I am a beginner of STAN and am reading User’s Guide and Reference Manual now.
I would like to ask any opinion or thoughts to extend the multivariate outcomes model at Chapter 9 to multivariate ‘multilevel’ model.

For example, I had a data set which had two scores (math and science) at level-1 nested students at level-2, middle school at level-3, and high school at level-4.
Because I will not add error term at level-1 because the level-1 error is omitted from the model because level one is only used to partition the multivariate outcomes’ structure into the model (Hox, 2010; Goldstein, 2010).
Therefore, it will be actually 3-level multivariate model as shown below.
Y = (Gamma1 + u11 + u21 + r1)*dummy1 + (Gamma2 + u12 + u22 + r2)*dummy2
(Actually, the middle schools are not purely nested to high school thus it is cross-classified, but as I search in the forum, it would not require specific codes to set the characteristic. If I am wrong, please let me know.)

However, I am struggle with extending the multivariate model to multivariate multilevel model because I cannot find the example to make 3-level multilevel model with multiple responses…
Could you give any your thoughts or suggestion about it? Or, is there any study or example I can follow to understand how the code can be set up to multivariate multilevel model?

Here is the code what I saw at the Manual (Chapter 9.15 Multivariate outcomes)

data {
int<lower=1> K;
int<lower=1> J;
int<lower=0> N;
vector[J] x[N];
vector[K] y[N];
parameters {
matrix[K, J] beta;
cholesky_factor_corr[K] L_Omega;
vector<lower=0>[K] L_sigma;
model {
vector[K] mu[N];
matrix[K, K] L_Sigma;
for (n in 1:N)
mu[n] = beta * x[n];
L_Sigma = diag_pre_multiply(L_sigma, L_Omega);
to_vector(beta) ~ normal(0, 5);
L_Omega ~ lkj_corr_cholesky(4);
L_sigma ~ cauchy(0, 2.5);
y ~ multi_normal_cholesky(mu, L_Sigma);

Thank you for your time in advance.

If you’re a beginner try to start with a univariate model. The majority of what’s hard about a multivariate hierarchical model is the level of bookkeeping required so you don’t want to be doing that at the same time as you construct the hierarchical piece.

It’s elephants all the way down. Really, it’s just more of the same. You just code the multiple levels up directly; each level can involve regression predictors, etc. The tricky part’s just keeping track of all the predictors and when to apply them.

Any examples of a univariate multi-level model?


The documentation has an example of hierarchical logistic regression:

What’s also helpful is the case studies:

In particular, this one may be helpful, for multilevel (hierarchical) modeling: