Transforming a multinomial model into a dirichlet-multinomial

Err, shouldn’t that be multinomial(trans_output[i])? And actually, how do the shapes make sense here, why does first loop have 1:L when output is of length S and the outermost dimension of trans_output is G?

Dirichlet-multinomial has a closed-form PMF. Maybe you could use that.

functions {
  real dirichlet_multinomial_lpmf(int[] y, vector alpha) {
    real sum_alpha = sum(alpha);
    return lgamma(sum_alpha) - lgamma(sum(y) + sum_alpha)
           // + lgamma(sum(y)+1) - sum(lgamma(to_vector(y)+1) // constant, may omit
           + sum(lgamma(to_vector(y) + alpha)) - sum(lgamma(alpha));
  }
}
parameters {
  real<lower=0> theta[P]; // parameters for the ODE
  real<lower=0> phi;
}
model {
  ...
  for(i in 1:S) cat_data[i] ~ dirichlet_multinomial(phi*output[i]);
}
3 Likes