Dear BugsListers
I want to thank respondents to a previous question on modeling MPN. I am compiling responses and
follow-up questions to submit to the listserve in case anyone is interested.
I have a new problem in how to structure the data (or model):
Dataset description:
Each record in the dataset is an individual water sample. For each water sample, a maximum of 5
dilutions are run with 3 test tubes per dilution. In some cases, less than 5 dilutions are run
when previous steps in the test method indicate that 3 or 4 dilutions are sufficient. Each test
tube results in a positive (=1) or negative (=0) reading depending on whether growth occurred. The
data for each dilution is reported as the number of positive tubes.
The problem is how to represent samples that have less than 5 dilutions. I would like to work with
the existing dataset structure and create an index for each sample. The problem is that I cannot
use 0 or NA as the 'n' in the binomial distribution. Any ideas?
Cheers, Chris
===========================================================================================
The set up of the model would be:
t[i,k] = # of positive tubes per dilution
v[k] = dilution level of test k
n[k]= # of tubes per dilution = 3
lambda[i] - concentration of organisms in water sample [i]
t[i,k] ~ binomial (p[i,k], n) - the distribution for the test results per dilution level
p[i,k] ~ 1-exp(-lambda[i]*v[k]) - obtained from the complement of the first term in a Poisson.
lglambda <-log(lambda[i])
lglambda ~dnorm(mu[i],tau[i])
mu[i] would then be fit to a GLMM with priors on parameters
--------------------------------------------------------------------------------------------
And in Bugs notation, I would like to set it up as:
# i - index on a water sample; I=4014
# t1...t5 - the number of positive tubes for dilutions 1 to 5, respectively
# n1...n5 - the total number of tubes dilutions 1 to 5, respectively. normally each are 3;
sometimes not available
# v1...v5 - the volume of the dilution
model
{
for (i in 1:I)
{
# number of positive tubes have binomial distributions.
#**** HERE IS WHERE I NEED TO FIGURE OUT HOW TO REPRESENT CASES WHERE n4, n5, or both DO NOT
EXIST.****
t1[i]~dbin(p1[i],n1[i])
t2[i]~dbin(p2[i],n2[i])
t3[i]~dbin(p3[i],n3[i])
t4[i]~dbin(p4[i],n4[i])
t5[i]~dbin(p5[i],n5[i])
# here a common lambda is found for each sample and set of dilutions
p1[i] <- 1-exp(-lambda[i]*v1)
p2[i] <- 1-exp(-lambda[i]*v2)
p3[i] <- 1-exp(-lambda[i]*v3)
p4[i] <- 1-exp(-lambda[i]*v4)
p5[i] <- 1-exp(-lambda[i]*v5)
# lambda with a lnorm distribution
lglambda <-log(lambda[i])
lglambda ~dnorm(mu[i],tau[i])
# mu[i] would then be fit to a GLMM
...
__________________________________
Yahoo! Mail - PC Magazine Editors' Choice 2005
http://mail.yahoo.com
-------------------------------------------------------------------
This list is for discussion of modelling issues and the BUGS software.
For help with crashes and error messages, first mail [log in to unmask]
To mail the BUGS list, mail to [log in to unmask]
Before mailing, please check the archive at www.jiscmail.ac.uk/lists/bugs.html
Please do not mail attachments to the list.
To leave the BUGS list, send LEAVE BUGS to [log in to unmask]
If this fails, mail [log in to unmask], NOT the whole list
|