Clinician's
corner

Back to main page

Programmer's
corner
So, you use WinBUGS a lot? Want more?
Patrick Blisle
Division of Clinical Epidemiology
McGill University Health Center
Montreal, Quebec CANADA
patrick.belisle@rimuhc.ca

Last modification: 11 sep 2015















Version 2.1 (November 2015)
WriteWinBUGSScript
The writing of WinBUGS Scripts made easier
If you have worked with WinBUGS for some time, you certainly know how using WinBUGS script can make life easier, especially if you run WinBUGS repeatedly!
If you have to write lots of WinBUGS scripts, however, you may feel the need for an automated and faster way to write them.

[ WriteWinBUGSScript is a Perl program that facilitates writing WinBUGS scripts. Once you prepare your model, data and initial values files, the rest of the task of writing a WinBUGS script is fully automated. WriteWinBUGSScript will examine these files, identify which nodes are data and which are stochastic and write a WinBUGS script that will load the model and appropriate data/inits files, compile the model, run WinBUGS for some burn-in iterations, monitor stochastic nodes for more iterations and then save their histories/traces and summary statistics to both odc and text output files. ]


Aside from the above basic task, when running WinBUGS models, one often wants to try out different initial values, different prior distributions, or even different data sets. A convenient way to do so is to write WinBUGS scripts; however, each different set-up (in terms of data analyzed, prior distributions used, etc.) must be specified by the mean of a different script, which involves writing a lot of WinBUGS script files being very similar to each other but each with its own peculiarities.
That's where WriteWinBUGSScript becomes interesting, since it can automate this process as well.

The present document introduces the reader to the syntax of WriteWinBUGSScript input files and gives a few examples that will hopefully make its use straightforward.
Top
Syntax

The input files that WriteWinBUGSScript will read are plain text files consisting of a series of commands (one command per line).
The exhaustive list of WriteWinBUGSScript commands is given below (and – in an attempt to make their interpretation more intuitive – illustrated with either file names or node names); each command is explained in detail in the next section. Click any command line to jump to the related detailed subsection.



Top
Details

Note that every WriteWinBUGSScript command line (above) except the first starts with a flag (that is, a letter or keyword preceded by a dash) that indicates the nature of the information on the corresponding line. Only the first line (path to the model to be compiled) and -s lines are mandatory.

The different pieces of information that can be passed to WriteWinBUGSScript (WriteWinBUGSScript ) are:

Flag Signification
-d
-i
-dnm
-m
-dnp
-exp
-inv.logit
-x
-coda
-coda>
-o
-s
Data file (path)
Initial values
Do Not Monitor
Monitor
Do Not Plot
report descriptive statistics for Exp() of listed nodes
report descriptive statistics for Inverse logit() of listed nodes
variables to drop from data files prior to running WinBUGS script
Monitor a node for further CODA or BOA analysis
CODA Output file location
WinBUGS Output file location
Script file location (WriteWinBUGSScript output file location)


Top
Model

The first line of each block of code must consist in the path to the model to be compiled, followed by a series of 2, 3 or 4 numbers specifying the number of iterations to run.

File paths
I personally always design a project directory the same way when it involves WinBUGS programming, making it easier to find the different relevant files. For a project – MyProject, say – I save the different models used in the subdirectory MyProject/WinBUGS/models, the inits files in MyProject/WinBUGS/inits, parameters and data files in MyProject/WinBUGS/data, WinBUGS output files in MyProject/WinBUGS/log, the scripts — which I write almost exclusively through WriteWinBUGSScript — in MyProject/WinBUGS/scripts and the WriteWinBUGSScript input files in MyProject/WinBUGS/scripts-src; the file paths used in WriteWinBUGSScript input files can be absolute or relative, which explains the presence of the ../ leading the file names given in the examples throughout this document.

Note that slashes (/) and backslashes (\) can be used equally in file paths and that file paths including spaces must be surrounded by quotes.

Relative paths
When used in a WriteWinBUGSScript input file, relative paths are relative to the directory where your WriteWinBUGSScript input file is found.

Dot-dot (..) refers to the directory one level above the one where the WriteWinBUGSScript input file is; you can use as many dot-dot's as necessary to point out the file you need, e.g. ../../AnotherProject/WinBUGS/data/mydata.txt.

Top
Number of burn-in and monitored iterations

The first line of code may be completed by the number of burn-in and monitored iterations, separated by a slash (/). Note that the k-notation, that is, 2k = 2000, 10k = 10,000, etc, can be used. If thinning is necessary, its value will be given in a third component, e.g. 1k/10k/20, which would mean: burn-in for 1000 iterations, monitor nodes values for 10,000 iterations, keeping values for each 20th value only.

