WinBUGS Code for Chapter 14, Binary Diagnostic Tests

## Exercise 14.5:  One test, one population with a gold-standard

Model 

model{
z1 ~ dbin(eta, nD)
z2 ~ dbin(theta, nDbar)
nD ~ dbin(pi, n)
PVP <- eta*pi / (eta*pi + (1-theta)*(1-pi))
PVN <- theta*(1-pi) / (theta*(1-pi) + (1-eta)*pi)
DLRp <- eta / (1-theta)
DLRn <- (1-eta) / theta
eta ~ dunif(0,1)
theta ~ dunif(0,1)
pi ~ dunif(0,1)
}
list(z1=815, z2 =327, nD=1023, nDbar=442, n=1465)   #Data
list(eta=0.5, theta=0.5, pi=0.5)  #Initial Values

## Exercise 14.6:  One test, one population with a gold-standard

Model 

model{ 
 y[1:2,1:2] ~ dmulti(p[1:2,1:2],n)
 p[1,1] <- pi*eta 
 p[1,2] <- (1-pi)*(1-theta)
 p[2,1] <- pi*(1-eta) 
 p[2,2] <- (1-pi)*theta
 eta ~ dbeta(a[1],b[1])
 theta ~ dbeta(a[2],b[2])
 pi ~ dbeta(a[3],b[3])
}
y[,1] y[,2]
815 115
208 327
END
list(n = 1465, a = c(1,1,1), b= c(1,1,1))
list(pi = 0.5, eta = 0.5, theta = 0.5)

## Example 14.2.2 and Exercise 14.8:  One test, one population with no gold-standard

Model 

model{
 y ~ dbin(PAP, n)
 PAP <- eta*pi + (1-theta)*(1-pi)
 eta ~ dbeta(19.34, 23.41)
 theta ~ dbeta(152.08, 2.53)
 pi ~ dbeta(1,1)
 #pi ~ dbeta(1,5)  # Alternate prior
 P <- 1-step(pi - 0.05)
}
list(y=0, n=100)   # Data
list(eta=0.5, theta=0.5, pi=0.5)  # Initial Values

Example 14.3.1 and Exercise 14.10: Two conditionally independent tests, two populations with no gold-standard.
Data are from two tests for N. salmonis in trout (Enoe et al. Pre Vet Med 2000).
Test 1 is ME and Test 2 is PCR.

Model 
model
{
y1[1:2, 1:2]  ~ dmulti(p1[1:2, 1:2], n1)
y2[1:2, 1:2]  ~ dmulti(p2[1:2, 1:2], n2)
p1[1,1] <- pi1*eta1*eta2 + (1-pi1)*(1-theta1)*(1-theta2)
p1[1,2] <- pi1*eta1*(1-eta2) + (1-pi1)*(1-theta1)*theta2
p1[2,1] <- pi1*(1-eta1)*eta2 + (1-pi1)*theta1*(1-theta2)
p1[2,2] <- pi1*(1-eta1)*(1-eta2) + (1-pi1)*theta1*theta2
p2[1,1] <- pi2*eta1*eta2 + (1-pi2)*(1-theta1)*(1-theta2)
p2[1,2] <- pi2*eta1*(1-eta2) + (1-pi2)*(1-theta1)*theta2
p2[2,1] <- pi2*(1-eta1)*eta2 + (1-pi2)*theta1*(1-theta2)
p2[2,2] <- pi2*(1-eta1)*(1-eta2) + (1-pi2)*theta1*theta2
eta1 ~ dbeta(2.82, 2.49)		
theta1 ~ dbeta(15.7, 1.30)	
eta2 ~ dbeta(8.29, 1.81)		
theta2 ~ dbeta(10.69, 2.71)  
pi1 ~ dbeta(1.27, 9.65)	     
pi2 ~ dbeta(1.73, 2.71)		 
}
list(n2=30, n1=132, y1=structure(.Data=c(0,0,3,129),.Dim=c(2,2)),
                                  y2=structure(.Data=c(3,0,24,3),.Dim=c(2,2)))
