tag:blogger.com,1999:blog-3240271627873788873.post4915029573384886721..comments2024-03-26T06:46:11.752-04:00Comments on Doing Bayesian Data Analysis: Negative Binomial ReparameterizationJohn K. Kruschkehttp://www.blogger.com/profile/17323153789716653784noreply@blogger.comBlogger5125tag:blogger.com,1999:blog-3240271627873788873.post-764029313455600672015-09-27T04:10:54.792-04:002015-09-27T04:10:54.792-04:00This seems to be working to here. Anything obvioul...This seems to be working to here. Anything obvioulsy wrong?<br /><br />y[i] ~ dnegbin( p[x[i]] , r[x[i]] )<br />}<br />for ( grpIdx in 1:2) {<br />p[grpIdx] <- r[grpIdx]/(r[grpIdx]+m[grpIdx])<br />m[grpIdx] ~ dgamma(0.01, 0.01)<br />r[grpIdx] ~ dgamma(0.01, 0.01)<br />v[grpIdx] <- r[grpIdx]*(1-p[grpIdx])/(p[grpIdx]*p[grpIdx])<br />}<br />}<br />"<br /><br />snail<-read.table('alive14-15.txt',header=T)<br />attach(snail)<br /><br /># Load the data:<br /> y = c( snail$alive14 , snail$alive15) # combine data into one vector<br /> x = c( rep(1,length(alive14)) , rep(2,length(alive15)) ) # create group membership code<br /> Ntotal = length(y)<br /> # Specify the data in a list, for later shipment to JAGS:<br />dataList <- list(<br /> x = x,<br /> y = y ,<br /> Ntotal = Ntotal<br /> )<br />length(y) #check<br />length(x)<br /><br />paramaters<-c("m","r","p","v") #collect these<br />jagsModel <- jags.model(textConnection(modelString), data=dataList, n.chains=3, n.adapt=1e3)<br /><br />codaSamples <- coda.samples(jagsModel, variable.names=parameters, n.iter=3000, thin=1)<br />m <- as.matrix(codaSamples)<br /><br />thanks Graeme<br />Anonymousnoreply@blogger.comtag:blogger.com,1999:blog-3240271627873788873.post-148782807658259082015-08-09T09:46:15.578-04:002015-08-09T09:46:15.578-04:00Graeme:
I won't have a chance to write a ful...Graeme: <br /><br />I won't have a chance to write a full script, but the key thing to do is use nested indexing for the group identifier. Row i of the data file contains the count y[i] and the group identifier x[i]. Then the model statement uses<br /><br />for ( i in 1:Nrows ) {<br /> y[i] ~ dnegbin( p[x[i]] , r[x[i]] )<br />}<br />for ( grpIdx in 1:2 ) {<br /> p[grpIdx] <- r[grpIdx]/(r[grpIdx]+m[grpIdx])<br /> m[grpIdx] ~ dgamma(0.01,0.01)<br /> # etc. ...<br />}<br />John K. Kruschkehttps://www.blogger.com/profile/17323153789716653784noreply@blogger.comtag:blogger.com,1999:blog-3240271627873788873.post-48493820033191458192015-08-09T01:23:48.905-04:002015-08-09T01:23:48.905-04:00G'day, What would a full script look like when...G'day, What would a full script look like when comparing two vectors both of which have a negbinom distribution? I have been trying to transpose this post with BEST.R function, but not very well.<br /><br />Graemenoreply@blogger.comtag:blogger.com,1999:blog-3240271627873788873.post-83575456188435658212013-12-25T19:42:02.330-05:002013-12-25T19:42:02.330-05:00Thanks, Professor Kruschke. You inspired me to try...Thanks, Professor Kruschke. You inspired me to try this -- and JAGS as a vehicle for learning BUGS -- on another problem, available for public inspection at: http://hypergeometric.wordpress.com/2013/12/23/jags-for-finding-highs-and-lows-in-a-week-of-wikipedia-accesses/Anonymousnoreply@blogger.comtag:blogger.com,1999:blog-3240271627873788873.post-3193078596235222012012-11-06T12:11:24.478-05:002012-11-06T12:11:24.478-05:00This is great. I'm going to try it with my dat...This is great. I'm going to try it with my data now. Thanks for the tip since this seems like a common problem I will encounter with my count data.<br /><br />Cheers,<br />DanDaniel Hockinghttps://www.blogger.com/profile/07743548705954290964noreply@blogger.com