问题描述
在我的模型中，我有 10 个选项（从 1 到 10）为每个主题和每个试验选择 (@H_404_1@expectancy)。我根据下图中的规则计算了每个选项的值，因此每个选项的值根据每次试验中 @H_404_1@shock 和 @H_404_1@v 之间的差异更新（乘以 @H_404_1@alpha） .然后，我使用 softmax 规则将每个选项的 @H_404_1@v 转换为具有相同函数的特定概率：JAGS errors: "Resolving undeclared variables" and "Invalid vector argument to exp"。
我想这里的问题是我不能让 jags 更新相同选择的值。
数据：@H_404_1@expectancy = 每次试验中 1-10 的数字。 @H_404_1@shock=每次试验中编号为 1 或 0。 （我在下面提供了示例数据）
第二个情节是如何在有 2 个选择/1 个主题的情况下完成此操作。
@H_404_1@RW_model <- function(){ # data for(i in 1:nsubjects) # for each person { # initial value for v v [i,1,expectancy[i,1]] <- 0 for (j in 2:ntrials) # for each trial { # expectancy chosen expectancy[i,j] ~ dcat(mu[i,j,1:10]) predk[i,1:10]) # softmax rule to calculate values of each expectancy for each subject # tau is the value sensitivity parameter mu[i,1:10] <- exp_v[i,1:10] / sum(exp_v[i,1:10]) exp_v[i,j-1]] <- exp(v[i,j-1]]/tau[i]) # prediction error: difference between Feedback and learned values of the chosen expectancy pe [i,j-1] <- shock [i,j-1] - v [i,j-1,j-1]] # value updating process for expectancy v [i,j-1]] <- v [i,j-1]] + alpha [i] * pe [i,j-1] } } # priors for (i in 1:nsubjects){ tau [i] ~ dunif (0,3) alpha [i] ~ dunif (0,1) } } # example data/ initial value/ parameters nsubjects <- 42 ntrials <- 14 shock <- matrix(c(0,0),nrow=42,ncol = 14,byrow = T) expectancy <- matrix(c(1,2,3,4,5,6,7,8,10,00),byrow = T) data <- list('shock','nsubjects','ntrials','expectancy') myinits <- list(list(tau = runif (42,3),alpha = runif (42,1))) parameters <- c("tau",'alpha','v','predk') # jags sampling samples <- jags(data,inits=myinits,parameters,model.file = RW_model,n.chains=1,n.iter=1000,n.burnin=500,n.thin=1,DIC=T)
解决方法
由于设置中没有随机节点，我不确定 MCMC 模拟会给您带来什么，但这里有一些有效的代码。据我所知，主要问题是当表显示 c(Vmax-Vn) 时，它没有使用 c() 作为连接运算符。 c 是定义为 0.3 的常数。
dat <- list(
v = c(0,rep(NA,8)),vn = rep(NA,8),vmax=1,c=0.3,Ntrials = 9
)
jmod <- "model{
for(i in 2:Ntrials){
v[i] <- c*(vmax - v[(i-1)]) + v[(i-1)]
}
}"
library(runjags)
out <- run.jags(jmod,data=dat,monitor="v")
编辑：回答已编辑的问题。
您必须查看这是否正在执行您想要的操作，但它会产生结果。据我所知，这是将 Stan 代码直接扩展到多个主题和每个试验的多个选择移植到 JAGS。一、模型：
RW_model <- function(){
for(i in 1:Nsubjects){
for(j in 1:Ntrials){
expectancy[i,j] ~ dcat(p[i,1:Nposs,j])
mu[i,j] <- tau[i] * v[i,j]
for(k in 1:Nposs){
q[i,k,j] <- exp(mu[i,j])
p[i,j] <- q[i,j]/sum(q[i,j])
}
pe[i,j] <- shock[i,j] - v[i,expectancy[i,j],j]
for(k in 1:Nposs){
v[i,j+1] <- v[i,j] + ifelse(expectancy[i,j] == k,alpha[i] * pe[i,0)
}
}
}
for (i in 1:Nsubjects){
tau [i] ~ dunif (0,3)
alpha [i] ~ dunif (0,1)
}
}
接下来，我们可以制作一些数据。此处的数据来自 20 名受试者，expectancy 中有 5 种可能的选择和 14 次试验。
set.seed(40120)
# v has to be a Nsubjects x Nchoices x Ntrials+1 matrix
v <- array(NA,dim=c(20,5,15))
# The first trial of v is initialized to 0 fo all subjects and choices
v[,1] <- matrix(0,nrow=20,ncol=5)
# expectancy and shock are both Nsubjects x Ntrials matrices
expectancy <- matrix(sample(1:5,20*14,replace=TRUE),ncol=14)
shock <- matrix(sample(c(0,1),ncol=14)
dlist <- list(
Nsubjects = 20,Ntrials = 14,Nposs = 5,expectancy = expectancy,shock = shock,v=v
)
最后，我们可以指定初始值并运行模型。
myinits <- list(list(tau = runif (nrow(expectancy),3),alpha = runif (nrow(expectancy)0,1)))
parameters <- c("tau",'alpha')
library(R2jags)
# jags sampling
samples <- jags(dlist,inits=myinits,parameters,model.file = RW_model,n.chains=1,n.iter=1000,n.burnin=500,n.thin=1,DIC=T)
您可以从输出中看到它生成了样本并按主题估计了不同的 tau 和 alpha 参数。
samples$BUGSoutput
Inference for Bugs model at "/var/folders/55/q9y1hbcx13b5g50kks_p0mb00000gn/T//Rtmp32zeG4/model1726bdee0c99.txt",fit using jags,1 chains,each with 1000 iterations (first 500 discarded)
n.sims = 500 iterations saved
mean sd 2.5% 25% 50% 75% 97.5%
alpha[1] 0.6 0.3 0.1 0.4 0.7 0.9 1.0
alpha[2] 0.4 0.3 0.0 0.1 0.3 0.7 0.9
alpha[3] 0.3 0.2 0.0 0.1 0.2 0.4 0.9
alpha[4] 0.4 0.3 0.0 0.1 0.3 0.6 0.9
alpha[5] 0.4 0.3 0.0 0.1 0.3 0.6 1.0
alpha[6] 0.3 0.3 0.0 0.1 0.2 0.5 0.9
alpha[7] 0.4 0.3 0.0 0.2 0.4 0.6 0.9
alpha[8] 0.3 0.2 0.0 0.1 0.2 0.5 0.9
alpha[9] 0.3 0.3 0.0 0.1 0.3 0.5 0.9
alpha[10] 0.4 0.2 0.0 0.2 0.4 0.5 0.9
alpha[11] 0.3 0.3 0.0 0.1 0.2 0.5 0.9
alpha[12] 0.3 0.3 0.0 0.1 0.2 0.5 0.9
alpha[13] 0.3 0.3 0.0 0.1 0.3 0.5 0.9
alpha[14] 0.4 0.3 0.0 0.2 0.4 0.6 0.9
alpha[15] 0.5 0.3 0.0 0.2 0.5 0.7 1.0
alpha[16] 0.4 0.3 0.0 0.1 0.3 0.7 0.9
alpha[17] 0.4 0.3 0.0 0.2 0.4 0.6 0.9
alpha[18] 0.4 0.3 0.0 0.2 0.4 0.7 1.0
alpha[19] 0.4 0.3 0.0 0.2 0.4 0.7 1.0
alpha[20] 0.3 0.3 0.0 0.0 0.2 0.5 0.9
deviance 909.3 5.1 901.1 905.8 908.8 912.5 919.9
tau[1] 1.1 0.6 0.2 0.7 1.1 1.5 2.3
tau[2] 0.8 0.7 0.0 0.3 0.6 1.2 2.7
tau[3] 1.0 0.8 0.0 0.3 0.8 1.6 2.8
tau[4] 1.1 0.7 0.0 0.4 0.9 1.6 2.7
tau[5] 0.9 0.7 0.0 0.3 0.7 1.3 2.7
tau[6] 1.0 0.8 0.0 0.3 0.9 1.6 2.8
tau[7] 1.0 0.7 0.1 0.5 0.9 1.5 2.7
tau[8] 1.1 0.8 0.0 0.4 0.9 1.7 2.9
tau[9] 1.0 0.8 0.0 0.3 0.8 1.6 2.7
tau[10] 1.6 0.8 0.1 0.9 1.7 2.3 2.9
tau[11] 1.1 0.9 0.0 0.3 0.9 1.8 2.8
tau[12] 1.0 0.8 0.0 0.4 0.8 1.6 2.9
tau[13] 0.9 0.8 0.0 0.3 0.6 1.4 2.7
tau[14] 1.1 0.8 0.0 0.5 0.9 1.6 2.7
tau[15] 1.2 0.7 0.1 0.7 1.1 1.6 2.7
tau[16] 1.0 0.8 0.1 0.4 0.9 1.5 2.8
tau[17] 1.1 0.8 0.1 0.4 0.8 1.6 2.9
tau[18] 1.2 0.8 0.1 0.5 1.1 1.7 2.7
tau[19] 1.2 0.8 0.1 0.5 1.0 1.8 2.8
tau[20] 1.0 0.9 0.0 0.3 0.8 1.5 2.9
DIC info (using the rule,pD = var(deviance)/2)
pD = 12.9 and DIC = 922.1
DIC is an estimate of expected predictive error (lower deviance is better).
因为我不知道对这个模型有什么期望，所以您必须验证它是否真的按照您的预期运行，但似乎这是一个很好的起点。