如何表示分类预测变量 rstan？答案

【问题标题】：How to represent a categorical predictor rstan?如何表示分类预测变量 rstan？
【发布时间】：2015-05-24 20:50:51
【问题描述】：

格式化分类预测器以在 STAN 中使用的正确方法是什么？我似乎无法将分类预测变量作为正常因子变量输入，那么转换正常分类变量以使 Stan 可以接受的最快方法是什么？

例如，假设我有一个连续预测器和一个分类预测器

your_dataset = data.frame(income = c(62085.59, 60806.33, 60527.27, 67112.64, 57675.92, 58128.44, 60822.47, 55805.80, 63982.99, 64555.45),
country = c("England", "England", "England", "USA", "USA", "USA", "South Africa", "South Africa", "South Africa", "Belgium"))

看起来像这样：

     income      country
1  62085.59      England
2  60806.33      England
3  60527.27      England
4  67112.64          USA
5  57675.92          USA
6  58128.44          USA
7  60822.47 South Africa
8  55805.80 South Africa
9  63982.99 South Africa
10 64555.45      Belgium

我将如何准备将其输入rstan？

【问题讨论】：

标签： r dataframe stan rstan

【解决方案1】：

Stan 只输入实数或整数变量是正确的。在这种情况下，您希望将分类预测变量转换为虚拟变量（可能不包括参考类别）。在 R 中，您可以执行类似的操作

dummy_variables <- model.matrix(~ country, data = your_dataset)

看起来像这样

   (Intercept) countryEngland countrySouth Africa countryUSA
1            1              1                   0          0
2            1              1                   0          0
3            1              1                   0          0
4            1              0                   0          1
5            1              0                   0          1
6            1              0                   0          1
7            1              0                   1          0
8            1              0                   1          0
9            1              0                   1          0
10           1              0                   0          0
attr(,"assign")
[1] 0 1 1 1
attr(,"contrasts")
attr(,"contrasts")$country
[1] "contr.treatment"

但是，如果您在其他一些变量上存在未建模的缺失，则可能无法得出正确数量的观察结果。通过输入整个模型公式，这种方法可以更进一步，例如

X <- model.matrix(outcome ~ predictor1 + predictor2 ..., data = your_dataset)

现在，您有了一个完整的预测变量设计矩阵，可以在带有线性代数的 .stan 程序中使用，例如

data {
  int<lower=1> N;
  int<lower=1> K;
  matrix[N,K]  X;
  vector[N]    y;
}
parameters {
  vector[K] beta;
  real<lower=0> sigma;
}
model {
  y ~ normal(X * beta, sigma); // likelihood
  // priors
}

建议使用设计矩阵，因为它使您的 .stan 程序可重复用于同一模型的不同变体甚至不同数据集。

【讨论】：

您好，这是一种“指标变量”方法。您将如何为“索引变量”执行此操作？指标变量：是k-1个虚拟变量，指标变量：1个k级虚拟变量
@Lefty ：查看替代答案。

【解决方案2】：

另一种方法是使用索引变量，在这种情况下，Stan 程序看起来像

data {
  int<lower = 1> N; // observations
  int<lower = 1> J; // levels
  int<lower = 1, upper = J> x[N];
  vector[N] y;      // outcomes
}
parameters {
  vector[J] beta;
  real<lower = 0> sigma;
}
model {
  y ~ normal(beta[x], sigma); // likelihood
  // priors 
}

你会像这样将数据从 R 传递给 Stan

list(N = nrow(my_dataset),
     J = nlevels(my_dataset$x),
     x = as.integer(my_dataset$x),
     y = my_dataset$y)

【讨论】：