6

假设我想在两个级别上对随机效应进行建模，即我有两个级别的嵌套：父组中的个人和祖父组中的父组。我知道如何从这些示例中为单个随机效应（如下）编写基本模型，但我不知道如何编写等价于

lmer(resp ~ (1|a/b), data = DAT)

在 lmer。

单个 RE 的 STAN 代码。问题是，如何嵌套a在更高级别b？

data{
  int<lower=0> N;
  int<lower=0> K;
  matrix[N,K] X;
  vector[N] price;
  int J;
  int<lower=1,upper=J> re[N];
}
parameters{
  vector[J] a;
  real mu_a;
  real tau;
  real<lower=0> sigma_a;
  real<lower=0> sigma;
  vector[K] beta;
}
transformed parameters{
  vector[N] mu_hat;
  for(i in 1:N)
    mu_hat[i] <- a[re[i]];
}
model {
  mu_a ~ normal(0,10);
  tau ~ cauchy(0,5);
  a ~ normal(mu_a,sigma_a);
  for(i in 1:N)
    price[i] ~ normal(X[i]*beta + mu_hat[i], sigma);
}
"
4

2 回答 2

4

我不确定 lmer 中的 a/b 表示法是什么，但如果你想要嵌套层次多层深，那么使用预测器很容易。假设您有一个 IRT 模型，其中学生（j in 1:J）嵌套在学校（school[j] in 1:S）和学校嵌套在城市（city[s] in 1:C）。

[2017 年 4 月 14 日更新]

您现在可以矢量化所有内容。所以而不是这个：

for (j in 1:J) 
  theta[j] ~ normal(alpha[school[j]], sigma_theta);
for (s in 1:S)
  alpha[s] ~ normal(beta[city[s]], sigma_alpha);
beta ~ normal(0, 5);

你可以有

theta ~ normal(alpha[school], sigma_theta);
alpha ~ normal(beta[city], sigma_alpha);
beta ~ normal(0, 5);
于 2015-06-16T19:32:59.310 回答
0

如果您的模型很简单，则该brms软件包值得一看。它将您的公式编译为 stan 并运行模型。它还具有从 lmer 借来的表达语法。我喜欢的是，您可以将模型编译为 stan 文件，然后在已有lmer公式的情况下在其上构建

当然，它还有一个额外的好处（来自 stan），即围绕“固定效应”与“随机效应”估计的令人困惑的差异消失了，并且两者基本上都被估计为具有后验分布的参数。

于 2017-05-10T17:32:37.917 回答