Random variable is nan: Poisson-gamma model for grouped data

Modeling
1

Hello everyone,

I am trying to build a model in Rstan for grouped data, but the sampler does not work because of the error: “gamma_lpdf: Random variable is nan, but must not be nan!”. The Stan code is as the following:

data {

  int<lower=1> I; // num of rows or groups of data

  int<lower=1> J; //num of data point in each group 

  int<lower=1> Y[I,J]; //the response variable
  
  // the initial data do not have equal size subset, use interpolation to generate more data points in smaller group.
}
parameters {

 real<lower=0> alpha[J];

  real<lower=0> beta[J];
} 

model {

   matrix[I,J] lambda; 

   for(i in 1:I){

     for (j in 1:J){

         target +=cauchy_lpdf(alpha[j]|0,5);// weakly informative

         target +=cauchy_lpdf(beta[j]|0,5);

         target +=gamma_lpdf(lambda[i,j]|alpha[j],beta[j]);

         target +=poisson_lpmf(Y[i,j]|lambda[i,j]);    //line 26
    }
   }
}

My data come in this way: 5 groups of data(J=5) with 10 data points in each group.

782 777 787 844 498
96  64  90 109  66
144 133 189 189 124
117  92 181 138  96
170 138 221 324 151
239 184 315 394 189
13   9  16  13   1
37  32  34  54  18
489 408 543 688 351
223 162 268 328 185

R code:

Y.c=c(782,	96,	144,	117,	170,	239,	13,	37,	489,	223,
      777,	64,	133,	92,	138,	184,	9,	32,	408,	162,
      787,	90,	189,	181,	221,	315,	16,	34,	543,	268,
      844,	109,	189,	138,	324,	394,	13,	54,	688,	328,
      498,	66,	124,	96,	151,	189,	1,	18,	351,	185)
Poisson_HB_kernel_dat=list(Y=array(Y.c,dim=c(10,5)),I=10,J=5)
n_sam=2000
n_warmup=n_sam/2
n_chain=4

fit=stan(file ='PoissonBayesiankernel_groupedData.stan', data = Poisson_HB_kernel_dat,iter =(n_sam+n_warmup), warmup=n_warmup,chains = n_chain)

When I ran the code, the sampler didn’t work and the error messages were:

Rejecting initial value:
Error evaluating the log probability at the initial value.
_**Exception: gamma_lpdf: Random variable is nan, but must not be nan!  (in 'model1dfc6d1b537f_PoissonBayesiankernel_groupedData' at line 26)   **_

Initialization between (-2, 2) failed after 100 attempts. 
 Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.
[1] "Error in sampler$call_sampler(args_list[[i]]) : Initialization failed."
[1] "error occurred during calling the sampler; sampling not done"

Line 26 of the code is

target +=poisson_lpmf(Y[i,j]|lambda[i,j]); //line 26

Could anyone suggest a possible cause for the problem “Random variable is nan”? I read one similar post, which can be found here (http://discourse.mc-stan.org/t/random-variable-is-nan-but-must-not-be-nan/1536), but I don’t think the reason is the same.

Thank you.

Julian

2

Variables in Stan are default NaNs unless you set them to anything else, so in the block:

model {
  matrix[I,J] lambda; // <-- lambda will be NaNs unless you set them equal to something
  ...
}

target +=gamma_lpdf(lambda[i,j]|alpha[j],beta[j]); or equivalently lambda[i, j] ~ gamma(alpha[j],beta[j]); don’t do any assignment to lambda.

Do you mean for lambda to be a parameter that gets sampled? If so, define it in the parameters block.

3

Thank you!

Yes, lambda is the parameter of interest. If I move it back to the parameter block, I will receive error massage:

Rejecting initial value:
Error evaluating the log probability at the initial value.
Exception: poisson_lpmf: Rate parameter is -1.5266, but must be >= 0!  (in 'model1dfc3c1f37b1_PoissonBayesiankernel_groupedData' at line 36)

Rejecting initial value:
Error evaluating the log probability at the initial value.
Exception: poisson_lpmf: Rate parameter is -0.0981493, but must be >= 0!  (in 'model1dfc3c1f37b1_PoissonBayesiankernel_groupedData' at line 36)


Initialization between (-2, 2) failed after 100 attempts. 
 Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.
[1] "Error in sampler$call_sampler(args_list[[i]]) : Initialization failed."
[1] "error occurred during calling the sampler; sampling not done"

I tried to bound lambda using

  if(lambda[i,j]>1e-6)
  {
   lambda[i,j] = lambda[i,j];
  }
 else
  {
   lambda[i,j]=1e-6;
  }

but this operation is not allowed for parameter.

4

If you want to put a lower bound on a parameter, you do that when you define it. Something like

matrix<lower = 0.0>[I, J] lambda;

should work.

5

This is the relevant line. You are using lambda before defining its contents. See the first entry of

6

Thanks for your help! I apologize for trying to bound the parameter in such an awkward way without noticing it.

7

No worries. Just keep in mind that you can’t assign to parameters in Stan, it’s a common source of errors.