1
$\begingroup$

I am writing a hierarchical BUGS model that involves both linear and angle variables. I want the hyper-parameters to be normally distributed, which is straight-forward for the linear variables, but I'm not sure what to do about the angle variables. I would gladly use a von Mises distribution or a truncated normal, but neither is available in the BUGS language. My angle variables range from $0$ to $\pi$.

Here is my model:

model {

   # parameters
   for ( i in 1:P ) {
      Y[i] ~ dnorm(Yhat[i], eachtau[plot[i]])
      Yhat[i] <- C[plot[i]] + a[plot[i]] * exp(-v[plot[i]] * year[i]) * cos(w[plot[i]] * year[i] - phi[plot[i]])
   }

   # hyperparameters
   for ( i in 1:M ) {
      C[i] ~ dnorm(mu.C, tau.C)
      a[i] ~ dnorm(mu.z, tau.z)
      v[i] ~ dnorm(mu.v, tau.v)

      # what do we do with these angles?!
      w[i] ~ ?
      phi[i] ~ ?

      eachtau[i] ~ dgamma(0.001, 0.001)
   }

   # priors
   mu.C ~ dnorm(0.0, 0.001)
   tau.C ~ dgamma(0.001, 0.001)
   mu.z ~ dnorm(0.0, 0.001)
   tau.z ~ dgamma(0.001, 0.001)
   mu.v ~ dnorm(0.0, 0.001)
   tau.v ~ dgamma(0.001, 0.001)
}

If it's important, I'm actually using JAGS (rather than OpenBUGS or WinBUGS), using the rjags package in R, but I believe the model syntax is the same for both.

What's the simplest way to model the angle variables with something like a normal in BUGS?

$\endgroup$
  • $\begingroup$ What you're actually doing is providing priors on them. You can always construct your own using the zeroes trick or the ones trick (users.aims.ac.za/~mackay/BUGS/Manuals/Tricks.html), or you could, for example, set them equal to $\pi$ times beta variates. The latter has the disadvantage of a discontinuity in the pdf at $\pi$, but maybe this doesn't matter. $\endgroup$ – jbowman Feb 3 '12 at 19:36
  • $\begingroup$ Do you mean that the way I've written the model makes the normal distributions priors? Because that's not what I'm intending to do. My data consists of M plots and at each plot I have a value that is measured over time (P=M*#time-points). I want to fit a damped oscillator (the Yhat line) to each plot; I suspect that the damped oscillator parameters (amplitude, damping coefficient, etc.) are similar from plot to plot. So I want to have hyper-parameters on the damped oscillator parameters that are normally distributed. I am really interested in the mean and variance of these hyper-parameters. $\endgroup$ – mkosmala Feb 3 '12 at 22:53
0
$\begingroup$

If you are willing/able to try a different MCMC package, you might consider Stan. There are rstan, pystan, and other interfaces, or you could use the command-line version CommandStan. It does have a von miss distribution.

Two caveats:

  1. I've indirectly used the von miles distribution, through an R package called brms, which uses rstan under the hood. I found it to be a bit tricky in my particular problem, but this may be just me.

  2. Stan, unlike BUGS/JAGS is more procedural. For example, a for loop in BUGS/JAGS is a concise way to describe plates in your model. It's not a literal loop that's executed at run time. In Stan, a for loop is literally executed at run time. There are other consequences of Stan's similar-looking-but-different language, such as the order of statements within a block mattering.

So it may be too large of a leap for you to switch to Stan for just a single distribution, but I wanted to suggest it. Stan is compiled to C++ code, which is compiled to an executable and it's a very efficient sampler. If you are using R and can set up your model as a regression, I can highly recommend the brms package or rstanarm, both of which use rstan (and hence Stan) under the hood. (You do need to have a C++ compiler on your system to use brms, which is the default under Linux, and straightforward on the Mac.)

| cite | improve this answer | |
$\endgroup$
0
$\begingroup$

Here's an example of estimating amplitude and phase parameters in a JAGS model: http://doingbayesiandataanalysis.blogspot.com/2012/07/its-getting-warmer-in-wisconsin.html That model is extended in a subsequent post: http://doingbayesiandataanalysis.blogspot.com/2012/10/bayesian-estimation-of-trend-with-auto.html Apologies if I misunderstood what you were asking; but I hope that helps.

| cite | improve this answer | |
$\endgroup$

Not the answer you're looking for? Browse other questions tagged or ask your own question.