# INLA functions (continued)

I have polished up one of the two functions I’ve thought of implementing for INLA and it’s now available in the development version of the R package. So if you’ve got INLA installed in your R version (see how you can do it here, if you don’t), you can update it to the development version by typing

`inla.upgrade(testing=TRUE)`

This will add the new function inla.contrib.sd which can be used to express the uncertainty in the structured (“random”) effects of your INLA model in terms of the standard deviation, instead of the precision (which INLA gives by default).

In the help, I consider a simple example: I simulate $N=12$ Binomial trials. For each, first I simulate the sample size to be a number between 80 and 100 (subjects)

```
n=12
Ntrials = sample(c(80:100), size=n, replace=TRUE)
```

Then I simulate the actual observations $y$ as Binomial variables depending on the sample size Ntrials and a probability prob which in turn depends on a variable $\eta\sim\mbox{Normal}(0,0.5)$, defined as $$\mbox{prob} = \frac{\exp(\eta)}{1+\exp(\eta)}$$

which I do with the code

```
eta = rnorm(n,0,0.5)
prob = exp(eta)/(1 + exp(eta))
y = rbinom(n, size=Ntrials, prob = prob)
```

At this point, I define a **very** simple hierarchical model where each trial represents a clustering unit (_ie _I assume a trial-specific effect). You do this in INLA by using the command

```
data=data.frame(y=y,z=1:n)
formula=y~f(z,model="iid")
m=inla(formula,data=data,family="binomial",Ntrials=Ntrials)
summary(m)
```

This produces an INLA object m and the standard summary is

```
Call:
"inla(formula = formula, family = \"binomial\", data = data, Ntrials = Ntrials)"
Time used:
Pre-processing Running inla Post-processing Total
0.20580840 0.02721024 0.78909874 1.02211738
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant kld
(Intercept) 0.1815634 0.1606374 -0.1371095 0.1815601 0.5003695 5.379434e-05
Random effects:
Name Model Max KLD
z IID model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant
Precision for z 4.506 2.300 1.508 4.032 10.290
Expected number of effective parameters(std dev): 10.16(0.5718)
Number of equivalent replicates : 1.181
Marginal Likelihood:,"56.18"
```

Now, while the summary of the posterior distribution for the *precision* of the structured effect $z$ is of course just as informative, it is in general (more) difficult to interpret, because it is on the scale of 1/standard deviation. On the other hand, you can’t just take take the reciprocal of the (square-rooted) summaries to obtain info about the posterior distribution of the standard deviation, because the transformation is not linear and thus it does not directly apply to the moments and quantiles.

One easy way out is to sample from the posterior marginal of the precision, transform the draws from this distribution to draws of the resulting distribution for $\sigma = \frac{1}{\sqrt{\mbox{precision}}}$ and then summarise *that* posterior marginal distributions. That’s what inla.contrib.sd does

`s <- inla.contrib.sd(m)`

[You need to specify the INLA model, m in this case, and you can also input the number of draws you require from the posteriors $-$ the default is nsamples=1000].

The resulting object s contains two elements:

`s$hyper`

returns a table

```
mean sd 2.5% 97.5%
sd for z 0.5096883 0.125651 0.3132221 0.7948447
```

with the summary on the standard deviation scale. The other element s$samples contains the simulated values from the posterior, which can be used for example to draw a histogram of the distribution

`hist(s$samples,xlab="standard deviation for z",main="")`

which in this case produces the following graph.