data { int N; int J; int county[N]; vector[N] y; } parameters { vector[J] a; real mu_a; real sigma_a; real sigma_y; } model { vector[N] y_hat; for (i in 1:N) y_hat[i] = a[county[i]]; mu_a ~ normal(0, 1); a ~ normal (mu_a, sigma_a); y ~ normal(y_hat, sigma_y); }