# Sampling variables and calculating likelihood in WinBUGS/OpenBUGS

I am trying to read some WinBUGS/OpenBUGS examples to figure out how to specify models. I can't seem to understand where the probabilistic dnorm, dunif etc. functions are sampling from the distribution and where they are calculating likelihoods. Also, does the order of the statements matter, because it seems that the order is reversed. The statements seem to use variables not defined yet. For example, on page 50 of the book "Introduction to WinBUGS for Ecologists", is a small example:

model(
# priors
pop.mean ~ dunif(0,5000)
prec <- 1/pop.variance
pop.variance <- pop.sd * pop.sd
pop.sd ~ dunif(0,100)

# likelihood
for (i 1:nobs) {
m[i] ~ dnorm(pop.mean, prec)
}
)


It seems that pop.mean ~ dunif(0,5000) samples pop.mean from a uniform distribution while m[i] ~ dnorm(pop.mean, prec) does not sample m[i] but calculates a part of the likelihood that is used to reject/accept the move. But there is no difference in the syntax! Also prec is calculated using pop.variance before pop.variance has been defined. What if we had written it the other way round? Would it be the same model?

-

1) The probabilistic dnorm, dunif etc. functions are just describing the probability distribution which the variable on the left hand side (lhs) is assumed to have. If the variable is a parameter, then it's a prior distribution. If the variable is data, then it's p(data | parameters).

2) The distribution is not what is used for sampling. Don't think that pop.mean ~ dunif(0,5000) means that pop.mean is actually drawn from a Uniform(0,5000) distribution when the code is executed! The sampling is done using various algorithms inside WinBUGS / JAGS. The code describes a statistical model, not an algorithm to be executed.

3) As Tomas points out, the order of statements doesn't matter.

prec <- 1/pop.variance
pop.variance <- pop.sd * pop.sd
pop.sd ~ dunif(0,100)


merely tells the underlying sampler that 1) the prior on pop.sd is a U(0,100) distribution, 2) pop.variance = pop.sd*pop.sd, and 3) prec = 1/pop.variance. You can rearrange those statements however you like without changing their meaning.

Pretty much everything that happens goes on under the hood of WinBUGS / JAGS. The software takes the statistical model, as described by the code, and processes it using whatever algorithms it determines are "best". I'll repeat: the code describes the model, it does not specify an algorithm for processing the model + data. This is the hardest thing to get about WinBUGS / JAGS, or functional programming languages in general (at least, it was for me,) and once you get it, programming becomes much easier.

-
+1 But I think it's more accurate to call BUGS a declarative language, not a functional one. – Jack Tanner Jan 5 '12 at 15:39
Good point; upon reflection, I think BUGS more resembles (in some loose sense) Prolog than Lisp, but it doesn't seem right to call it a logic programming language either. Maybe it should, though. – jbowman Jan 5 '12 at 18:15

There are two categories of programming languages (which WinBUGS definitely is):

• procedural (or imperative) languages - they specify how the thing should be done - they specify the algorithm, the sequence of steps

• declarative languages - they specify what should be done - they just declare the principle

The WinBUGS model code is a declarative language and your "problem" is that you are looking at it as if it was a procedural language. The genial thing on WinBUGS is that you just write formulas as if you were writing it in a paper :-). You just specify the model, you specify what, and you don't solve how the thing is computed (this is some MCMC magic you don't need to understand). So don't bother about calculating likelihoods - you just specify the facts about variables! This will make your head a lot lighter :-)

• the order of commands does not matter (but the indices do - it is different to write m from m[i])
• there is no difference in syntax because the principle is in fact the same! The variable m[i] or pop.mean is having the particular distribution of given parameters and it doesn't matter if they are derived from another distribution (like in hierarchical models) or they are fixed constants. This is the beauty of the language - don't be scare of it, enjoy it :-)