Gamma distribution samples nan at initialization with apparently reasonable parameters

I’m trying to run a hierarchical model in which duration of a hypothetical process is drawn from a gamma distribution; the parameters of the gamma vary between individuals, and are drawn from a population distribution. The model fails to initialize; the error message points at the line with the gamma sampling, stating that it produces nan. When I print the parameters of the gamma, they look reasonable (mostly between 1 and 10). My code and data are attached.
IMSalience6_demo1.stan (2.8 KB)
IMSalience6Demo.R (572 Bytes)
IMSalience6_DemoData.RData (853.5 KB)

Is any of your data nan/negative/zero?

Looks like the problem is line 103

    t1 ~ gamma(Tshape, Trate);

Here t1 is neither a parameter nor data but a local variable declared in the model block. Despite being called a “sampling statement” this does not draw and assign a new value to the left-hand side. Instead, it computes the likelihood of the lefthand-side. In fact it is basically

   target += gamma_lpdf(t1 | Tshape, Trate);

This works only if t1 is already initialized. Since it’s unknown it must be declared in the parameters block.

parameters {
  ...
  vector<lower=0>[N] t1;
}

Another thing that caught my attention is line 115-116

    logLik += sum(Lik);
    target += logLik;

logLik has not been initialized so incrementing just gives NaN, and if logLik is a cumulative quantity then maybe target += is supposed to be ouside the loop?
Also, it’s strange to go from Lik to logLik without taking a logarithm. If it’s mean to be log(sum(Lik)) then note that log_sum_exp(log(x) + y) is usually more stable than log(sum(x*exp(y))).

1 Like

Thanks a lot, this solved my problem!