Inferring very small parameters when implementing non-centered hierarchical model

Modeling
1

Hi!

I’m trying to fit a model to some data in stan. I’ve got some weakly informative priors already based off of some early data analysis and have a good idea of what the parameters in the model should be, but I’m running into some issues with my non-centered hierarchical model. I know exactly where the issue is, but I’m looking for some guidance on how to try and get past it. I imagine I’ll need some sort of transformation or rescaling somewhere but I’m not sure how best to do this without drastically changing the model itself. I’ve recreated the crux of the problem:

data{
    int<lower = 1> N; //Sample Size
    int<lower = 1> n_rep; //Number of replicates 
    array[n_rep, N] real y; // data
    array[N] real x; //Predictor array, all positive
}


parameters{
    real<lower=0> mu_theta; 
    real<lower = 0> tau; 
    array[n_rep] real theta_tilde;
    real<lower = 0> mu_sigma; 
    real<lower = 0> xi; 
    array[n_rep] real sigma_tilde; 
} 

transformed parameters{
    array[n_rep] real theta;
    array[n_rep] real<lower = 0> sigma;

    for (n in 1:n_rep){
        theta[n] = mu_theta + tau * theta_tilde[n];
        sigma[n] = mu_sigma + xi * sigma_tilde[n];
    }
}


model{
    mu_theta ~ normal(0.1,0.1);
    tau ~ normal(0,0.1);
    mu_sigma ~ normal(0,0.1);
    xi ~ normal(0,0.1);

    for (n in 1:n_rep){
        theta_tilde[n] ~ std_normal();
        sigma_tilde[n] ~ std_normal();
        for(j in 1:N){
        y[n,j] ~ normal(theta[n]*x, sigma[n]);
        }
    }
}

I know that the values for theta for each replicate should be around 0.1 and that the variance is relatively small here too. The issue is that during the step in the transformed parameters block, since, as I understand it, theta_tilde should not be constrained to be positive, this often comes out to be a negative value, the location parameter comes out as negative and I get the “location parameter is nan, but should be finite!” error, presumably as a result of the log of this negative value being taken.

I’ve implemented a near identical model for another dataset with larger parameters and had no issues in the past, so I’d be grateful for any input or suggestions!

2

Hi, @sproodle and welcome to the Stan community!

Do you mean more data? Are you sure that the priors you’re laying down now are consistent with the new data? Just check that the posterior’s in the central 99% interval of the prior.

In terms of Stan, you have a problem with sigma, which is required to be positive, but is defined as follows.

real<lower = 0> mu_sigma; 
real<lower = 0> xi; 
array[n_rep] real sigma_tilde; 

sigma[n] = mu_sigma + xi * sigma_tilde[n];

sigma_tilde[n] ~ std_normal();

If you add a <lower=0> constraint on sigma_tilde, you will ensure the constraint is satisfied.

Having said that, it’s easier to just do this on the unconstrained scale, which is much easier for hierarchical models. Using the offset and multiplier gives this a non-centered parameterization on the log scale.

parameters {
  real mu_log_sigma;
  real<lower=0> sigma_log_sigma;
  vector<offset=mu_sigma, multiplier=sigma_sigma>[n_rep] log_sigma;

model {
  vector[n_rep] sigma = exp(log_sigma);
  log_sigma ~ normal(mu_log_sigma, sigma_log_sigma);
  mu_log_sigma ~ normal(...);
  sigma_log_sigma ~ lognormal(...);
1 Like