0

我有一个模型，用于rho从. 每个项目都有一个固定的效果（平均向量），但我还想为每个人添加一个随机效果。我不是 100% 确定如何做到这一点，但我有一个猜测，如果有人能确认或更正它，我将不胜感激。我只需要将最后一行更改为以下内容：N_itemsN_subjectsmu

y[i]' ~ multi_normal(mu + gamma[i],Sigma)

gamma[i]person 的随机效应在哪里i？（此外，在块中声明实数列向量，parameters然后在model block.

顺便说一句，如果有人对提高这个模型的效率有任何建议，我将永远感激不尽。

data {
  int N_subjects;
  int N_items;
  matrix[N_subjects,N_items] y;
}

parameters {
  vector[N_items] mu;
  real<lower=0> sigma;
  real<lower=0,upper=1> rho;
}

transformed parameters {
  cov_matrix[N_items] Sigma;
  for (j in 1:N_items)
    for (k in 1:N_items)
      Sigma[j,k] <- pow(sigma,2)*pow(rho,step(abs(j-k)-0.5));
}

model {
  sigma ~ uniform(0,100);
  rho ~ uniform(0,1);
  for (i in 1:N_items)
    mu[i] ~ normal(0,100);
  for (i in 1:N_subjects)
    y[i]' ~ multi_normal(mu,Sigma);
}
4

2 回答 2

2

在计算效率的问题上，更改 for (i in 1:N_items) mu[i] ~ normal(0,100); 为， mu ~ normal(0,100); 因为这会将 mu 的内存分配数量从 N_items 减少到 1。类似地， 如果您在块 中 声明为数组，则可以 for (i in 1:N_subjects) y[i]' ~ multi_normal(mu,Sigma); 替换 为 另外，如果您先前的信念是均匀分布在 0 和 100 之间，那么在块中声明这些边界并在块 中 注释掉它的先验在计算 上更有效这样做更有效，因为它避免生成大于 100 的建议，这将被自动拒绝并且没有哈密顿蒙特卡罗所需的定义明确的偏导数。 y ~ multi_normal(mu,Sigma); # or mu + gamma yrow_vectordata row_vector[N_items] y[N_subjects]; sigmaparameters real<lower=0,upper=100> sigma; model // sigma ~ uniform(0,100); sigma

于 2015-05-22T15:10:37.220 回答
2

在模型规范问题上，您可以根据需要指定多元正态似然的期望值。在您的情况下，最好先声明y为row_vectorlike 的数组， row_vector[N_items] y[N_subjects]; 然后在parameters块 vector[N_subjects] unit; row_vector[N_items] item; 中，然后在model块中，建立mu以用于可能性 row_vector[N_items] mu[N_subjects]; for (i in 1:N_subjects) mu[i] <- unit[i] + item; y ~ multi_normal(mu, Sigma); 您将需要适当的先验unit和item。

于 2015-05-23T16:37:46.900 回答