Bernoulli_lpmf: Probability parameter is nan, but must be in the interval [0, 1]

Modeling
1

Hello Stan community! I have been struggling with this issue for a few days, and I suspect it may be due to a problem with matrix multiplication:

data {
  int<lower=1> N;
  int<lower=0> D;
  matrix[N,N] T_dis;
  matrix[N,2] d;
  vector[N] s;
  int<lower =0> y[N];
  real tau;
}

parameters {
  real<lower=0> a0;
  real<lower=0> a1;
  simplex[D] theta[N];
  vector[D] G[N];
}

transformed parameters{
  matrix[N, D] Z;
  vector[N] p;
  vector[N] ones_vector = rep_vector(1, N);
  for(n in 1:N){
    Z[n] = softmax((log(theta[n]) + G[n])/tau)';
  }
  for(n in 1:N){
    //p[n] = 100
    for(j in 1:N){
        if(j!=n){
          p[n] += dot_product(Z[n], Z[j])*pow(T_dis[n,j], -2) + (1-dot_product(Z[n],Z[j]))*pow(distance((d'*Z*(Z[j])')/(ones_vector'*Z*(Z[j])'),d[n]),-2);
          }
        }
    p[n] = 1-exp(-(a0+a1*s[n])*p[n]);
  }
}

model {
  a0 ~ gamma(3,1);
  a1 ~ gamma(3,1);
  for(n in 1:N){
    theta[n] ~ uniform(0,1);
    G[n] ~ gumbel(0,1);
    target += bernoulli_lpmf(y[n]|p[n]);
  }
}

The error message looks like:

Chain 1: Rejecting initial value:
Chain 1: Error evaluating the log probability at the initial value.
Chain 1: Exception: bernoulli_lpmf: Probability parameter is nan, but must be in the interval [0, 1] (in ‘anon_model’, line 39, column 4 to column 40)
Chain 1: Initialization between (-2, 2) failed after 100 attempts.
Chain 1: Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.

Any thought or suggestions would be helpful! Thank you!

2

Welcome to the Stan community.

It looks as if you never set p to an initial value. Consequently, your loop adds some value to nothing, i.e. NaN += any number equals NaN. The issue might be resolved just by assigning 0 to all elements of p at the start.

vector[N] p = rep_vector(0, N);

As a side-note for the future, it looks like you have a lot of repeated calculations. You could probably speed up the model a lot by doing them once–once per iteration for calculations that use parameters and once per chain for those that just use data. For example,

transformed data{
    matrix[N,N] T_dist_neg2;

    for(i in 1:N){
        for(j in 1:N){
            T_dist_neg2[i,j] = pow(T_dist[i,j], -2);
        }
    }
}
transformed parameters {
    matrix[N,N] Z_dotprod;

    for(i in 1:N)
        for(j in 1:N){
            Z_dotpred[i,j] = dotp_product(Z[i], Z[j]);
        }
    }

    for(n in 1:N){
        for(j in 1:N){
            if(j != n){
                p[n] += Z_dotpred[n,j] * T_dist_neg2[i,j] ...;
            }
        }
    }
}
2 Likes
3

Hi Simon! Thank you for the answer! Setting initial values for p[n] did the trick, and relocating the matrix computations to the transformed parameters block has noticeably sped up the entire process. Many thanks!

1 Like