Nested Model in STAN?

Question

Say I want to model a random effect at two levels, i.e. I have two levels of nesting: individuals within a parent group and parent groups within a grandparent group. I know how to write a basic model for a single random effect (below) from examples like these but I don't know how to write the equivalent to

lmer(resp ~ (1|a/b), data = DAT)

in lmer.

STAN code for single RE. Question is, how to nest a within a higher level b?

data{
  int<lower=0> N;
  int<lower=0> K;
  matrix[N,K] X;
  vector[N] price;
  int J;
  int<lower=1,upper=J> re[N];
}
parameters{
  vector[J] a;
  real mu_a;
  real tau;
  real<lower=0> sigma_a;
  real<lower=0> sigma;
  vector[K] beta;
}
transformed parameters{
  vector[N] mu_hat;
  for(i in 1:N)
    mu_hat[i] <- a[re[i]];
}
model {
  mu_a ~ normal(0,10);
  tau ~ cauchy(0,5);
  a ~ normal(mu_a,sigma_a);
  for(i in 1:N)
    price[i] ~ normal(X[i]*beta + mu_hat[i], sigma);
}
"

Answer 1

I'm not sure what the a/b notation is in lmer, but if you want nested levels multiple layers deep, then it's easy with a predictor. Say you have an IRT model with students (j in 1:J) nested in schools (school[j] in 1:S) and schools nested in cities (city[s] in 1:C).

[Update 14 April 2017]

You can now vectorize everything. So rather than this:

for (j in 1:J) 
  theta[j] ~ normal(alpha[school[j]], sigma_theta);
for (s in 1:S)
  alpha[s] ~ normal(beta[city[s]], sigma_alpha);
beta ~ normal(0, 5);

you can have

theta ~ normal(alpha[school], sigma_theta);
alpha ~ normal(beta[city], sigma_alpha);
beta ~ normal(0, 5);