Memory clear

When running a program that monitors a large number of nodes and/or runs for a very large number of iterations, WinBUGS sometimes fails due to a lack of available memory, leading to a Trap error message. One work-around is to clear the memory after a pre-specified number of iterations. For example, we may decide to clear memory after each 1000 or 2000 iterations, etc.

To specify a number of iterations to run before a memory clear, a fourth number (see the above notation) can be added.
The first example below would clear the memory after each block of 2000 iterations, with a thinning factor of 2.
In the second example, the two consecutive slashes (and hence the missing third argument) would indicate that no thinning is to be done, but the memory is to be cleared after each block of 2500 iterations.

../models/model.txt 1k/20k/2/2k
../models/model.txt 1k/20k//2500


Top
Data files

-d ../data/data.txt ../data/hyperparameters1.txt

List data files to be loaded in WinBUGS.
I often use data files to specify hyperparameters' values (rather than explicitely using them in the WinBUGS program), as we often want to compare results obtained with different sets of hyperparameters' values. Doing so, we only need one copy of the WinBUGS program (and not one program per set of hyperparameters' values); WriteWinBUGSScript will write one script per set, which is all we need (and we'll be careful to save the results in different
output files).

WinBUGS input data files can be either rectangular data structures or lists. Multiple data files can be listed on the same -d command line or on multiple lines, as long as each line starts with the -d flag.

Top
Initial values

Initial values can be given in text files or as a series values for each node.

-i b.age=0 b.weight=0 mu=2.5

When given node values directly (rather than through an inits file), WriteWinBUGSScript saves initial values to a temporary file. See
example 3 for details and Removing temporary inits files section for inits files removal.

Sometimes, the number of components of a vector will depend on the data and may differ from run to run (e.g. when different data subsets are analyzed [female or male only], or when data are dropped if some explanatory variables are missing [age, weight, etc.]); that can be easily worked around through the use of the data-dependent vector length indicator, given between curly brackets.

The line

-i b=-2{J}

would initialize the vector b as b = c(-2, -2, -2, -2, ..., -2), of length J, where J would be defined in a data file referenced to by a -d flag. Not only does this save you from typing the value -2 a possibly large number of times, but it most interestingly makes your WriteWinBUGSScript input file readily reusable. Indeed, if the length if the vector b changes in an (often unexpected but often happening) data update, you need not worry about updating your inits files as well; simply resubmit your WriteWinBUGSScript input file and a new temporary inits file (including the right number of b components) will be rewritten as well as the WinBUGS script (it will link to the updated inits file).

(Remember that WriteWinBUGSScript reads your data file(s) before writing the script(s): thus the length of b in the initial values file will correspond to the value of J (in the data file) when WriteWinBUGSScript was run: if data changes, scripts written by WriteWinBUGSScript should be rewritten.)

Finally, it is possible to define initial values for a vector in which the first component is fixed: suppose that in your model a vector called b.race – of length 4 – models the differences between categories 2, 3 and 4 and the reference category race=1, as in
  b.race[1] <- 0 # reference category

  for (i in 2:4)
  {
    b.race[i] ~ dnorm(0, 0.0001)
  }
The following WriteWinBUGSScript command line

-i b.race=0{4*}

would initially set b.race=c(NA, 0, 0, 0).

It is also possible to let WriteWinBUGSScript randomly select the initial values from a Bernoulli distribution. In this case, the probability of sampling the value 1 can be specified directly through the Bernoulli probability p argument, or alternatively through its logit value, via the logit.p argument, as shown in code below.

If the node thus initialized is a chance node, it is necessary to specify its length between curly brackets after the dbern() statement. On the other hand, if the variable is data with missing values for which we want to impute 0/1 data with some specified probability p of sampling a 1, then specifying its length is not necessary. See last line of the examples below, with node a.variable.with.missing.values being imputed when NAs occur in the data.

-i healthy=dbern(p=0.6){4}
-i healthy=dbern(p=0.7){C}
-i healthy=dbern(logit.p=-0.35){C}
-i a.variable.with.missing.values=dbern(logit.p=-1.2)


Top
Stochastic nodes not to monitor

-dnm node1 node2
-dnm tau.*
-dnm mu{N\10}

This options allows you to give the list of stochastic nodes that are not of interest: neither summary statistics nor history/density plots will be saved for nodes listed in -dnm instruction line(s).

The *-notation (as in tau.*, above) can be used to indicate a series of nodes names beginning with the same letters.

