So I’ve made some changes but I’m still running into problems with passing the variables into the mult_norm_par function. My code now looks like this:

```
functions {
vector mult_norm_par(vector mu_L, vector theta, real[] x_r, int[] x_i) {
int NvarsY = x_i[2];
vector[NvarsY] mu = mu_L[1];
matrix[NvarsY, NvarsY] L_Sigma = mu_L[2];
real lp = multi_normal_cholesky_lpdf(x_r | mu, L_Sigma);
return [lp]';
}
}
data{
int<lower=0> NvarsX; // num independent variables
int<lower=0> NvarsY; // num dependent variables
int<lower=0> N; //Total number of rows
matrix[N, NvarsX] x;
vector[NvarsY] y [N];
}
parameters{
cholesky_factor_corr[NvarsY] L_Omega; //For the correlations between outcome variables
vector<lower=0>[NvarsY] L_sigma; //Residual SD of outcome variables
vector[NvarsY] Beta0; // Intercepts
vector[NvarsX] Beta [NvarsY]; // Slopes
}
transformed parameters{
vector[NvarsY] mu[N];
matrix[NvarsY, NvarsY] L_Sigma;
L_Sigma = diag_pre_multiply(L_sigma, L_Omega);
for (i in 1:N){
for (dv in 1:NvarsY){
mu[i, dv] = Beta0[dv] + dot_product(to_vector(Beta[dv, 1:NvarsX]), x[i, 1:NvarsX]);
}
}
// vars for map_rect
vector theta[N];
real x_r[N, NvarsX];
int x_i[N, 2];
{
for (n in 1:N) {
x_r[n] = to_array_1d(y[n]);
x_i[n, 1] = n;
x_i[n, 2] = NvarsY;
}
}
}
model{
//Priors
L_Omega ~ lkj_corr_cholesky(1);
L_sigma ~ cauchy(0,5);
//for (i in 1:N){
// y[i, 1:NvarsY] ~ multi_normal_cholesky(mu[i, 1:NvarsY], L_Sigma);
//}
target +=sum(map_rect( mult_norm_par, append_row(to_vector(mu[, 1:NvarsY]), L_Sigma)), theta, x_r, x_i) );
}
generated quantities{
matrix[NvarsY, NvarsY] Omega;
matrix[NvarsY, NvarsY] Sigma;
Omega = multiply_lower_tri_self_transpose(L_Omega);
Sigma = quad_form_diag(L_Omega, L_sigma);
}
```

When I run this I’m getting the error:

variable definition base type mismatch, variable declared as base type: vector variable definition has base: real error in ‘examples/map_rect/linear_multivariate_map_rect.stan’ at line 4, column 36

```
2: vector mult_norm_par(vector mu_L, vector theta, real[] x_r, int[] x_i) {
3: int NvarsY = x_i[2];
4: vector[NvarsY] mu = mu_L[1];
^
5: matrix[NvarsY, NvarsY] L_Sigma = mu_L[2];
```

I’'ve tried a few variations, but each time I’m running into a problem unpacking elements of `mu_L`