# Assistance to convert OpenBugs Code into Stan Code

Hi all,

I was introduced to OpenBUGS by a colleague and has been using it since. Recently I came across STAN and thought maybe I could give it a try. However there seems to be some error. Is there anyone who is kind enough to vet through my code below and enlightened me. Many thanks.

``````model{
mu ~ dnorm(0,1.0E-6)
beta ~ dgamma(1.0E-5,1.0E-5)

for (i in 1:6){
delta[i] ~ dnorm(mu,beta)

prec[i] <- 1/(u[i] * u[i])
x[i] ~ dnorm(delta[i],prec[i])

pred[i] ~ dnorm(mu,prec[i])
DOE[i] <- x[i] - pred[i]
}
}
``````
``````list(
x = c(34.30,32.90,34.53,32.42,31.90,35.80),
u = c(1.03,0.69,0.83,0.29,0.40,0.38)
)
``````
``````list(beta=1)
``````

I am trying to convert the above into:

``````data {
int<lower=0> N;
real x[N];
real<lower=0> u[N];
}

parameters {
real mu;
real <lower=0> beta;
real <lower=0> delta[N];
real <lower=0> pred[N];
}

transformed parameters{
real uu[N];
uu[N] = u[N] * u[N];
}

model {
mu ~ normal(0,1E-6);  //prior for mu
beta ~ gamma(1E-5,1E-5);  //prior for beta
delta[N]~normal(mu,beta);
x[N] ~ normal(delta[N],uu[N]);
pred[N]~normal(mu,uu[N]);
}

generated quantities{
real DOE[N];
DOE[N] = x[N] -pred[N];
}

I saved the above as STAN extension and calling it in R with the following:
```R script
library(rstan)
library(shinystan)

options(mc.cores = 4)
data <-list(N=6, y = c(34.3,32.9,34.53,32.42,31.9,35.8), u=c(1.03,0.69,0.83,0.29,0.4,0.38))

model1 <- stan(file='Openbugs convert HB REM1.stan', data=data)
print(model1)``````

Hi Benny, welcome to Stan!

There are a few things here. Firstly, I’d highly recommend reading the User’s Guide chapter on transitioning from BUGS: https://mc-stan.org/docs/2_25/stan-users-guide/stan-for-bugs-appendix.html

### Vectors

When working with multiple values for a single parameter, it is more efficient to use vectors, rather than arrays of `real`, see this section of the manual for more information on the various data types: https://mc-stan.org/docs/2_25/reference-manual/data-types-chapter.html.

``````real<lower=0> u[N];
``````

Use:

``````vector<lower=0>[N] u;
``````

### Indexing

Something to be very careful with your models is how you use indexes. Take the following statement for example:

``````delta[N]~normal(mu,beta);
``````

This is not putting a prior on the whole vector of `delta` values, only on the sixth value. This is because `N` is passed as the single integer 6, which means that statement is equivalent to:

``````delta~normal(mu,beta);
``````

If you want the prior (or any function) to apply to the whole vector, then you can just specify it without any indexing:

``````delta ~ normal(mu,beta);
``````

### Normal Distribution

Just a quick one. Stan parameterises the normal distribution using the mean and standard deviation, not the variance, so you don’t need to square the SD parameter (`u`)

### Predicted Values

To generate the model-predicted values, you should use the `_rng` function of your desired distribution in the `generated quantities` block, rather than declaring them as a new parameter.

So rather than having in the `model` block:

``````pred[N]~normal(mu,uu[N]);
``````

You would instead have in the generated quantities:

``````pred = normal_rng(mu, u);
``````

### The Fixed Model

I’ll post a corrected version of your model below, but make sure you take some time to see the changes I’ve made and why, and I’d also recommend reading the Stan documentation here: https://mc-stan.org/users/documentation/ to get more familiar with the language.

``````data {
int<lower=0> N;
vector[N] x;
vector<lower=0>[N] u;
}

parameters {
real mu;
real<lower=0> beta;
vector<lower=0>[N] delta;
}

model {
mu ~ normal(0,1E-6);  //prior for mu
beta ~ gamma(1E-5,1E-5);  //prior for beta
delta ~ normal(mu,beta);
x ~ normal(delta,u);
}

generated quantities{
real pred[N];
vector[N] DOE;

pred = normal_rng(mu, u);

// The vectorised rng returns a real[],
// need to convert to vector for subtraction
DOE = x - to_vector(pred);
}
``````
5 Likes

That also affects the prior parametrization. It should be

``````  mu ~ normal(0,1000);  //prior for mu
``````
2 Likes

Oh good catch, I missed that! Except the original prior was `normal(0,1E-6)`, so the SD prior would be `normal(0, 0.001)`

No, OpenBUGS uses precision instead of variance. You can see this in the original model where it says

``````prec[i] <- 1/(u[i] * u[i])
x[i] ~ dnorm(delta[i],prec[i])
``````
1 Like

Oh right, I was assuming that Benny had already put it in the variance scale, like with the `uu` parameter, but I should have checked that

1 Like

Hi @andrjohns @nhuurre, Great! Many thanks for taking your time to understand my OpenBUGS code and converting to STAN (also at the same time guiding me along!). To be honest, I am a chemist by training, I picked up the stats while on the job and recently tried exploring Bayesian stats. Hope I can learn from you all on this forum.
Thanks again, greatly appreciate your assistance.

3 Likes