Sometimes, it is of interest — for the sake of safety — to monitor a node for a limited number of its components only. Consider the example where a node mu is defined for a large number of individuals N (N=10,000, say); monitoring the nodes mu[1], mu[2], ..., mu[N] will most likely be of little interest — not accounting for the additional work load it would impose to WinBUGS, possibly making it substantially slower — but monitoring a limited (and randomly selected) number of them may be a good idea — at least in the first stages of the WinBUGS program development – to assess the good behavior of the program (in terms of mu estimation). The -dnm command variant above (where 10 can be replaced by any number less than the dimension of the node, which can be stated in terms of a constant read through data files – N, in this example – or by a specified number, e.g., -dnm mu{200\10}) would monitor 10 randomly selected items from the vector mu.

Top
Force node monitoring

The -m flag forces the monitoring of the following nodes declared, no matter what was said through the -dnm option.

-dnm *[i]
-m somevar

This option shows useful when all variables indexed by [i] (for example) are not to be monitored but for somevar[i].

Top
Stochastic nodes not to plot (history/density)

The -dnp flag allows to give a list of stochastic nodes for which neither estimated densities nor histories are to be plotted (in .odc output file); summary statistics for these nodes will still be reported, though (unless they already appear in a -dnm instruction line, in which case it is not necessary to list them again in a -dnp instruction line).

-dnp some.binary.node


Top
Exp()

It is sometimes interesting to define exp() of some nodes and monitor these exponentiated values, e.g. in computation of odds ratios. While models involving such nodes will often be based on log odds (log.OR, say), the values of interest are the actual values of the odds ratios. It seems (and often is) natural to simply define the node of interest (e.g., OR <- exp(log.OR)), but in practice it will sometimes cause WinBUGS to crash, especially if log.OR can take — even with low probability — high values: taking the exp() of high values in WinBUGS can lead to numerical problems and indeed cause it to crash.

-exp log.OR
-exp b.*

