An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20111110/8ec24d3e/attachment.pl>
Gibbs sampler
2 messages · Gyanendra Pokharel, Duncan Murdoch
On 11-11-10 1:28 PM, Gyanendra Pokharel wrote:
I have the following code,
gibbs<-function(m,theta = 0.25, lambda =0.55, n =1){
alpha<- 1.5
beta<- 1.5
gamma<- 1.5
x<- array(0,c(m+1, 3))
x[1,1]<- theta
x[1,2]<- lambda
x[1,3]<- n
for(t in 2:(m+1)){
x[t,1]<- rbinom(1, x[t-1,3], x[t-1,1])
x[t,2]<-rbeta(1, x[t-1,1] + alpha, x[t-1,3] - x[t-1,1] + beta)
Are you certain that x[t-1,3] - x[t-1,1] + beta will always be positive? If it is negative, rbeta will return NaN. Duncan Murdoch
x[t,3]<- rpois(1,(1 - x[t-1,1])*gamma)
}
x
}
gibbs(100)
it returns only 1 or 2 values of theta, some times NaN, this may be if any
theta is greater than 1, which is used as the probability for the next
rbinom(), so returns NaN. Can some one suggest to solve this problem?
[[alternative HTML version deleted]]
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.