list(pi1=0.03, pi2=0.30, eta1=0.55, theta1=0.98, eta2=0.90, theta2=0.85)  

Model with U(0,1) Priors 
Uniform priors
model
{
y1[1:2, 1:2]  ~ dmulti(p1[1:2, 1:2], n1)
y2[1:2, 1:2]  ~ dmulti(p2[1:2, 1:2], n2)
p1[1,1] <- pi1*eta1*eta2 + (1-pi1)*(1-theta1)*(1-theta2)
p1[1,2] <- pi1*eta1*(1-eta2) + (1-pi1)*(1-theta1)*theta2
p1[2,1] <- pi1*(1-eta1)*eta2 + (1-pi1)*theta1*(1-theta2)
p1[2,2] <- pi1*(1-eta1)*(1-eta2) + (1-pi1)*theta1*theta2
p2[1,1] <- pi2*eta1*eta2 + (1-pi2)*(1-theta1)*(1-theta2)
p2[1,2] <- pi2*eta1*(1-eta2) + (1-pi2)*(1-theta1)*theta2
p2[2,1] <- pi2*(1-eta1)*eta2 + (1-pi2)*theta1*(1-theta2)
p2[2,2] <- pi2*(1-eta1)*(1-eta2) + (1-pi2)*theta1*theta2
eta1 ~ dbeta(1,1)		
theta1 ~ dbeta(1,1)	
eta2 ~ dbeta(1,1)	
theta2 ~ dbeta(1,1) 
pi1 ~ dbeta(1,1)	     
pi2 ~ dbeta(1,1)	 
}
list(n2=30, n1=132, y1=structure(.Data=c(0,0,3,129),.Dim=c(2,2)),
                                  y2=structure(.Data=c(3,0,24,3),.Dim=c(2,2)))
list(pi1=0.03, pi2=0.30, eta1=0.55, theta1=0.98, eta2=0.90, theta2=0.85)  

Exercise 14.14:  Estimating the prevalence distribution of OPP in the US

Model 
model{
 for(i in 1:k){
   y[i] ~ dbin(p[i],n[i])
   p[i] <- eta*pi[i] + (1-theta)*(1-pi[i])
   pi[i] ~ dbeta(alpha,beta) 
 }
 PPD ~ dbeta(alpha,beta) # Generates
 alpha <- mu*psi         # pred prev dens
 beta <- psi*(1-mu)
 mu ~ dbeta(2.64,7.55)
 psi ~ dgamma(4 , 0.2)
 eta ~ dbeta(71,3)
 theta ~ dbeta(32,4) 
}

Inits 
list(eta=0.98,theta=0.99,mu=0.20,psi=20, PPD=0.2,
pi=c(0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 
0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2))

Exercise 14.15: Estimating the prevalence distribution of hypopneumonia in swine
Model, data, and inits 
model{
 for(i in 1:k){
   y[i] ~ dbin(p[i],n[i])
   p[i] <- pi[i]*eta + (1-pi[i])*(1-theta)
   pi[i] ~ dbeta(alpha,beta)I(0.0001,0.9999)
 }
 alpha <- mu*psi 
 beta <- (1-mu)*psi
 PPD ~ dbeta(alpha,beta) # This will generate
                                             # the pred prev density
 mu ~ dbeta(1, 10) 
 psi ~ dgamma(1, 1)
 eta ~ dbeta(24,3) 
 theta ~ dbeta(38,2) 
}
list(pi=c(0.0001,0.7,0.7,0.7,0.7,0.7,0.7,0.7,0.7,
              0.7,0.7,0.7,0.7,0.7,0.7,0.7,0.7,0.7,0.7), 
     eta=0.9, theta=0.95, mu=0.7, psi=10, PPD=0.7)
list(k=19, n=c(2128,20,20,19,20,12,20,19,29,10,21,20,20,21,20,20,21,21,20),
      y=c(5,15,13,4,14,10,17,12,10,9,20,8,16,20,13,6,12,20,17))