Dear Bugs list,
I was wondering if there is an established method for integrating the
predictive distribution of a model over a subset of the predictive
distributions support? I came across this problem whilst writing the
following programe. It is comparing the bivariate cost-effect data between
old and young patients, what I want to estimate is the probability of their
difference laying each of the four quadrants of a 2-dimensional plane. It's
very simple to sample from the predicitive distribution, and then define a
dichotomous variable by the indicator for each quadrant. However I don't
really want to use the resulting means, as it seems to me that they are
just plug-in estimators. But, whenever I add the last two lines, so as to
define them as multinomial variables, I get the error message informing me
that p[] has multiple definitions. It's easy to see that this approach will
fail whenever I want to integrate the predictive distribution. Is there a
better way to do it?
Simon Bond.
model
{
for( i in 1 :M ) {
old[i , 1:2] ~ dmnorm(mu.old[], R.old[ , ])
young[i ,1:2]~ dmnorm(mu.young[], R.young[,])
}
mu.old[1:2] ~ dmnorm(mean[],prec[ , ])
R.old[1:2 , 1:2] ~ dwish(Omega[ , ], 2)
mu.young[1:2]~dmnorm(mean[],prec[,])
R.young[1:2,1:2]~dwish(Omega[,], 2)
old.predict[1:2]~dmnorm(mu.old[], R.old[,])
young.predict[1:2]~dmnorm(mu.young[], R.young[,])
effect<-old.predict[1]-young.predict[1]
cost<-old.predict[2]-young.predict[2]
p[1]<-step(effect)*step(cost)
p[2]<-step(-effect)*step(cost)
p[3]<-step(-effect)*step(-cost)
p[4]<-step(effect)*step(-cost)
p[1:4]~dmulti(q[],4)
q[1:4]~ddirch(alpha[])
}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|