A possible workaround in such cases is to not compute and monitor the actual values of exp(log.OR) but simply take the exp() of the summary statistics for log.OR in the WinBUGS output, which ( — for median, 2.5% and 97.5% percentiles, not for mean and s.d.!! — will give some useful info on OR. The instructions line -exp log.OR would add a comment to the script written by WriteWinBUGSScript, as in

display(log)
check('c:/users/pbelisle/My Documents/Home/SomeProject/WinBUGS/models/model.txt')
compile(1)
gen.inits()
update(1000)
set(log.OR)
update(1000)
stats(*)
save('c:/users/pbelisle/My Documents/Home/SomeProject/WinBUGS/log/OddsRatio.txt')
save('c:/users/pbelisle/My Documents/Home/SomeProject/WinBUGS/log/OddsRatio.odc')
quit()
# when done, compute stats for exp() of log.OR


This comment will of course be ignored by WinBUGS if you use it 'as is', but the program
RunWinBUGSScript will take good note of its presence, run the script and, after completion, read the .txt WinBUGS output file, compute the exp() of the summary statistics for the nodes listed in the -exp instructions line(s) and add them to the Node statistics section.

Exp() of nodes multiplied by a numeric constant can be computed:

-exp b.age*5 b.avgrent/10

Is is also possible to label the nodes computed through the -exp command above through the "→" notation, as in

-exp b.age*10->decade.effect

which would add the label decade.effect to the node exp(b.age*10), as illustrated below.

Node statistics
node mean sd MC error 2.5% median 97.5% start sample
b.age 0.01753 0.3143 0.002933 -0.6043 0.0185 0.6362 501 10000
/odds ratio/
decade.effect -> -> -> 0.002374 1.2032 579.4040 501 10000 exp(10*b.age)


Top
Inverse logit()

In our experience, monitoring the values of a binomial proportion when it is extremely close to 0 or 1 will cause WinBUGS to crash when it tries to write the corresponding node statistics to output file. A work-around is not to monitor it but to define its inverse logit through a WinBUGS code line such as logit.p <- logit(p) and monitor the latter through the WriteWinBUGSScript command:

-inv.logit logit.p

The WinBUGS .txt output file will contain descriptive statistics for both logit.p and p (through inv.logit(logit.p)).

As illustrated above for -exp, nodes for which descriptive statistics are requested through the -inv.logit command can be labeled through the "→" notation, as in

-inv.logit logit.sensitivity->sensitivity


Top
Dropping variables from data files

When developping WinBUGS models, one sometimes starts with models less complex than the final model and increase its complexity step by step, after some validation of its good behavior. In that process, all data may not be used at every step; a hierarchical model, for example, may be built by adding independent explanatory variables one at a time. In theory, that implies preparing new data files at each step, as WinBUGS will not allow data variables not to be used in a model.

-x datavar1 datavar2

The -x option adds a comment to the WinBUGS script generated, which will read

check('c:/users/pbelisle/My Documents/Home/MyProject/WinBUGS/models/SomeModel.txt')
data('c:/users/pbelisle/My Documents/Home/MyProject/WinBUGS/data/list.txt')
data('c:/users/pbelisle/My Documents/Home/MyProject/WinBUGS/data/rectangle.txt')
# RunWinBUGSScript will drop the following variable(s) from file(s) above:
# datavar1 datavar2
As is the case for the comments generated by the
-exp command, this will have no effect if you run the generated script in WinBUGS directly: RunWinBUGSScript, however, will read the comment and drop the listed variables from data files at run time. Data files will remain unchanged on disk, as temporary copies of the files that are dropped one or more variables will be used.

Top
CODA

Some nodes may need further analysis through
CODA or BOA after your WinBUGS program has finished; that will be possible only if you saved the values for the node(s) in question for each iteration, which is done through this option.

-coda node1.to.monitor node2.to.monitor
-coda> ../log/coda/mynodevalues.txt


Note that the -coda flag must be used jointly with a -coda> flag, which gives the (relative) path of the file where the values of the coda-monitored nodes will be saved.

If only one node is monitored through the -coda command, node values will be saved in files with the name specified on the -coda> command line augmented with the coda Index and 1 suffixes — ../log/coda/mynodevaluesIndex.txt and ../log/coda/mynodevalues1.txt in the example above.

If more than one node is monitored, two CODA files will be written for each node, and the file names given above will be augmented by a litteral extension for each node (_a, _b, _c and so on).

It is possible to save values for every stochastic node to the CODA output file through the command

-coda *


Top
WinBUGS Output file location

The -o command line specifies the WinBUGS output file location.

-o output.txt

Whether you use a file name with .txt or .odc file extension – or even without a file extension – is irrelevant ; both .txt and .odc output files will be created by the WinBUGS script written by WriteWinBUGSScript; actually, the output files created will be – in this example – output-WinBUGSlog.txt and output.odc.

Top
Script file location

Specify where the script written by WriteWinBUGSScript will be saved through the -s command.

-s script.txt


Top
Exclude a few nodes from the node statistics section

Specify a few nodes to exclude from node statistics section.

-drns b.x b.y b.z

This option can be useful if the model includes, for example, the definition of a few Odds Ratios, e.g.,
  b.x ~ dnorm(0, 0.0001)
  b.y ~ dnorm(0, 0.0001)
  b.z ~ dnorm(0, 0.0001)
  
  exp.b.x <- exp(b.x)
  exp.b.y <- exp(b.y)
  exp.b.z <- exp(b.z)
where one would plot the density and the trace (or history) of the regression parameters b.x, b.y and b.z but would not really be interested in the corresponding node statistics but only in the node statistics of their exponential counterpart, which whould be done through the options below.

-dnp exp.b.x exp.b.y exp.b.y
-dnrs b.x b.y b.z


Top
Notes

  • Multiple scripts can be written through a single WriteWinBUGSScript input file; simply insert one or more blank lines between each script's description;
  • a script can be ignored by typing the line "# skip" on the line before its description;
  • comments are allowed (everything after a sharp sign (#) is ignored — except for # skip, above);
  • flags must appear in first column;
  • file names (either model, data or inits files) with spaces must be surrounded by double quotes (as illustrated in logreg example);
  • after completion, a list of scripts written can be found in WriteWinBUGSScript.log.


)) WriteWinBUGSScript was written to fill my own personal needs in terms of WinBUGS scripting: hence it may not fill all of your needs. It certainly does not exploit all the possibilities of WinBUGS scripting, and thus certainly should not be seen as a complete alternative to learning the basics of WinBUGS scripting. Nonetheless, the experienced (and frequent) WinBUGS user may still appreciate the time saved by using WriteWinBUGSScript.

If you come to use WriteWinBUGSScript on a regular basis but often need to modify/complete the WinBUGS scripts obtained, please do not hesitate to contact me if you feel that a larger audience may benefit from an upgraded version — along the lines that would meet your needs — of WriteWinBUGSScript. If I agree with you, I may very well consider the possibility of implementing new characteristics in versions to come.


Top
Examples

Example 1

Let's consider a logistic regression model, presented on Lawrence Joseph webpage. Model, data and initial values are reproduced below.

# Data (saved to file c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/data/fractures.txt)

list(sex=c(1, 1, 1, 0, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 1, 0, 0, 0,
1, 1, 1, 0, 0, 1, 1, 0, 1, 1, 1, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 1, 0, 0,
0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1,
0, 0, 1, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 0, 1, 0, 0, 1, 1, 1, 0, 1, 0, 1,
1, 1, 0, 1, 1, 1, 1, 1),
age= c(69, 57, 61, 60, 69, 74, 63, 68, 64, 53, 60, 58, 79, 56, 53, 74,
56, 76, 72, 56, 66, 52, 77, 70, 69, 76, 72, 53, 69, 59, 73, 77, 55, 77,
68, 62, 56, 68, 70, 60, 65, 55, 64, 75, 60, 67, 61, 69, 75, 68, 72, 71,
54, 52, 54, 50, 75, 59, 65, 60, 60, 57, 51, 51, 63, 57, 80, 52, 65, 72,
80, 73, 76, 79, 66, 51, 76, 75, 66, 75, 78, 70, 67, 51, 70, 71, 71, 74,
74, 60, 58, 55, 61, 65, 52, 68, 75, 52, 53, 70),
frac=c(1, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 0,
1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 0, 1, 0, 1,
0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 1, 0, 1, 0, 0, 1,
1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1,
0, 0, 1, 0, 0, 1),
n=100)
# Model (saved to file c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/models/age and sex.txt)

model
{
  age.bar <- mean(age[])

    for (i in 1:n) {
          logit(p[i]) <- alpha + b.sex*sex[i] + b.age*(age[i] - age.bar)
          frac[i]   ~ dbern(p[i])
    }
    alpha    ~ dnorm(0.0,1.0E-4)
    b.sex     ~ dnorm(0.0,1.0E-4)
    b.age     ~ dnorm(0.0,1.0E-4)
}
# Inits (saved to file c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/inits/logreg_parms.txt)

list(alpha=0, b.sex=1, b.age=1)
A WriteWinBUGSScript input file for this problem would be:

"../models/age and sex.txt" 250/5k
-d ../data/fractures.txt
-i ../inits/logreg_parms.txt
-o "../log/age and sex"
-s "../scripts/age and sex.txt"

The code above was saved to file c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/scripts-src/pines.txt;
as the model used is found in file c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/models/age and sex.txt, the first line of code could have consisted in the absolute path

"c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/models/age and sex.txt" 250/5k

but I personally prefer to use the shorter relative paths, as in

"../models/age and sex.txt" 250/5k


The WriteWinBUGSScript output file (a WinBUGS script) would look like:
display(log)
check('c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/models/age and sex.txt')
data('c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/data/fractures.txt')
compile(1)
inits(1, 'c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/inits/logreg_parms.txt')
gen.inits()
update(250)
set(alpha)
set(b.age)
set(b.sex)
set(p)
update(5000)
stats(*)
density(alpha)
density(b.age)
density(b.sex)
density(p)
history(alpha)
history(b.age)
history(b.sex)
history(p)
save('c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/log/age and sex.odc')
save('c:/users/pbelisle/My Documents/Home/logreg/WinBUGS/log/age and sex-WinBUGSlog.txt')
quit()
Note that the node age.bar, defined in the model, does not appear in the list of monitored nodes since WriteWinBUGSScript recognizes it as being deterministic, as it is the mean of other deterministic nodes.

Since list of nodes for which initial values were given is small (three nodes only), one could prefer to specify initial values directly in WriteWinBUGSScript input file, as shown below

"../models/age and sex.txt" 250/5k
-d ../data/fractures.txt
-i alpha=0 b.sex=1 b.age=1
-o "../log/age and sex"
-s "../scripts/age and sex.txt"

rather than through actually writing and saving an inits file and giving its path (through the same -i option).

Finally, the nodes p[1]--p[n] (if judged of no practical interest) could have been dropped from the monitoring list with the WriteWinBUGSScript input file below:

"../models/age and sex.txt" 250/5k
-d ../data/fractures.txt
-i alpha=0 b.sex=1 b.age=1
-dnm p
-o "../log/age and sex"
-s "../scripts/age and sex.txt"

This will potentially make the WinBUGS script faster to run and the WinBUGS odc output file smaller.

Example 2

Let's consider a more complex model – that is, a program with numerous nodes – where the idea of writing a WinBUGS script by hand is much less attractive. Consider the pines example, presented on The BUGS project webpage, for which models, data and initial values are reproduced below.

# MODEL (was saved in c:/users/pbelisle/My Documents/Home/pines/WinBUGS/models/pines.txt)
#
# Bayes factors using Carlin and Chib method.  
# For full description see Page 47 of Classic BUGS examples Vol 2.

model
{
# standardise data

   for(i in 1:N){
       Ys[i] <- (Y[i] - mean(Y[]))/sd(Y[])
       xs[i] <- (x[i] - mean(x[]))/sd(x[])
       zs[i] <- (z[i] - mean(z[]))/sd(z[])
    }
 
# model node

       j  ~ dcat(p[])
   p[1] <- 0.9995
  p[2] <- 1-p[1] # use for joint modelling
     pM2 <- step(j - 1.5)


# model structure

   for(i in 1:N){
       mu[1,i] <- alpha + beta *xs[i]
       mu[2,i] <- gamma + delta*zs[i]
       Ys[i]    ~ dnorm(mu[j,i],tau[j])
   }

# Model 1 

   alpha  ~ dnorm(mu.alpha[j],tau.alpha[j])
   beta   ~ dnorm(mu.beta[j],tau.beta[j])
   tau[1] ~ dgamma(r1[j],l1[j])


# Model 2

   gamma  ~ dnorm(mu.gamma[j],tau.gamma[j])
   delta  ~ dnorm(mu.delta[j],tau.delta[j])
   tau[2] ~ dgamma(r2[j],l2[j])

# Code below appears uncommented in the original WinBUGS code;
# I personally prefer to specify these values through data files 
# rather than through WinBUGS coding. 
# (Furthermore, doing so will allow WriteWinBUGSScript 
# to recognize the nodes below as deterministic!)

## estimation priors
#
#   mu.alpha[1]<- 0
#   tau.alpha[1] <- 1.0E-6
#   mu.beta[1] <- 0
#   tau.beta[1]  <- 1.0E-4
#   r1[1]      <- 0.0001
#   l1[1] <- 0.0001
#
## pseudo-priors 
#
#   mu.alpha[2]<- 0
#   tau.alpha[2] <- 256
#   mu.beta[2] <- 1
#   tau.beta[2]  <- 256
#   r1[2]      <- 30    
#   l1[2] <- 4.5
#
## pseudo-priors
#
#   mu.gamma[1] <- 0
#   tau.gamma[1] <- 400
#   mu.delta[1] <- 1
#   tau.delta[1] <- 400
#   r2[1]       <- 46      
#   l2[1] <- 4.5
#
## estimation priors
#
#   mu.gamma[2] <- 0
#   tau.gamma[2] <- 1.0E-6
#   mu.delta[2] <- 0
#   tau.delta[2] <- 1.0E-4 
#   r2[2]       <- 0.0001
#   l2[2] <- 0.0001
}
# DATA (were saved in c:/users/pbelisle/My Documents/Home/pines/WinBUGS/data/pines.txt)
#

list(N=42,
Y = c(3040, 2470, 3610, 3480, 3810, 2330, 1800, 3110, 3160, 2310,
           4360, 1880, 3670, 1740, 2250, 2650, 4970, 2620, 2900, 1670,
           2540, 3840, 3800, 4600, 1900, 2530, 2920, 4990, 1670, 3310,
           3450, 3600, 2850, 1590, 3770, 3850, 2480, 3570, 2620, 1890,
           3030,3030),
     x = c(29.2, 24.7, 32.3, 31.3, 31.5, 24.5, 19.9, 27.3, 27.1, 24.0,
           33.8, 21.5, 32.2, 22.5, 27.5, 25.6, 34.5, 26.2, 26.7, 21.1,
           24.1, 30.7, 32.7, 32.6, 22.1, 25.3, 30.8, 38.9, 22.1, 29.2,
           30.1, 31.4, 26.7, 22.1, 30.3, 32.0, 23.2, 30.3, 29.9, 20.8,
           33.2, 28.2),

     z = c(25.4, 22.2, 32.2, 31.0, 30.9, 23.9, 19.2, 27.2, 26.3, 23.9, 
           33.2, 21.0, 29.0, 22.0, 23.8, 25.3, 34.2, 25.7, 26.4, 20.0,
           23.9, 30.7, 32.6, 32.5, 20.8, 23.1, 29.8, 38.1, 21.3, 28.5,
           29.2, 31.4, 25.9, 21.4, 29.8, 30.6, 22.6, 30.3, 23.8, 18.4,
           29.4, 28.2))
# DATA [parameters list] 
# (were saved in c:/users/pbelisle/My Documents/Home/pines/WinBUGS/data/parameters.txt)
#
list(mu.alpha  = c(0, 0),
     mu.beta   = c(0, 1),
     mu.gamma  = c(0, 0),
     mu.delta  = c(1, 0),

     tau.alpha = c(1.0E-6, 256),
     tau.beta  = c(1.0E-4, 256),
     tau.gamma = c(400, 1.0E-6),
     tau.delta = c(400, 1.0E-4),

     l1 = c(0.0001, 4.5),
     l2 = c(4.5, 0.0001),
     r1 = c(0.0001, 30),
     r2 = c(46, 0.0001)
)
# INITS (were saved in file c:/users/pbelisle/My Documents/Home/pines/WinBUGS/inits/pines1.txt)
#

list(j = 2, tau = c(1,1), alpha = 0, beta = 0, gamma = 0, delta = 0)
A WriteWinBUGSScript input file could be as follows:

../models/pines.txt 1k/10k
-d ../data/pines.txt ../data/parameters.txt
-i ../inits/inits1.txt
-dnm j mu p xs Ys zs
-dnp pM2
-o ../log/pines
-s ../scripts/pines.txt

In this WriteWinBUGSScript input file, we listed a bunch of nodes not to be monitored (-dnm nodes lists), for diverse reasons:
  • deterministic nodes (xs, Ys and zs): if you omit them from this list, the WinBUGS script written will still be valid but less esthetic as it will cause WinBUGS warnings (monitor could not be set) when run;
  • nodes of little interest, e.g. mu (be it in terms of posterior densities or nodes statistics);
  • nodes that are only interesting after some computation (e.g. j is ignored in favor of pM2);
  • since pM2 is a binary variable, its density plot and trace are not interesting, hence it is listed in a 'do not plot' (-dnp) instruction line.

The result to the above WriteWinBUGSScript input file will look like below:
display(log)
check('c:/users/pbelisle/My Documents/Home/pines/WinBUGS/models/pines.txt')
data('c:/users/pbelisle/My Documents/Home/pines/WinBUGS/data/pines.txt')
data('c:/users/pbelisle/My Documents/Home/pines/WinBUGS/data/parameters.txt')
compile(1)
inits(1, 'c:/users/pbelisle/My Documents/Home/pines/WinBUGS/inits/inits1.txt')
gen.inits()
update(1000)
set(alpha)
set(beta)
set(delta)
set(gamma)
set(pM2)
set(tau)
update(10000)
stats(*)
density(alpha)
density(beta)
density(delta)
density(gamma)
density(tau)
history(alpha)
history(beta)
history(delta)
history(gamma)
history(tau)
save('c:/users/pbelisle/My Documents/Home/pines/WinBUGS/log/pines.odc')
save('c:/users/pbelisle/My Documents/Home/pines/WinBUGS/log/pines-WinBUGSlog.txt')
quit()

Example 3

Let's consider a last example, the pig weights, also presented on The BUGS project webpage, for which model and data are reproduced below.

# Model (was saved in c:\users\pbelisle\My Documents\Home\pigweights\WinBUGS\models\pigweights.txt)

# Histogram smoothing adapted from Example 5.9 from Congdon (2001), p 180.  
#  Illustrates a structured precision matrix for a multivariate normal prior
#
# Pig Weight Gain
	
model{ 
	y[1:s] ~ dmulti(th[1 : s] , n)
	sum.g <- sum(g[])

# smoothed frequencies

	for (i in 1 : s) {     
		Sm[i] <- n * th[i]
		g[i] <- exp(gam[i])    
		th[i]  <- g[i] / sum.g
	}

# prior on elements of AR Precision Matrix  

	rho ~ dunif(0, 1)
	tau ~ dunif(0.5, 10)

# MVN for logit parameters

	gam[1 : s] ~ dmnorm(mu[], T[ , ])
	for (j in 1:s) { 
		mu[j] <- -log(s)
	}

# Define Precision Matrix

	for (j in 2 : s - 1) {
		T[j, j] <- tau * (1 + pow(rho, 2))
	}

	T[1, 1] <- tau T[s, s] <- tau

	for (j in 1 : s-1) { 
		T[j, j + 1] <- -tau * rho
		T[j + 1, j] <- T[j, j + 1]
	}

	for (i in 1 : s - 1) {
		for (j in 2 + i: s) {
			T[i, j] <- 0
			T[j, i] <- 0 
		}
	}

}
# Data (were saved in c:\users\pbelisle\My Documents\Home\pigweights\WinBUGS\data\pigweights.txt)

list(y=c(1,1,0,7,5,10,30,30,41,48,66,72,56,46,45,22,24,12,5,0,1),n=522,s=21)
A WriteWinBUGSScript input file could be as follows:

../models/pigweights.txt 500/10k
-d ../data/pigweights.txt
-i gam=-3{s}
-o ../log/pigweights
-s ../scripts/pigweights.txt

where each component of the vector gam (of length s, defined through data file) is initialized with the value -3.
display(log)
check('c:/users/pbelisle/My Documents/Home/pigweights/WinBUGS/models/pigweights.txt')
data('c:/users/pbelisle/My Documents/Home/pigweights/WinBUGS/data/pigweights.txt')
compile(1)
inits(1, 'C:/DOCUME~1/patrick/LOCALS~1/Temp/WBugs/__wbinits20071008i.txt')
gen.inits()
update(500)
set(Sm)
set(T)
set(g)
set(gam)
set(mu)
set(rho)
set(sum.g)
set(tau)
set(th)
update(10000)
stats(*)
density(Sm)
density(T)
density(g)
density(gam)
density(mu)
density(rho)
density(sum.g)
density(tau)
density(th)
history(Sm)
history(T)
history(g)
history(gam)
history(mu)
history(rho)
history(sum.g)
history(tau)
history(th)
save('c:/users/pbelisle/My Documents/Home/pigweights/WinBUGS/log/pigweights.odc')
save('c:/users/pbelisle/My Documents/Home/pigweights/WinBUGS/log/pigweights-WinBUGSlog.txt')
quit()
Note the line giving the path to the file where initial values were temporarily written, C:/DOCUME~1/patrick/LOCALS~1/Temp/WBugs/__wbinits20071008i.txt in this case.

By default, temporary files (including initials values) are written to your system temporary directory, which WriteWinBUGSScript will find through system environment variables. However, if you do not like that location, you can choose your own by editing the program WriteWinBUGSScript.pl and change the line
# $TMPDIR = "c:/my computer/my tmp dir";
for, e.g.,
$TMPDIR = "c:/myname/tmp";
(Do not forget to remove the leading sharp sign (#), which is used to comment out lines in Perl programs.)

Note that temporary inits files are not removed automatically, in case you'd wish to resubmit a WriteWinBUGSScript-generated WinBUGS script involving one of them in the future: however, WriteWinBUGSScript can easily be used to remove the tmp inits file it has written in the past upon a simple request.

Top
Preparing multiple scripts

Very often, we want to run an analysis on a data set with slightly different models, different sets of parameters, different sets of initial values, etc. Preparing a WinBUGS script file for each desired scenario can be done by writing the relevant code to the WriteWinBUGSScript source file, separating the code for each script by a blank line.

The code below, for example, will create two scripts (saved to files ../scripts/age and sex-1.txt and ../scripts/age and sex-2.txt) when submitted to WriteWinBUGSScript, one for the initial values found in file ../inits/inits-1.txt and another for the initial values found in file ../inits/inits-1.txt, respectively.

"../models/age and sex.txt" 250/5k
-d ../data/fractures.txt
-i ../inits/inits-1.txt
-o "../log/age and sex-1"
-s "../scripts/age and sex-1.txt"

"../models/age and sex.txt" 250/5k
-d ../data/fractures.txt
-i ../inits/inits-2.txt
-o "../log/age and sex-2"
-s "../scripts/age and sex-2.txt"

The writing of such source files can be made easier through the use of loop, a tool to repeat code with slight modifications each time the code is repeated. The code below, for example, when submitted to loop, would write a WriteWinBUGSScript source file consisting of the code to prepare 24 WinBUGS scripts, one for each combination of the two models, the three data sets and the four sets of initial values compared, listed in the upper part of the file.

"age and sex scripts.txt"
{
    M     1 2
    PART     hip knee femur
    I     1 2 3 4
}

"../models/model-<M>.txt" 250/5k
-d ../data/<PART>.txt
-i ../inits/inits-<I>.txt
-o "../log/age and sex-<M>-<PART>-<I>"
-s "../scripts/age and sex-<M>-<PART>-<I>.txt"

See the page on loop Perl program for details.

Top
Download WriteWinBUGSScript

WriteWinBUGSScript is a free program: version 2.1 is directly available from this page. Save and unzip it under the name WriteWinBUGSScript.pl and read the section below.

Top
How to use WriteWinBUGSScript

WriteWinBUGSScript is a program written in Perl. Please refer to my generic page on running Perl programs for instructions.

You can use WinBUGS scripts written by WriteWinBUGSScript as you would for any manually typed scripts but some features — like computing exp() and computing inverse logits — will only be handled by RunWinBUGSScript.

Removing temporary inits files

Temporary inits files written by WriteWinBUGSScript can be easily removed: at the DOS-prompt, type the path to Perl followed by the path to WriteWinBUGSScript.pl with the -clean flag, as in, e.g.,

c:\perl\bin\perl.exe c:\users\pbelisle\My Documents\Home\bin\WriteWinBUGSScript.pl -clean

Log file

Success or error messages are written to the log file WriteWinBUGSScript.log, found in the directory where you saved WriteWinBUGSScript.pl.