Error while using categorical ditribution; Number of Categories is 0

[Queodel]

So I am pretty new to Stan. I am posting the code for my model below. The JAGS version of this model runs fine. But I am getting an error when trying to use the categorical distribution in Stan. If anyone can point out where I am going wrong I shall be grateful.

The model

data{
  real<lower=0> rt[1058];
  int<lower=1,upper=2> ttype[1058];
  int<lower=8, upper=20> setSize[1058];
}
parameters{
  real<lower=0> saillance[2];
  real sigma;
  real tempsItems;
  real t0;
}
model{
  // Local variables

  vector[21]K[2,20];
  real Temp[2,20,21,21];
  real P[2,20];
  real nItemsInspect[1058];
  int nII[1058];
  int z[3];
  z[1]<-8;
  z[2]<-14;
  z[3]<-20;
 
  //Prior
for (s in 1:2){
    saillance[s] ~ exponential(1);
    }
sigma ~ uniform(0,1000);
tempsItems ~ normal(88, .01);
t0 ~ uniform(10, 1000);

for (c in 1:2){
    for(N in 1:20) { // 20 is the max number of items in this experiment.
      P[c, N] <- saillance[c] / (N-1+saillance[c]);
    }
  }  

for (c in 1:2){
    for (N in z[1]:z[3]){
      K[c, N, 1]  <- P[c, N];
      for (k in 2:N){
        for (i in 1:(k-1)) {
          Temp[c, N, k, i] <- 1-P[c, N-i+1];
        }
        for (j in k:20){
          Temp[c, N, k, j] <- 1;
        }
        K[c, N, k]  <- P[c, N-k+1] * prod(Temp[c,N,k,]);
      }
      for (l in (N+1):21){
        K[c, N, l]  <- 0;
      }
    }
  }
 
for (i in 1:1058){
      nII[i]  ~  categorical(K[ttype[i],setSize[i]]);
      nItemsInspect[i] <- nII[i] * tempsItems +t0;
      rt[i] ~ normal(nItemsInspect[i], sigma);  //check
    }
  }

and the consequent error which I get is-

Exception thrown at line 56: stan::math::categorical_log: Number of categories is 0, but must be between (1, 21)

line 56 being-

nII[i]  ~  categorical(K[ttype[i],setSize[i]]);


[Bob Carpenter]

That is a confusing error message. We're going to take
a pass at all of our distribution error messages.

> vector[21] K[2,20];
> int nII[1058];

> nII[i] ~ categorical(K[ttype[i], setSize[i]]);

Each nII value must be between 1 and 21, but there was an i for
which nII[i] == 0. See the definition of categorical() in the manual.

- Bob

> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

[Jonah Gabry]

A few other tips on your model: 

• sigma needs to be constrained with <lower=0> since it's a standard deviation parameter. We also don't recommend priors like sigma ~ uniform(0,1000). If you really do want to use it then you need to declare sigma with <lower=0,upper=1000> instead of just <lower=0>, but I would try something without a hard upper bound, e.g. student_t(3, 0, scale), where scale is something reasonable (informative is good if possible!). With a <lower=0> constraint in the parameters block Stan will know to use a half-student_t distribution. 
  
• Similarly, if t0 going to have a uniform(10,1000) prior then it needs <lower=10,upper=1000> when it's declared in the parameters block. But again, unless constraining it to be between 10 and 1000 is a logical necessity, its better to use soft constraints. 
  
• What does tempsItems represent? The prior you have is  tempsItems ~ normal(88, .01), which is very strong. This prior is essentially saying that there's probability 0 it will be less than 87.95 or more than 88.05. Are you sure that makes sense? 
• For efficiency I think you should be able to vectorize many of the calls to the probability distributions in the model block. If speed isn't a problem then don't worry about that, but if it's a concern then you can probably speed things up quite a bit. So, for example, 
  

for (i in 1:1058){

  nItemsInspect[i] <- nII[i] * tempsItems +t0;
  rt[i] ~ normal(nItemsInspect[i], sigma);
}

will be faster as 

for (i in 1:1058){

  nItemsInspect[i] <- nII[i] * tempsItems + t0;
}

rt ~ normal(nItemsInspect, sigma); // vectorized

Jonah

[Queodel]

Bob and Jonah,

Thanks for your replies. So I read that Stan does not support discrete parameters. Now the argument for the categorical function K[ttype[i],setSize[i]] is a vector representing a probability mass function. So will this be a problem. As you said one of the values in nII[i] is ==0, I was trying to print all the values and see where this was going wrong but my console is full off errors and warning and scrolling up does not show any printed values. Thanks again.

[Bob Carpenter]

There's a bug in printing that when there's an exception in
the density, nothing gets printed. So you'll have to comment
things out to see where the problem is. Sorry compilation takes
so long---it doesn't make this easy enough.

Parameters can't be categorical, but we have the categorical
distributions. You declare nII as a local variable, and then
you have to set it before using it and you should be good to go.

If you know one of the nII[i] is 0, that's a bug in the input.

- Bob