To get myself used to stan, so that I can code up a complex model I’m trying to fit in the future, I’m trying to get a basic linear mixed effects model running in stan. However I’ve run into trouble when I try to extract the random effects that relate to each individual. I’ve included the code so far below, but this will not run yet.

What I am trying to do is create an array of length `N`

of vectors of length `2`

in `b2`

. I then want to extract, using an integer id number held in integer array `idnum`

, the relevant part of `b2`

. The array `idnum`

will be of length `M`

and will take unique values `1...N`

.

When I try to do this using code:

`linpred_y = X_l*beta1 + Z2_l*b2[idnum,];`

I get the following error message (with extra lines deleted):

```
No matches for:
matrix * vector[ ]
...
No matches for:
vector + ill-formed
```

I’d really appreciate any advice how to extract vectors from a vector array, using id values held in an integer array, such that I can multiply a design matrix, and the extracted vector, to make up part of a linear predictor. Code below:

```
functions{
}
data {
//dimensions of data
int<lower=1> M; //number of longitudinal observations
int<lower=1> N; //number of individuals
//id variables
int<lower=1> idnum[M]; //long id variable
//values for priors
real<lower=0> sigmaBeta;
real<lower=0> AU;
real<lower=0> Aeps1;
vector[M] y1; //longitudinal outcome
matrix[M,3] X_l; //design matrix long fixed
matrix[M,2] Z2_l; //design matrix long indrand
}
transformed data{
vector[2] zeroVec;
zeroVec = rep_vector(0,2); //used for mean vector for random effects
}
parameters{
vector[3] beta1; //fixed effects for longitudinal component
vector[2] b2[N]; //individual specific random effects
real<lower=0> sigmaEps1;
corr_matrix[2] Omega_Dl;
vector<lower=0>[2] tau_Dl;
}
transformed parameters{
vector[M] linpred_y;
matrix[2,2] Dl;
linpred_y = X_l*beta1 + Z2_l*b2[idnum,];
Dl = quad_form_diag(Omega_Dl, tau_Dl);
}
model{
//Priors
tau_Dl ~ student_t(3,0,AU);
Omega_Dl~ lkj_corr(1);
beta1 ~ normal(0,sigmaBeta);
sigmaEps1 ~ student_t(3,0,Aeps1);
//Likelihood
y1 ~ normal(linpred_y,sigmaEps1);
b2 ~ multi_normal(zeroVec,Dl);
}
generated quantities{
}
```

Thanks in advance