大都会黑斯廷斯 MCMC 与 R
Metropolis-Hastings MCMC with R
我正在尝试使用 MH 算法和 R 实现一个简单的 MCMC 问题是我得到了这个错误(我试图计算 alpha 但它不是 NA 问题)
Error in if (runif(1) <= alpha) { : missing value where TRUE/FALSE needed
这是我的函数,有人能发现问题吗?
PoissonMetropolisHastingRW = function(n=100000,lambda=10,p=0.5,x0=0){
x=rep(0,n); y=0; alpha = 0
x[1]=x0
for(i in 2:n){
if (x[i-1] == 0){
y = sample(c(0,1),1, prob=c(0.5,0.5))
alpha = min(1,((lambda^y)*x[i-1]*p)/((lambda^x[i-1])*y*(1-p)))
#alpha = min(1, ( ((lambda^y)*x[i-1])/( (lambda^x[i-1])*y) )*(p/(1-p)) ))
if(runif(1)<=alpha) {x[i]=y}
else {x[i]= x[i-1]}
}
if (x[i-1] > 0){
y = sample(c(x[i-1]-1,x[i-1]+1), 1, prob=c(1-p,p))
alpha = min(1,((lambda^y)*x[i-1]*p)/((lambda^x[i-1])*y*(1-p)))
#alpha = min(1, (((lambda^y)*x[i-1]/((lambda^x[i-1])*y))*(p/(1-p))))
if(runif(1) <= alpha) {x[i]=y}
else {x[i]= x[i-1]}
}
}
return(x)
}
如果 y
恰好为 0(并且每次迭代的概率为 0.5,这肯定会发生),则 alpha
为 0 / 0(因为 x[i-1] == 0
)。它给你 NaN
。条件 something <= NaN
提供 NA
.
我正在尝试使用 MH 算法和 R 实现一个简单的 MCMC 问题是我得到了这个错误(我试图计算 alpha 但它不是 NA 问题)
Error in if (runif(1) <= alpha) { : missing value where TRUE/FALSE needed
这是我的函数,有人能发现问题吗?
PoissonMetropolisHastingRW = function(n=100000,lambda=10,p=0.5,x0=0){
x=rep(0,n); y=0; alpha = 0
x[1]=x0
for(i in 2:n){
if (x[i-1] == 0){
y = sample(c(0,1),1, prob=c(0.5,0.5))
alpha = min(1,((lambda^y)*x[i-1]*p)/((lambda^x[i-1])*y*(1-p)))
#alpha = min(1, ( ((lambda^y)*x[i-1])/( (lambda^x[i-1])*y) )*(p/(1-p)) ))
if(runif(1)<=alpha) {x[i]=y}
else {x[i]= x[i-1]}
}
if (x[i-1] > 0){
y = sample(c(x[i-1]-1,x[i-1]+1), 1, prob=c(1-p,p))
alpha = min(1,((lambda^y)*x[i-1]*p)/((lambda^x[i-1])*y*(1-p)))
#alpha = min(1, (((lambda^y)*x[i-1]/((lambda^x[i-1])*y))*(p/(1-p))))
if(runif(1) <= alpha) {x[i]=y}
else {x[i]= x[i-1]}
}
}
return(x)
}
如果 y
恰好为 0(并且每次迭代的概率为 0.5,这肯定会发生),则 alpha
为 0 / 0(因为 x[i-1] == 0
)。它给你 NaN
。条件 something <= NaN
提供 NA
.