Errors when trying to sample from a hierarchical linear model

Hi - new stan user here.

I’m trying to build a hierarchical linear regression model similar to the one described in section 1.13 of the Stan manual.

Unlike the multivar regression example presented there, my model does not have group-level predictors. Instead, I have a mean \mu for the regression coefficients \beta. In math-speak:

y_n \sim {\rm Normal}(x_n \cdot \beta_{jj[n]}, \sigma) \\ \beta_j \sim {\rm MVNormal}(\mu_j, \Sigma) \\ \mu_j \sim {\rm Normal}(0, 2.) \\ \Sigma := {\rm diag}(\tau) \cdot L_\Omega L_\Omega^T \cdot {\rm diag}(\tau) \\ \tau \sim {\rm Cauchy}(2.5) \\ L_\Omega \sim {\rm LKJCholesky}(3)
// adapted from
// https://mc-stan.org/docs/2_23/stan-users-guide/multivariate-hierarchical-priors-section.html

data {
  int<lower=0> N;              // num samples
  int<lower=1> K;              // num predictors
  int<lower=1> J;              // num groups (stations)
  int<lower=1,upper=J> jj[N];  // group for individual samples
  matrix[N, K] x;              // individual predictors
  vector[N] y;                 // outcomes
}

parameters {
  matrix[J, K] mu;
  matrix[K, J] z;
  cholesky_factor_corr[K] L_Omega;
  vector<lower=0,upper=pi()/2>[K] tau_unif;
  real<lower=0> sigma;                       // prediction error scale
}

transformed parameters {
  matrix[J, K] beta;
  vector<lower=0>[K] tau;     // prior scale
  for (k in 1:K) tau[k] = 2.5 * tan(tau_unif[k]);
  beta = mu + (diag_pre_multiply(tau, L_Omega) * z)';
}

model {
  sigma ~ cauchy(0, 2.5);
  to_vector(z) ~ std_normal();
  L_Omega ~ lkj_corr_cholesky(3);
  to_vector(mu) ~ normal(0, 2.);
  y ~ normal(rows_dot_product(beta[jj], x), sigma);
}

The data for the predictors and target look OK - I’ve cleaned these up and standardized everything before inference:

             count          mean       std       min       25%       50%       75%       max
PREDICTOR1  7300.0  6.793297e-09  1.000068 -4.911924 -0.649995 -0.061916  0.573672  5.153406
PREDICTOR2  7300.0  1.254147e-08  1.000068 -4.760993 -0.647560 -0.059847  0.578347  4.476428
(...etc...)

At runtime, I get a bunch of

Rejecting initial value:
  Log probability evaluates to log(0), i.e. negative infinity.
  Stan can't start sampling from this initial value.

Initialization between (-2, 2) failed after 100 attempts. 
 Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.

The model eventually crashes with:

RuntimeError: Initialization failed.

Clearly my model is mis-specified somehow, but how exactly? And what exactly is initialized between (-2,2)? Thank you!

Full disclosure: I’m using pystan 2.19.1.1.

Hi,

This means the set of initial values used led to an infinite log probability, so sampling couldn’t start. This is sometimes due to initial parameter values being outside reasonable bounds, e.g. a negative initial value for the sigma parameter in a Normal distribution.

By default, Stan takes a random value between -2 and 2 for each parameter at initialization. If this leads to a log probability of Inf, it tries again, with new values between -2 and 2 for a certain number of times, eventually gives up, and outputs the error message you got.

Specifying the bounds of your parameters within the model doesn’t affect which initial value is used by Stan when it initializes, if I’m not mistaken.

You may want to try setting reasonable initial values for your parameters (especially the ones that have boundaries) and see if that solves your problem.

Best,
D

I’ve been staring at the code for quite a while now but I can’t figure out what parameter may be the culprit for the negative infinity log-prob. How can I access the stan initial guesses? I suppose those would tell me what parameters in my model are out of bounds, if any. Thanks!

This is the reference manual for initialization. For bounded parameters, random initial values are in fact constrained appropriately by using the same transform as for the parameter itself.

Maybe you can use the print(); statement to print your parameter values at the first iteration? Otherwise I don’t know how to retrieve initial values.

I would still try setting initial values that are appropriate and see if that fixes the issue.