Hi,

I am building a model Y=X\beta + \epsilon, where both X and Y are matrices. I want to model the correlation between the Y.

As you can see, I have some pretty big parameter matrices, so I want to declare them in the model block so that they can be discarded. But I got error, saying “std_normal_lpdf: Random variable[1] is nan, but must not be nan!”

However, if I move the declaration to the transformed parameters blocks, then it runs fine, although very slow.

Can you advise where it went wrong?

Also beside the error message, I feel that there are many places that the code can be improved, can you also advise on that?

Thanks a lot!

```
data {
int<lower=0> nSamples; //number of samples
int<lower=0> nBio; //number of y
int<lower=0> nX; // number of x
matrix[nSamples, nBio] y; // the multivariate outcome matrix
matrix[nSamples, nX] X; // predictor matrix
}
parameters {
matrix[nX, nBio] beta; // betas from N(0,1)
vector<lower=0>[nBio] sigma_eps; // sd of the residual
vector<lower=0>[nBio] tau; // prior scale
cholesky_factor_corr[nBio] L;
}
transformed parameters {
}
model {
matrix[nSamples, nBio] z;
matrix[nSamples, nBio] mu;
mu = X * beta + z * (diag_pre_multiply(tau,L))';
sigma_eps ~ exponential(1);
tau ~ cauchy(0, 2.5);
L ~ lkj_corr_cholesky(1);
for(i in 1:nX){
to_vector(beta[i]) ~ student_t(10, 0, 10);
}
for(i in 1:nSamples){
to_vector(z[i]) ~ std_normal();
}
for(i in 1:nSamples){
for(j in 1:(nBio)){
y[i][j] ~ normal(mu[i][j], sigma_eps[j]);
}
}
}
```