In a Hierarchical Bayesian Model, how can we sample and see how a prior distribution looks like if it contains hyperparameters with hyperpriors?

I have a Bayesian Hierarchical Model that looks like:

$$\begin{equation} Y_i \sim N(\mu, \sigma^2) \\ \mu \sim N(\mu_0, \sigma_0^2) \\ \sigma^2 \sim Gamma(1,1) \\ \mu_0 \sim N(0,1) \\ \sigma_0^2 \sim Gamma(2,2)\\ \end{equation}$$

In this set-up, I have found the posterior of $$\mu$$ and am trying to compare it to how the prior of $$\mu$$ looks like on a density plot to view shrinkage.

I am not sure how to construct the prior density of $$\mu$$ via sampling since it doesn't have fixed hyperparameters.

I thought about the following code in R:

mu0 <- rnorm(100, 0, 1)
sigma02 <- rgamma(100,2,2)
mu <- rnorm(100,mu0,sigma02)


But am not sure if I am properly doing the average, as it seems very spread out. Instead, should I take the mean of mu0 and sigma02 and plug them in like so?

mu0 <- rnorm(100, 0, 1)
sigma02 <- rgamma(100,2,2)
mu <- rnorm(100,mean(mu0),mean(sigma02))


I lose variability here, but I am not sure if either works.

• No the first solution is correct, it shows the diffusive effect of using an hyperprior. – Xi'an Apr 5 at 14:49