stan 中具有单纯形列的矩阵

Question

有没有办法在 Stan 中构造一个包含单纯形列的矩阵？我要构建的模型类似于以下，我的模型算作狄利克雷多项式：

data {
  int g;
  int c;
  int<lower=0> counts[g, c];
}

parameters {
  simplex [g] p;
}

model {
  for (j in 1:c) {
    p ~ dirichlet(rep_vector(1.0, g));
    counts[, j] ~ multinomial(p);
  }
}

但是，我想使用潜在的[g, c] 矩阵来实现类似于以下的分层模型的更多层：

parameters {
  // simplex_matrix would have columns which are each a simplex.
  simplex_matrix[g, c] p;
}
model {
  for (j in 1:c) {
    p[, j] ~ dirichlet(rep_vector(1.0, g));
    counts[, j] ~ multinomial(p[, j]);
  }
}

如果有另一种方法来构造这个潜在变量，那当然也很棒！我对 stan 只实现了几个层次模型并不十分熟悉。

Answer 1

要回答您提出的问题，您可以在 Stan 程序的参数块中声明一个单纯形数组，并使用它们填充一个矩阵。例如，

parameters {
  simplex[g] p[c];
}
model {
  matrix[g, c] col_stochastic_matrix;
  for (i in 1:c) col_stochastic_matrix[,c] = p[c];
}

但是，您实际上不需要在您给出的示例中形成列随机矩阵，因为您可以通过索引一组单纯形来实现多项式狄利克雷模型，例如

data {
  int g;
  int c;
  int<lower=0> counts[g, c];
}
parameters {
  simplex [g] p[c];
}
model {
  for (j in 1:c) {
    p[j] ~ dirichlet(rep_vector(1.0, g));
    counts[, j] ~ multinomial(p[j]);
  }
}

最后，您实际上根本不需要声明一个单纯形数组，因为它们可以从后验分布中积分并在 Stan 程序的生成量块中恢复。详情见wikipedia，但其本质是由这个Stan函数给出的

functions {
  real DM_lpmf(int [] n, vector alpha) {
    int N = sum(n);
    real A = sum(alpha);
    return lgamma(A) - lgamma(N + A) 
           + sum(lgamma(to_vector(n) + alpha)
           - sum(lgamma(alpha));
  }
}