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.