// lda.stan

data {
  int<lower=2> m;
  int<lower=2> V;
  int<lower=1> n;
  int<lower=1> pdot;               // sum of p_i
  array[pdot] int<lower=1,upper=V> x;
  array[pdot] int<lower=1,upper=n> doc;  // doc ID for word l
  vector<lower=0>[V] alpha;
  vector<lower=0>[m] gamma;
}
parameters {
  array[n] simplex[m] w;          // topic dist for document i
  array[m] simplex[V] theta;            // word dist for topic j
}
model {
  for (i in 1:n)
    w[i] ~ dirichlet(gamma);  // prior
  for (j in 1:m)
    theta[j] ~ dirichlet(alpha);     // prior
  for (l in 1:pdot) {
    real A[m];
    for (j in 1:m)
      A[j] = log(w[doc[l], j]) + log(theta[j, x[l]]);
    target += log_sum_exp(A);  // likelihood;
  }
}