0
$\begingroup$

Bayesian Inference for the Normal Distribution, I use the following r code to obtain the posterior distribution. Let's say the data, $X\sim N(\mu, \sigma^{2})$ and $\mu \sim N(0,10)$ and $\sigma \sim \text{discrete uniform} (0,10)$. I used rethinking package in r to build the above model. Also, I simulated the $10000$ samples from the posterior distribution. It can be done as follows.

install.packages(c("mvtnorm","loo","coda"), repos="https://cloud.r-project.org/",dependencies=TRUE)
options(repos=c(getOption('repos'), rethinking='http://xcelab.net/R'))
install.packages('rethinking',type='source')
library(rethinking)

set.seed(6)
x <- data.frame(x = rt(100,3))
xbar <- mean(x$x)
fit <- rethinking::map(
  alist(
    x ~ dnorm(mu, sigma),
    mu ~ dnorm(1, 10),
    sigma ~ dunif(0, 50)
  ),
  data=x)

precis(fit, corr=TRUE)

sim_post <- extract.samples(fit) 
dim(sim_post)
post_par <- apply(sim_post, 2, mean)
post_par

quantile(sim_post$mu ,  c(.05, .95))
quantile(sim_post$sigma, c(0.05, 0.95))

level <- 0.95

CI_low <- qnorm((1-level)/2, sim_post[,1], sim_post[,2])
CI_up <- qnorm((1+level)/2, sim_post[,1], sim_post[,2])

cover_prob <- mean((CI_low < 0) & (0 < CI_up))
cover_prob

I graphed both posterior parameters below. enter image description here

How can I obtain the coverage probability for mu? Why am I getting $1$ always? Are there any problem in r code? Are there any ways to obtain the coverage probability for this analysis?

$\endgroup$
  • $\begingroup$ It's not clear what kind of object sim_post might be and therefore makes it difficult to determine how you are trying to compute the interval. Nevertheless, I believe the answers to your questions will immediately become apparent if you were to visualize the posterior distribution: why not plot its histogram, for instance? $\endgroup$ – whuber Sep 26 '18 at 20:49
  • $\begingroup$ @whuber, the dim(sim_post) is [1] 10000 2, it contains simulated values of the posterior parameters. $\endgroup$ – score324 Sep 26 '18 at 20:51

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.