4

R(R2OpenBUGS)を介してOpenBUGSを使用して、観測「時間」を確率変数としてモデル化しようとしています。すべての観測時間が利用可能である(NAがない)場合、すべてが機能しますが、いずれかの時間をNAに設定すると、何も起こりません。WinBUGSで同じコードをテストしたところ、トラップエラー「NIL逆参照(読み取り)」が発生しました。だから私の質問は、私のコードに本当に何か問題があるのか​​、それとも私のモデルがバグに対してあまりにも奇妙なのかということです。

私のモデルは次のようなものです。

model{
 for(i in 1:k){
  obs[i] ~ dbern(p) #is the observation done at time 1 or 2?
  y[(i-1)*2 + obs[i]+1] <- x[i]
 }    
 for(i in 1:n){    
   y[i] ~ dnorm(mu,tau) 
 }    
 mu ~ dnorm(0,0.0001)
 tau~ dgamma(0.001,0.001)  
 p ~ dunif(0,1) 
}

そして、Rコードは次のようになります。

library(R2OpenBUGS)
x<-obs<-rep(NA,5)
for(i in 1:k)
{
  obs[i]<-sample(c(0,1),1) #observation time of ith observation
  x[i]<-rnorm(1) #observed values
}

obs[2]<-NA #one of the sampling times is missing
INITS <- list(list(tau=1,mu=0,p=0.5))
DATA  <- list(x=x,n=n,k=k,obs=obs)

ob <- bugs(
  data=DATA,
  inits=INITS,
  parameters.to.save=c("tau","mu","p","y"),
  model.file="BUGSModel.R",
  n.chains=1,
  n.iter=50,
  n.burnin=10,
  n.thin=1,    
  DIC=FALSE)
4

1 に答える 1

4

私があなたの質問をよく理解しているなら、あなたはこの表現かどうか尋ねています

obs[i] ~ dbern(p)

Win / OpenBUGSの場合は奇妙なので、欠落している値を処理しません。いいえ、そうは思いません。バグはこの方法で欠落している値を処理することができ、事後分布でそれらを代入することさえできます。

しかし、私はそれを強く疑っています

y[(i-1)*2 + obs[i]+1] <- x[i]

本当に変です!obs[i]これにより、 nullの観測値を使用してインデックスを計算するように強制されるため、バグに問題が発生する可能性があります。これは本当に奇妙です、あなたはそれをする別の方法を見つけようとするべきです。最初にモデルを単純化してこのルールをスキップしてみてください。問題がなくなることは間違いありません。

于 2012-02-27T22:05:59.403 に答える