OpenBUGS:valor ausente na distribuição de Bernoulli
-
13-11-2019 - |
Pergunta
Estou tentando modelar o "tempo" de observação como variável aleatória com OpenBUGS via R (R2OpenBUGS).Se todos os horários de observação estiverem disponíveis (sem NA's) tudo funciona, mas se eu definir um dos horários para NA, nada acontece.Testei o mesmo código com WinBUGS e recebo o erro de trap 'NIL dereference (read)'.Então, minha pergunta é: há algo realmente errado no meu código ou meu modelo é muito estranho para BUGS?
Meu modelo é assim:
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)
}
E o código R fica assim:
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)
Solução
Se entendi bem sua pergunta, você está perguntando se esta expressão
obs[i] ~ dbern(p)
é estranho para Win/OpenBUGS, pois ele não irá lidar com o valor ausente.Não, eu não penso assim;bugs é capaz de lidar com valores ausentes dessa forma e até mesmo imputá-los - com distribuição posterior.
Mas tenho uma forte suspeita de que
y[(i-1)*2 + obs[i]+1] <- x[i]
é realmente estranho!Isso pode causar problemas de bugs à medida que você força o cálculo do índice usando a observação obs[i]
que é nulo.Isso é realmente estranho, você deveria tentar encontrar outra maneira de fazer isso.Primeiro tente simplificar o modelo para pular esta regra e aposto que o problema desaparece.