JAGS:可变数量的簇
JAGS: variable number of clusters
我正在尝试 运行 一个贝叶斯聚类模型,其中聚类的数量是随机的,具有二项分布。
这是我的 Jags 模型:
model{
for(i in 1:n){
y[ i ,1:M] ~ dmnorm( mu[z[i] , 1:M] , I[1:M, 1:M])
z[i] ~ dcat(omega[1:M])
}
for(j in 1:M){
mu[j,1:M] ~ dmnorm( mu_input[j,1:M] , I[1:M, 1:M] )
}
M ~ dbin(p, Mmax)
omega ~ ddirich(rep(1,Mmax))
}
到运行它,我们需要定义参数和变量的初始值,这是在这个R脚本中完成的
Mmax=10
y = matrix(0,100,Mmax)
I = diag(Mmax)
y[1:50,] = mvrnorm(50, rep(0,Mmax), I)
y[51:100,] = mvrnorm(50, rep(5,Mmax), I)
plot(y[,1:2])
z = 1*((1:100)>50) + 1
n = dim(y)[1]
M=2
mu=matrix(rnorm(Mmax^2),nrow=Mmax)
mu_input=matrix(2.5,Mmax,Mmax) ### prior mean
p=0.5
omega=rep(1,Mmax)/Mmax
data = list(y = y, I = I, n = n, mu_input=mu_input, Mmax = Mmax, p = p)
inits = function() {list(mu=mu,
M=M,
omega = omega) }
require(rjags)
modelRegress=jags.model("cluster_variabile.txt",data=data,inits=inits,n.adapt=1000,n.chains=1)
但是,运行最后一个命令,得到
Error in jags.model("cluster_variabile.txt", data = data, inits = inits,
: RUNTIME ERROR: Compilation error on line 6.
Unknown variable M Either supply values
for this variable with the data or define it on the left hand side of a relation.
这对我来说毫无意义,因为错误出现在第 6 行,即使 M 已经出现在模型的第 4 行! 运行这个脚本的实际问题是什么?
所以,我认为主要问题是您无法更改正在更新的随机节点的维数。这似乎是可逆跳跃 MCMC 的问题,但我认为你不能在 JAGS 中做到这一点。
所以 JAGS 不像 R 或其他编程过程语言,因为它实际上并不 运行 逐行,它是一种声明性语言,这意味着命令的顺序实际上并不重要至少在关于错误如何弹出的条款。所以仅仅因为它没有在第 4 行抛出错误并不意味着那里也没有错误。我不是肯定的,但我相信错误的发生是因为 JAGS 尝试在输入值之前先构建数组,所以 M 在这个阶段实际上没有定义,但你最终无能为力。
除此之外,应该有一个相当简单的解决方法,只是效率较低。不是从 1:M
开始循环,而是从 1:MMax
开始循环,这样尺寸实际上不会改变,它始终是 MMax x MMax。然后第 7 行只是将这些位置的 1:M 分配给一个值。这样做的缺点是需要你在模型拟合后做一些处理。因此,在每次迭代中,您需要提取采样的 M 并将矩阵 mu 过滤为 M x M,但这应该不会太难。如果您需要更多帮助,请告诉我。
我正在尝试 运行 一个贝叶斯聚类模型,其中聚类的数量是随机的,具有二项分布。 这是我的 Jags 模型:
model{
for(i in 1:n){
y[ i ,1:M] ~ dmnorm( mu[z[i] , 1:M] , I[1:M, 1:M])
z[i] ~ dcat(omega[1:M])
}
for(j in 1:M){
mu[j,1:M] ~ dmnorm( mu_input[j,1:M] , I[1:M, 1:M] )
}
M ~ dbin(p, Mmax)
omega ~ ddirich(rep(1,Mmax))
}
到运行它,我们需要定义参数和变量的初始值,这是在这个R脚本中完成的
Mmax=10
y = matrix(0,100,Mmax)
I = diag(Mmax)
y[1:50,] = mvrnorm(50, rep(0,Mmax), I)
y[51:100,] = mvrnorm(50, rep(5,Mmax), I)
plot(y[,1:2])
z = 1*((1:100)>50) + 1
n = dim(y)[1]
M=2
mu=matrix(rnorm(Mmax^2),nrow=Mmax)
mu_input=matrix(2.5,Mmax,Mmax) ### prior mean
p=0.5
omega=rep(1,Mmax)/Mmax
data = list(y = y, I = I, n = n, mu_input=mu_input, Mmax = Mmax, p = p)
inits = function() {list(mu=mu,
M=M,
omega = omega) }
require(rjags)
modelRegress=jags.model("cluster_variabile.txt",data=data,inits=inits,n.adapt=1000,n.chains=1)
但是,运行最后一个命令,得到
Error in jags.model("cluster_variabile.txt", data = data, inits = inits,
: RUNTIME ERROR: Compilation error on line 6.
Unknown variable M Either supply values
for this variable with the data or define it on the left hand side of a relation.
这对我来说毫无意义,因为错误出现在第 6 行,即使 M 已经出现在模型的第 4 行! 运行这个脚本的实际问题是什么?
所以,我认为主要问题是您无法更改正在更新的随机节点的维数。这似乎是可逆跳跃 MCMC 的问题,但我认为你不能在 JAGS 中做到这一点。
所以 JAGS 不像 R 或其他编程过程语言,因为它实际上并不 运行 逐行,它是一种声明性语言,这意味着命令的顺序实际上并不重要至少在关于错误如何弹出的条款。所以仅仅因为它没有在第 4 行抛出错误并不意味着那里也没有错误。我不是肯定的,但我相信错误的发生是因为 JAGS 尝试在输入值之前先构建数组,所以 M 在这个阶段实际上没有定义,但你最终无能为力。
除此之外,应该有一个相当简单的解决方法,只是效率较低。不是从 1:M
开始循环,而是从 1:MMax
开始循环,这样尺寸实际上不会改变,它始终是 MMax x MMax。然后第 7 行只是将这些位置的 1:M 分配给一个值。这样做的缺点是需要你在模型拟合后做一些处理。因此,在每次迭代中,您需要提取采样的 M 并将矩阵 mu 过滤为 M x M,但这应该不会太难。如果您需要更多帮助,请告诉我。