# NOT RUN {
## ===================================================== ##
## 2-TEST EXAMPLE: Campylobacter ##
## ----------------------------------------------------- ##
## Two tests were performed on 656 chicken meat samples ##
## -> T1 = enrichment culture ##
## -> T2 = direct plating ##
## The following assumption were made: ##
## -> TP is larger than 45% and smaller than 80% ##
## -> SE1 must lie within 24% and 50% ##
## -> SP1 and SP2 both equal 100% ##
## -> beta(30, 12) describes P(T2+|D+,T1+) ##
## The following results were obtained: ##
## -> 113 samples T1+,T2+ ##
## -> 46 samples T1+,T2- ##
## -> 156 samples T1-,T2+ ##
## -> 341 samples T1-,T2- ##
## ===================================================== ##
## how is the 2-test model defined?
define_x(2)
define_prior(2)
## fit campylobacter 2-test model
campy <-
truePrevMulti(
x = c(113, 46, 156, 341),
n = 656,
prior = {
theta[1] ~ dunif(0.45, 0.80)
theta[2] ~ dunif(0.24, 0.50)
theta[3] <- 1
theta[4] ~ dbeta(30, 12)
theta[5] ~ dbeta(1, 1)
theta[6] <- 1
theta[7] <- 1
}
)
## fit same model using 'list-style'
campy <-
truePrevMulti(
x = c(113, 46, 156, 341),
n = 656,
prior =
list(
theta1 = list(dist = "uniform", min = 0.45, max = 0.80),
theta2 = list(dist = "uniform", min = 0.24, max = 0.50),
theta3 = 1,
theta4 = list(dist = "beta", alpha = 30, beta = 12),
theta5 = list(dist = "beta", alpha = 1, beta = 1),
theta6 = 1,
theta7 = 1
)
)
## show model results
campy
## explore model structure
str(campy) # overall structure
str(campy@par) # structure of slot 'par'
str(campy@mcmc) # structure of slot 'mcmc'
campy@model # fitted model
campy@diagnostics # DIC, BGR and Bayes-P values
## standard methods
print(campy)
summary(campy)
par(mfrow = c(2, 2))
plot(campy) # shows plots of TP by default
plot(campy, "SE1") # same plots for SE1
plot(campy, "SE2") # same plots for SE2
## coda plots of TP, SE1, SE2
par(mfrow = c(1, 3))
densplot(campy, col = "red")
traceplot(campy)
gelman.plot(campy)
autocorr.plot(campy)
## ===================================================== ##
## 3-TEST EXAMPLE: Giardia ##
## ----------------------------------------------------- ##
## Three tests were performed on stools from 272 dogs ##
## -> T1 = immunofluorescence assay ##
## -> T2 = direct microscopy ##
## -> T3 = SNAP immunochromatography ##
## The following assumption were made: ##
## -> TP is smaller than 20% ##
## -> SE1 must be higher than 80% ##
## -> SP1 must be higher than 90% ##
## The following results were obtained: ##
## -> 6 samples T1+,T2+,T3+ ##
## -> 4 samples T1+,T2+,T3- ##
## -> 12 samples T1+,T2-,T3+ ##
## -> 12 samples T1+,T2-,T3- ##
## -> 1 sample T1-,T2+,T3+ ##
## -> 14 samples T1-,T2+,T3- ##
## -> 3 samples T1-,T2-,T3+ ##
## -> 220 samples T1-,T2-,T3- ##
## ===================================================== ##
## how is the 3-test model defined?
define_x(3)
define_prior(3)
## fit giardia 3-test model
giardia <-
truePrevMulti(
x = c(6, 4, 12, 12, 1, 14, 3, 220),
n = 272,
prior = {
theta[1] ~ dunif(0.00, 0.20)
theta[2] ~ dunif(0.90, 1.00)
theta[3] ~ dunif(0.80, 1.00)
theta[4] ~ dbeta(1, 1)
theta[5] ~ dbeta(1, 1)
theta[6] ~ dbeta(1, 1)
theta[7] ~ dbeta(1, 1)
theta[8] ~ dbeta(1, 1)
theta[9] ~ dbeta(1, 1)
theta[10] ~ dbeta(1, 1)
theta[11] ~ dbeta(1, 1)
theta[12] ~ dbeta(1, 1)
theta[13] ~ dbeta(1, 1)
theta[14] ~ dbeta(1, 1)
theta[15] ~ dbeta(1, 1)
}
)
## show model results
giardia
## coda densplots
par(mfcol = c(2, 4))
densplot(giardia, col = "red")
# }
Run the code above in your browser using DataLab