The following file contains an example of running the program tt2.gibbs.
Other programs are run in a similar fashion.

To obtain the output for the AJE paper, I ran the program listed below
in R.  The program was run 4 times using 4 different starting
values, to check for convergence, although only one run is shown here.
Sensitivity runs for the prior distribution may also be important,
depending on your application.  Answers will not be identical from run
to run, due to the Monte Carlo variance.

R program (for 2 by 2 table, other situations are similar):

First run the command

> gibbs.sampler.out<-tt2.gibbs(38,87,2,35,5,5,2,30,0.5,0.5,0.5,0.5,
    0.5,1,1,21.96,5.49,4.1,1.76,4.44,13.31,71.25,3.75,20500)

to perform the Gibbs sampler.

As described in the tt2.doc file, the first 4 parameters are the data,
38,87,2,35 from the 2 by 2 table (see AJE paper for source of this
data).  The next 4 parameters, 5,5,2,30, are arbitrary starting values
for the number of true positive subjects in each cell.  Other starting
values can also be used, provided they are in the correct range. For
example, all are numbers between 0 and the cell size, eg, 38 for the
first entry, etc.  The next 5 parameters, 0.5,0.5,0.5,0.5,0.5 are
starting values for the two sets of sensitivities and specificities,
and the prevalence.  Again, these are arbitrary values, but must be
within the range 0 to 1.  Although they are theoretically arbitrary,
these starting values may affect rates of convergence.  Therefore in
practice, you may wish to use your "best guess" values for each
parameter. The next 10 parameters,
1,1,21.96,5.49,4.1,1.76,4.44,13.31,71.25,3.75, are the 10 parameters of
the 5 Beta prior distributions for the five parameters of interest, the
prevalence of disease, and the sensitivities and specificities of each
test.  For example, the beta prior parameters for the prevalence here
was a Beta(1,1) distribution, ie, a uniform prior.
Similarly, the next two parameters indicate a Beta(21.96,5.49) prior
for the sensitivity of test 1.  Please refer to my AJE paper for the
source of these prior distributions, and for references on how to
derive prior densities for your own particular situation.  The built-in
subroutines beta.to.mu and mu.to.beta may be helpful in this regard.
For example:

> beta.to.mu(21.96,5.49)
$mean:
[1] 0.8

$sd:
[1] 0.07499268

>  mu.to.beta(0.8,0.07499268)
$alpha:
[1] 21.96

$beta:
[1] 5.49

This allows conversion between beta parameters and their means and
standard deviations, and vice versa, as indicated by the example.

The last parameter, 20500, indicates the number of iterations of the
Gibbs sampler to run.  Depending on your computer setup, this number
may or may not be feasible in one run.  I would suggest starting with
small values, say 100, and gradually increasing the number until you
find the limits of your computer setup.

The line above creates output which can be summarized by next running the
command:

> tt2.sum(gibbs.sampler.out,500,1)

The first parameter is the file named in the first command above, which
stores the vectors.  The second command is the number of iterations
before convergence.  Inference is then based on all of the output after
the first 500 iterations, since the skip parameter is 1.  The output
you should obtain is:

$size:
[1] 20500

$qprev:
      2.5%      5.0%     25.0%     50.0%     75.0%    95.0%     97.5%
 0.5260121 0.5750092 0.7097868 0.7751134 0.830716 0.9004895 0.9230528

$qsens1:
      2.5%      5.0%     25.0%     50.0%     75.0%     95.0%     97.5%
 0.7900652 0.8067705 0.8568488 0.8874314 0.9144861 0.9458535 0.9544222

$qspec1:
      2.5%      5.0%     25.0%     50.0%     75.0%     95.0%     97.5%
 0.3731874 0.4141996 0.5772343 0.7041374 0.8210112 0.9343055 0.9570939

$qppv1:
      2.5%      5.0%    25.0%     50.0%     75.0%    95.0%    97.5%
 0.6285226 0.6900756 0.8487309 0.9189233 0.9603374 0.9887083 0.9931422

$qnpv1:
     2.5%      5.0%     25.0%     50.0%     75.0%     95.0%     97.5%
 0.273309 0.3407579 0.5270505 0.6340492 0.7269634 0.8297648 0.8571152

$qsens2:
      2.5%      5.0%     25.0%     50.0%     75.0%     95.0%     97.5%
 0.2226895 0.2347676 0.2741186 0.3048363 0.3390376 0.399007 0.4248368

$qspec2:
      2.5%      5.0%     25.0%     50.0%     75.0%     95.0%    97.5%
 0.9071047 0.9173849 0.9453828 0.9605303 0.9726511 0.9854504 0.9886054

$qppv2:
      2.5%      5.0%    25.0%     50.0%     75.0%     95.0%     97.5%
 0.8735122 0.8978902 0.9447286 0.9647505 0.9785255 0.9902885 0.9929387

$qnpv2:
      2.5%      5.0%     25.0%     50.0%     75.0%     95.0%     97.5%
 0.09755736 0.1269252 0.2156531 0.2844137 0.3663822 0.529226 0.5874898

These are all interpreted as quantiles for the various parameters.
For example, the output line

$qprev:
      2.5%      5.0%     25.0%     50.0%     75.0%    95.0%     97.5%
 0.5260121 0.5750092 0.7097868 0.7751134 0.830716 0.9004895 0.9230528

says that the 2.5% quantile of the prevalence is estimated to be
0.5260121, and the uper 97.5% quantile is estimated to be 0.9230528.
Therefore a 95% credible interval is ( 0.5260121, 0.9230528 ). The
median (50 percentile) is 0.7751134, and the inter-quartile range is
( 0.7097868, 0.830716 ).


The output is given in the order:

$size:    the number of gibbs iterations used for inference
$qprev:   prevalence of the condition
$qsens1:  sensitivity of test 1
$qspec1:  specificity of test 1
$qppv1:   positive predictive value of test 1
$qnpv1:   negative predictive value for test 1
$qsens2:  sensitivity of test 2
$qspec2:  specificity of test 2
$qppv2:   positive predictive value for test 2
$qnpv2:   negative predictive value for test 2