How to model individual specific variance of error-term in multilevel model

Modeling
1

Does anyone know how to model an individual specific variance of the error-term in the simple Radon example? So, the model is

y_i = \alpha_{j[i]} + \beta x_{i} + \epsilon_i, with \epsilon_{i} \sim N(0,\sigma^2_{i}). The standard stan code with \epsilon_{i} \sim N(0,\sigma^2) equals:

data {
  int<lower=0> J; 
  int<lower=0> N; 
  int<lower=1,upper=J> county[N];
  vector[N] x;
  vector[N] y;
} 
parameters {
  vector[J] a;
  real b;
  real mu_a;
  real<lower=0,upper=100> sigma_a;
  real<lower=0,upper=100> sigma;
} 
transformed parameters {

  vector[N] y_hat;

  for (i in 1:N)
    y_hat[i] <- a[county[i]] + x[i] * b;
}
model {
  sigma_a ~ uniform(0, 100);
  a ~ normal (mu_a, sigma_a);

  b ~ normal (0, 1);

  sigma ~ uniform(0, 100);
  y ~ normal(y_hat, sigma);
}

My next question would be how to do this in a three-level model, for example \epsilon would then follow \epsilon \sim N(0,\sigma^2_{ij}).

2

Just declare

vector<lower = 0>[N] sigma;

and then the vectorization will happen automatically. Unless there’s very strong hierarchical strucuture, you can’t let sigma vary by item as there’s no information to fit it!

The answer’s the same going to any number of levels.

1 Like
3

@Bob_Carpenter thanks for your reply!

Can you explain a little bit more on what you mean with: