##### Section 2:  WinBUGS model for ovarian cancer data

 model{ # Model for the outcomes for both tests
 for(i in 1:51){
  ty0CA199[i] <- log(y0CA199[i])   # Transform controls
  ty0CA125[i] <- log(y0CA125[i])
  ty0CA199[i] ~ dnorm(mu0199,tau0199) # Log normal model
  ty0CA125[i] ~ dnorm(mu0125,tau0125) # for controls 
 }
 for(i in 1:90){
  ty1CA199[i] <- log(y1CA199[i])  # Transform cases
  ty1CA125[i] <- log(y1CA125[i])
  ty1CA199[i] ~ dnorm(mu1199,tau1199)  # Log normal model
  ty1CA125[i] ~ dnorm(mu1125,tau1125)  # for cases
 }       # Diffuse priors
 mu0199~dnorm(0,0.01)   
 tau0199~dgamma(0.01,0.01)
 mu0125~dnorm(0,0.01)   
 tau0125~dgamma(0.01,0.01)
 mu1199~dnorm(0,0.01)   
 tau1199~dgamma(0.01,0.01)
 mu1125~dnorm(0,0.01)   
 tau1125~dgamma(0.01,0.01)
 AUC199 <- phi( (mu1199 - mu0199)/sqrt(1/tau0199 + 1/tau1199))
 AUC125 <- phi( (mu1125 - mu0125)/sqrt(1/tau0125 + 1/tau1125))
 for(i in 1:101){
# Get ROC values over the grid x[] = phi^{-1}(t[])
  ROC199[i]<- phi(sqrt(tau1199)*(mu1199-mu0199)+sqrt(tau1199/tau0199)*x[i])
  ROC125[i]<- phi(sqrt(tau1125)*(mu1125-mu0125)+sqrt(tau1125/tau0125)*x[i])
 }# Predictive probabilities of infection over the grid y[]
 for(i in 1:321){
  f1199[i] <- sqrt(tau1199)*exp((-tau1199/2)*pow(y[i]-mu1199,2))
  f0199[i] <- sqrt(tau0199)*exp((-tau0199/2)*pow(y[i]-mu0199,2))
  f1125[i] <- sqrt(tau1125)*exp((-tau1125/2)*pow(y[i]-mu1125,2))
  f0125[i] <- sqrt(tau0125)*exp((-tau0125/2)*pow(y[i]-mu0125,2))
  probMale125[i] <- f1125[i]*pm/(f1125[i]*pm + f0125[i]*(1-pm))
  probFemale125[i] <- f1125[i]*pf/(f1125[i]*pf + f0125[i]*(1-pf))
  probMale199[i] <- f1199[i]*pm/(f1199[i]*pm + f0199[i]*(1-pm))
  probFemale199[i] <- f1199[i]*pf/(f1199[i]*pf + f0199[i]*(1-pf))
 }  
} 



##### Section 3:  WinBUGS model, data, and inits for the regression example with continuous y

model{ 
 for(i in 1:n){ 
   z[i] ~ dbern(pi[i])
   y[i] ~ dnorm(mmu[i], ttau[i])
   mmu[i] <- z[i]*mu[1] + (1-z[i])*mu[2]
   ttau[i] <- z[i]*tau[1] + (1-z[i])*tau[2]
   logit(pi[i]) <- beta[1] + beta[2]*x[i]
 }
 tildepi[1] ~ dbeta(a[1],b[1])
 tildepi[2] ~ dbeta(a[2],b[2])
 beta[1] <- logit(tildepi[1])
 beta[2] <- logit(tildepi[2]) - logit(tildepi[1])
 for(i in 1:2){ 
   mu[i] ~ dnorm(d[i],0.001)
   tau[i] ~ dgamma(0.001,0.001)
 }
 v[1] ~ dnorm(mu[1],tau[1])  #Pred Dens Popn 1
 v[2] ~ dnorm(mu[2],tau[2])  #Pred Dens Popn 0
 # Get predictive probabilities of $D$
 for(i in 1:n){
   f1[i] <- sqrt(tau[1])*exp((-tau[1]/2)*pow(y[i]-mu[1],2))
   f2[i] <- sqrt(tau[2])*exp((-tau[2]/2)*pow(y[i]-mu[2],2))
   prob[i] <- f1[i]*pi[i]/(f1[i]*pi[i] + f2[i]*(1-pi[i]))
 }  
}
list(mu=c(55,45), tau=c(1,1), tildepi=c(.95,.05))
list(z=c(0, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1,
1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0))

list(n=300,a=c(19,1),b=c(1,19),d=c(20,10),
x = c(0, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,
1, 1, 1, 1, 1, 1,
1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0),
y = c(46.9688885299112, 52.3361408968752, 55.0587095330899, 67.063801419522,
68.8339157013614, 35.1997003210818, 43.8434567496886, 62.5678556871719,
51.7401497932518, 40.0210644491513, 55.5803936018858, 77.8540400706143,
53.2577006714074, 60.978096870618, 44.8963836121789, 45.8611825455399,
31.8055612025579, 44.0913994633267, 57.722105348682, 74.6109500714428,
52.7449019892315, 39.1355235642223, 46.0213939033561, 60.391334720546,
55.9977718278595, 44.1538782443230, 39.8951487485588, 57.8750840531612,
52.741948596222, 62.5174297594526, 39.7587406576789, 42.1347008412321,
29.9739391852203, 57.737985500304, 60.9027770410002, 40.8896197976564,
39.3549237757351, 42.9913777609811, 47.4636851909366, 40.19975650674,
46.8049462243584, 24.9840101346073, 53.913147855837, 31.7877258787095,
46.8465351493863, 44.6579914658539, 65.5782255615336, 47.3727177437165,
61.4277824354533, 41.3232437613357, 35.55101258649, 47.2046690388387,
54.7556206451842, 45.41255030363, 42.8147384529643, 22.7165793798118,
46.2008902700476, 55.2319677884315, 53.0094398983507, 48.3747036387213,
66.3509396059549, 43.739428059669, 83.3584841203558, 63.8976098188874,
42.9241764802072, 47.580661251016, 66.5911443421293, 57.8933498142899,
58.7485101903814, 41.1459589104188, 45.0476882605214, 48.17446507216,
41.1737178765139, 42.6121526136231, 48.7625758570532, 34.7531217958322,
47.8557167896222, 58.0750301395481, 32.8325393816697, 47.1903811877905,
40.163918721602, 45.2280622030011, 56.0167740512132, 43.8191866357756,
54.5691711543537, 37.7168944517151, 34.9276809305544, 62.3866423609381,
50.1037075520153, 58.4355197239564, 40.2671573484903, 60.7805867474096,
42.8168780864038, 55.4999551097814, 59.9444187794951, 44.5238375557945,
64.1243050143463, 24.5124935705296, 55.0657551062564, 58.9778393827595,
29.1429226326105, 31.4571459000624, 36.7874450804490, 34.5896554597068,
40.750392806842, 30.0467801951236, 40.7489409160689, 36.3809046657188,
58.2365214172514, 49.4283022810525, 33.9009170390867, 48.8656471648306,
33.020910338002, 41.4326955179173, 51.7548175724649, 39.979249177136,
55.3342775942271, 47.0222040653348, 59.6656425570817, 46.8155625838758,
52.9454435060603, 41.9878708909172, 26.4391187698536, 46.7718284534991,
46.8378047376953, 41.1191896022904, 42.2769904851212, 32.5460505660805,
55.2361044591861, 27.2078694980604, 35.9630078835694, 33.2292633602079,
41.6574609333793, 27.3834733152, 36.9821495866163, 54.1524350368548,
36.2217513370123, 31.8023503054429, 48.5910824045674, 36.3423467217617,
23.7152203775348, 51.896720662421, 49.0922660578311, 61.939437189927,
48.6708629284579, 43.7515213645761, 31.6912919345410, 42.6274611575775,
47.4228070397373, 44.5473638282546, 39.5870019393295, 50.8869434593182,
50.6223583666625, 40.5110515231160, 32.7256715242291, 35.2715751255559,
34.3933905356635, 32.004668642066, 41.4421841859785, 36.4075760228735,
43.3017576155405, 43.2893845266867, 36.035282852267, 42.4105732554219,
43.0190461497302, 32.312692090048, 28.2851310897362, 39.3906583601767,
37.9559659430498, 36.9724028427918, 28.5775671471887, 26.4029930319749,
25.8990848376682, 30.753554487765, 38.3609784614964, 35.7225535772337,
60.3281246520733, 41.9787824083943, 52.4280506016782, 40.3210195605944,
43.4956744280586, 46.2467042503358, 41.5944337046651, 28.0434248253336,
43.257898810863, 50.4003324736337, 37.4314642486052, 39.5386248666131,
29.8912644825252, 22.769206083817, 44.0710870114558, 28.0969782004725,
48.1169564193978, 45.7481824181866, 27.9146110824603, 39.0397759395597,
36.4960026167076, 40.1374871963424, 34.5431268184916, 46.7976471195881,
33.6964177698529, 58.0351050129309, 43.6274635621055, 26.5279249236766,
35.1883181962994, 47.6689047578843, 52.4641974360265, 49.7837014026475,
32.9192405411713, 57.4101669218539, 37.1894216284287, 34.8265181312403,
45.4235939524782, 23.0698515464671, 28.6360956987204, 45.3501290192816,
42.5864893361019, 47.2288184942705, 26.7834412894001, 50.1376097691537,
41.9241131506460, 34.9265683332031, 34.0555733161908, 42.1871614544994,
58.7520051047603, 20.6597945839066, 37.0059402377429, 55.1002441280777,
43.6682558700408, 44.0347669147137, 50.9780002317417, 37.4367349743465,
28.3676745207623, 48.809809867479, 27.7786770201964, 48.9781051215074,
29.4726186109028, 30.9009637589159, 39.6642946948769, 53.2769007089319,
41.0543914126919, 47.4113125630637, 41.3876752720242, 46.5457869044804,
34.5549470388404, 49.5040805196856, 43.0760514235605, 36.1391362732467,
32.6672931431525, 52.4161653503593, 52.2073058573022, 55.0187630400466,
60.1764373110256, 47.7763123062371, 43.5228153909299, 36.5087925071478,
42.2058517214953, 36.5612145931967, 28.8107075612686, 38.4372181214345,
45.9028943324997, 50.3490752334236, 53.8057938466392, 34.661751608125,
23.2465775481229, 50.9503934503450, 49.1443814623675, 39.0080927086418,
44.7646345585957, 21.2309912343723, 45.1648900374699, 38.164479174638,
48.523777652212, 25.7757424033611, 36.4532815990225, 46.0198156332465,
50.95234528742, 37.071586659518, 24.5320355627188, 19.1501760022900,
53.4340318219755, 23.8494224018, 35.0166384167196, 44.5539754278086,
29.5616357081007, 42.7567609118962, 41.7929758492584, 56.5246962019851,
39.5328891883512, 24.4217078546884, 38.5956517992541, 34.5088896020053,
48.6931763462509, 40.7392340133278, 57.8453253696438, 39.5016466623566,
14.5639898571131, 44.3014497066879, 28.1193738972429, 48.283793128681))










##### Section 4: WinBUGS model, data, and inits for the smoking cessation example

model {for (i in 1:N) {
y[i] ~ dbern(q[i]) # specification of the LR model with error
q[i] <- pi[i]*Se+(1-pi[i])*(1-Sp)
logit(pi[i]) <- beta[1] + x2[i]*beta[2] +
x3[i]*beta[3] + x4[i]*beta[4] + x5[i]*beta[5] + x6[i]*beta[6]
 }

Se ~ dbeta(99, 1) # the priors are specified for Se and Sp
Sp ~ dbeta(14, 2)
p[1] ~ dbeta(8, 15) # spec of the prior on 6 probs of smoking
p[2] ~ dbeta(10, 15) # cessation for 6 covariate combinations
p[3] ~ dbeta(3, 13)
p[4] ~ dbeta(8, 10)
p[5] ~ dbeta(4, 15)
p[6] ~ dbeta(6, 15)

# relate the regression coefficients to the
# probabilities on which prior was specified
# this results in an induced informative prior
# on the regression coefficients
beta[1] <- xinv[1,1]*logit(p[1]) + xinv[1,2]*logit(p[2])
                 + xinv[1,3]*logit(p[3]) + xinv[1,4]*logit(p[4])
                 + xinv[1,5]*logit(p[5]) + xinv[1,6]*logit(p[6])
beta[2] <- xinv[2,1]*logit(p[1]) + xinv[2,2]*logit(p[2])
                 + xinv[2,3]*logit(p[3]) + xinv[2,4]*logit(p[4])
                 + xinv[2,5]*logit(p[5]) + xinv[2,6]*logit(p[6])
beta[3] <- xinv[3,1]*logit(p[1]) + xinv[3,2]*logit(p[2])
                 + xinv[3,3]*logit(p[3]) + xinv[3,4]*logit(p[4])
                 + xinv[3,5]*logit(p[5]) + xinv[3,6]*logit(p[6])
beta[4] <- xinv[4,1]*logit(p[1]) + xinv[4,2]*logit(p[2])
                 + xinv[4,3]*logit(p[3]) + xinv[4,4]*logit(p[4])
                 + xinv[4,5]*logit(p[5]) + xinv[4,6]*logit(p[6])
beta[5] <- xinv[5,1]*logit(p[1]) + xinv[5,2]*logit(p[2])
                 + xinv[5,3]*logit(p[3]) + xinv[5,4]*logit(p[4])
                 + xinv[5,5]*logit(p[5]) + xinv[5,6]*logit(p[6])
beta[6] <- xinv[6,1]*logit(p[1]) + xinv[6,2]*logit(p[2])
                 + xinv[6,3]*logit(p[3]) + xinv[6,4]*logit(p[4])
                 + xinv[6,5]*logit(p[5]) + xinv[6,6]*logit(p[6])
}

# inits
list(Se = 0.95, Sp = 0.93, p = c(0.33, 0.39, 0.14, 0.44, 0.18, 0.26))
# data
list(y=c(0,1,0,1,0,0,0,0,0,1,0,0,0,0,1,1,0,0,0,0,0,0,1,0,0,1,0,0,0,1,0,
0,0,0,1,0,0,1,0,0,0,1,0,1,1,1,0,1,0,1,0,0,0,1,0,0,0,0,1,0,1,0,1,1,0,0,
0,0,0,0,1,0,1,0,1,0,1,1,1,1,0,0,0,1,0,1,0,0,0,0,0,0,1,0,1,1,0,1,1,0,1,
0,0,0,0,1,0,1,1,0,0,0,0,0,1,0,0,0,0,0,1,0,0,0,0,1,1,0,0,1,0,1,0,0,0,1,
0,0,0,1,0,0,0,0,1,0,1,0,1,0,0,0,0,0,0,0,0,0,1,1,0,0,1,0,0,0,1,0,0,0,0,
0,0,1,0,0,1,0,0,0,1,0,0,0,1,0,0,0,0,0,0,0,1,0,0,1,0,1,1,1,0,0,0,0,0,1,
1,0,0,1,0,0,1,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,1,1,0,0,0,0,0,
1,0,0,0,1,0,1,0,0,0,0,1,1,1,0,0,0,1,1,0,0,0,1,0,1,1,0,0,0,0,1,0,1,1,0,
0,0,1,1,0,1,0,1,1,0,0,1,0,0,1,1,0,0,1,1,0,1,0,0,0,0,0,0,1,1,1,0,1,0,0,
0,0,1,0,1,0,0,0,0,1,0,0,0,1,0,1,1,0,1,0,0,0,0,1,1,0,0,0,0,0,1,0,0,1,0,
0,0,1,0,0,0,0,0,0,0,0,1,0,1,1),
x2 = c(0,1,1,1,1,1,1,0,1,1,1,1,0,0,1,1,1,1,1,1,1,1,1,1,0,0,1,1,1,1,0,0,
0,1,1,1,0,1,1,1,1,1,1,0,0,1,1,0,1,1,1,1,0,0,1,1,0,0,1,0,1,1,0,1,1,1,1,
1,1,1,1,0,0,1,1,1,1,1,0,1,1,1,1,1,1,1,0,1,1,1,1,1,1,1,0,1,0,1,1,1,1,1,
0,0,0,0,1,0,1,1,1,1,1,1,0,1,1,1,1,0,0,0,1,1,1,1,1,1,1,1,0,1,1,1,1,1,1,
0,0,1,1,1,0,1,1,1,1,0,1,1,1,1,0,0,1,1,0,1,1,1,1,0,0,1,1,0,0,1,1,0,0,1,
1,1,1,0,0,1,1,0,1,0,1,1,0,1,1,1,1,0,1,1,1,0,1,0,1,0,0,1,1,1,1,1,1,0,1,
1,1,1,0,1,0,0,0,1,1,0,1,1,0,1,1,1,0,1,0,1,1,0,1,0,1,0,0,0,1,1,1,1,1,1,
1,1,1,1,0,1,1,1,1,0,1,1,1,0,0,1,1,1,1,1,1,1,0,1,1,0,1,1,0,1,0,0,1,0,1,
0,0,0,0,0,1,0,1,0,1,1,0,1,1,0,1,1,1,0,0,0,1,1,1,1,1,0,1,1,1,0,0,0,0,0,
0,1,0,0,1,1,0,1,0,1,1,1,0,0,1,1,1,1,1,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,
1,0,0,1,1,0,1,1,1,0,0,0,0,1),
x3 = c(1,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,1,0,
1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,1,0,1,0,0,0,0,0,0,0,
0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,
1,0,1,0,0,1,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,
1,1,0,0,0,0,0,0,0,0,1,0,0,0,0,1,0,0,0,1,0,0,0,0,0,1,0,0,1,1,0,0,1,1,0,
0,0,0,1,1,0,0,1,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,1,0,1,1,0,0,0,0,0,0,0,0,
0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,
0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,1,0,0,0,0,
1,1,0,1,0,0,1,0,1,0,0,0,0,0,1,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,1,0,1,
1,0,1,1,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
0,1,0,0,0,0,0,0,0,0,0,1,0,0),
x4 = c(0,1,1,1,0,1,1,1,0,1,1,0,1,0,0,1,0,0,0,0,1,1,1,0,0,0,1,1,1,1,1,0,
0,0,1,0,0,0,1,1,0,1,0,0,0,1,1,0,0,0,1,0,1,0,0,1,1,0,0,1,0,1,0,0,1,1,1,
0,1,1,0,0,0,1,1,1,1,1,0,1,0,0,0,1,1,0,0,0,0,0,0,1,1,0,0,1,1,0,0,1,1,1,
0,0,0,1,0,1,1,1,0,1,0,0,1,1,0,0,0,0,0,1,1,1,0,0,0,1,1,0,0,1,1,0,0,0,1,
1,1,1,1,0,0,0,1,1,0,1,0,0,0,1,1,0,0,1,0,0,1,0,0,0,0,0,0,0,1,1,1,0,1,0,
0,1,1,0,1,0,0,0,1,0,0,1,0,0,0,1,0,0,1,0,1,0,1,1,0,0,1,1,1,1,1,0,0,0,0,
1,0,1,0,0,0,0,0,1,0,0,0,1,0,0,1,1,0,1,0,1,0,0,0,0,1,0,0,0,0,1,0,1,1,1,
1,0,0,1,0,1,1,0,1,0,1,1,1,1,0,1,0,0,1,0,0,1,1,1,1,1,0,1,1,0,1,1,0,0,1,
1,0,0,1,1,0,0,0,0,0,1,1,0,0,0,1,0,1,0,0,1,0,1,0,1,0,1,1,1,1,0,1,1,0,0,
1,1,0,0,1,1,0,1,0,1,1,1,1,0,1,1,0,1,1,0,0,1,1,0,1,1,1,0,1,1,1,0,1,1,0,
1,0,1,1,1,1,1,1,0,0,0,1,0,1),
x5 = c(0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,1,1,1,1,0,0,0,0,1,1,0,0,0,0,0,0,
1,0,0,1,0,1,0,0,1,0,1,0,0,0,0,0,0,1,0,0,0,0,1,0,0,1,1,0,0,0,0,0,0,0,0,
1,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,1,1,1,0,0,1,1,0,0,1,1,0,0,0,
1,0,1,0,1,0,0,0,1,0,1,0,0,0,0,0,0,0,0,0,0,0,1,0,1,0,0,1,1,0,0,1,1,1,0,
0,0,0,0,0,0,1,0,0,1,0,1,0,1,0,0,0,0,0,1,1,0,1,0,0,1,1,1,1,0,0,0,0,0,0,
0,0,0,1,0,1,0,1,0,0,0,0,0,0,1,0,1,0,0,1,0,1,0,0,0,1,0,0,0,0,0,1,0,0,0,
0,0,0,0,1,0,0,1,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,
0,1,1,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,1,0,0,0,0,0,0,1,0,0,1,0,0,0,0,0,
0,1,0,0,0,0,1,1,1,0,0,0,1,1,1,0,1,0,0,1,0,0,0,1,0,1,0,0,0,0,0,0,0,0,1,
0,0,1,0,0,0,0,0,0,0,0,0,0,1,0,0,1,0,0,0,1,0,0,1,0,0,0,1,0,0,0,0,0,0,1,
0,1,0,0,0,0,0,0,0,0,0,0,0,0),
x6 = c(0,1,1,1,0,1,1,0,1,1,1,0,1,1,0,1,0,1,1,1,1,0,1,0,1,1,1,1,1,1,1,0,
0,0,1,1,1,1,1,1,1,1,0,1,1,1,1,1,1,1,1,0,1,1,1,1,0,0,1,0,0,0,1,1,0,1,0,
1,0,1,1,0,1,1,0,1,1,1,1,0,1,1,1,1,1,1,0,1,1,0,0,1,1,1,1,1,1,1,1,1,1,1,
1,1,0,0,1,0,1,0,1,1,1,0,1,0,1,0,1,1,1,1,1,0,0,1,1,1,1,1,1,1,1,1,1,1,1,
0,1,1,1,0,0,1,1,0,1,1,1,0,0,1,1,1,1,1,1,0,1,0,0,0,0,1,0,0,1,1,0,0,0,0,
0,0,1,1,1,1,0,0,1,0,0,1,1,0,1,1,0,0,1,1,1,1,1,0,0,0,1,1,1,0,1,0,1,1,1,
1,0,1,1,1,1,0,1,1,0,1,1,1,0,0,0,1,0,0,1,1,1,1,1,1,1,1,1,1,1,0,1,1,0,1,
0,1,1,1,0,1,1,1,1,1,0,1,1,0,1,1,0,1,1,1,1,1,1,1,1,0,1,1,0,1,0,0,1,0,0,
1,1,1,0,1,1,1,1,1,0,1,0,1,1,1,0,1,0,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,0,0,
1,1,1,1,0,0,1,0,1,1,0,0,1,1,1,1,1,1,1,0,1,1,1,1,1,1,1,1,1,1,0,1,1,0,0,
0,1,1,1,1,0,0,1,1,1,1,1,1,1),
xinv=structure(.Data=c(-1.0, 1.94E-16, 3.33E-16, 2.57E-16, 1.0, 1.0,
-7.77E-16, 2.50E-16, 1.0, -4.86E-17, 3.33E-16, 1.0,-1.0, 1.0, 1.0,
-3.75E-16, -3.33E-16, 1.0,1.0, -5.83E-16, -5.55E-16, -9.02E-17, -1.0,
8.88E-16,4.16E-16, -6.11E-16, -3.89E-16, 1.0, -1.0, 5.41E-16,1.0,
1.67E-16, -1.0, 3.33E-16, 5.55E-17, 7.77E-16),.Dim=c(6,6)), N=361)








##### Section 5: WinBUGS model, data, and inits for classical swine fever data

model{
x[1:4] ~ dmulti(p[1:4], n)
p[1] <- pi*(Sefat1*Sefat2+covDp) + (1-pi)*((1-Spfat1)*(1-Spfat2)+covDn)
p[2] <- pi*(Sefat1*(1-Sefat2)-covDp) + (1-pi)*((1-Spfat1)*Spfat2-covDn)
p[3] <- pi*((1-Sefat1)*Sefat2-covDp) + (1-pi)*(Spfat1*(1-Spfat2)-covDn)
p[4] <- pi*((1-Sefat1)*(1-Sefat2)+covDp) + (1-pi)*(Spfat1*Spfat2+covDn)
ls <- (Sefat1-1)*(1-Sefat2)
us <- min(Sefat1,Sefat2) - Sefat1*Sefat2
lc <- (Spfat1-1)*(1-Spfat2)
uc <- min(Spfat1,Spfat2) - Spfat1*Spfat2
pi ~ dbeta(13.322, 6.281) ### Mode=0.70, 95% sure > 0.50
Sefat1 ~ dbeta(9.628,3.876) ### Mode=0.75, 95% sure > 0.50
Spfat1 ~ dbeta(15.034, 2.559) ### Mode=0.90, 95% sure > 0.70
Sefat2 ~ dbeta(9.628, 3.876) ### Mode=0.75, 95% sure > 0.50
Spfat2 ~ dbeta(15.034, 2.559) ### Mode=0.90, 95% sure > 0.70
covDn ~ dunif(lc, uc)
covDp ~ dunif(ls, us)
rhoD <- covDp / sqrt(Sefat1*(1-Sefat1)*Sefat2*(1-Sefat2))
rhoDc <- covDn / sqrt(Spfat1*(1-Spfat1)*Spfat2*(1-Spfat2))
}
list(n=214, x=c(121,6,16,71))
list(pi=0.7, Sefat1=0.75, Spfat1=0.90, Sefat2=0.75, Spfat2=0.90)











##### Section 5: WinBUGS model, data, and inits for brucellosis data

model{
 for(i in 1:K){
  x[i, 1:3] ~ dmulti(p[i, 1:3], n[i])
  p[i,1] <- pi[i]*(Sebapa*Seriv+covDp) + (1-pi[i])*((1-Spbapa)*(1-Spriv)+covDn)
  p[i,2] <- pi[i]*(Sebapa*(1-Seriv)-covDp) + (1-pi[i])*((1-Spbapa)*Spriv-covDn)
  p[i,3] <- 1-p[i,1]-p[i,2]
  pi[i] ~ dbeta(alpha, beta)
 }
 Sebapa ~ dbeta(88.3,1.9) ## Mode=0.99, 95% sure > 0.95
 Spbapa ~ dbeta(13.3,6.3) ## Mode=0.70, 95% sure > 0.50
 Seriv ~ dbeta(130.7,15.4) ## Mode=0.90, 95% sure > 0.85
 Spriv ~ dbeta(99.7,6.2) ## Mode=0.95, 95% sure > 0.90
 alpha <- mu*psi
 beta <- psi*(1-mu)
 mu ~ dbeta(16.9,48.6) ## Mode=0.25; 95% sure < 0.35
 psi ~ dgamma(7.23, 1.28)
 ls <- (Sebapa-1)*(1-Seriv)
 lc <- (Spbapa-1)*(1-Spriv)
 us <- min(Sebapa,Seriv)-Sebapa*Seriv
 uc <- min(Spbapa, Spriv) - Spbapa*Spriv
 covDn ~ dunif(lc, uc)
 covDp ~ dunif(ls, us)
 rhoDp <- covDp / sqrt(Sebapa*(1-Sebapa)*Seriv*(1-Seriv))
 rhoDn <- covDn / sqrt(Spbapa*(1-Spbapa)*Spriv*(1-Spriv))
 pi21 ~ dbeta(alpha,beta)
 a[1] <- 1-step(pi21-0.50)
 a[2] <- 1-step(pi21-0.10)
 a[3] <- step(pi21-0.90)
 a[4] <- step(pi21-0.75)
}
list( Sebapa=0.99, Spbapa=0.70, Seriv=0.90, Spriv=0.95, mu=0.25, psi=7)
list(K=20, n=c(67,12,22,71,18,80,147,45,37,118,12,62,11,20,56,60,26,54,17,17),
x=structure(.Data=c(4,2,61,
1,1,10,
2,2,18,
2,0,69,
2,0,16,
2,0,78,
80,5,62,
4,5,36,
10,0,27,
10,12,96,
1,0,11,
1,1,60,
3,0,8,
3,7,10,
38,0,18,
1,0,59,
3,0,23,
13,3,38,
7,0,10,
1,1,15),.Dim=c(20,3)))