normGibbs {Bolstad2} | R Documentation |
normGibbs draws a Gibbs sample from the posterior distribution of the parameters given the data fron normal distribution with unknown mean and variance. The prior for μ given var is prior mean m0 and prior variance var/n0 . That means n0 is the 'equivalent sample size.' The prior distribution of the variance is s0 times an inverse chi-squared with kappa0 degrees of freedom. The joint prior is the product g(var)g(mu|var).
normGibbs(y, steps = 1000, type = 'ind', ...)
y |
A vector containing the data |
steps |
The number of iterations of Gibbs sampling to carry out |
type |
Either 'ind' for sampling from an independent conjugate prior or 'joint' for sampling from a joint conjugate prior. 'i' and 'j' can be used as compact notation |
... |
If type = 'ind' then the user can specify the prior for μ with a parameter priorMu which can either be a single number m0, or m0 and n0. if m0 and n0 are not specified then m0 and n0 are 0 by default. The user can also specify priorVar, which if given, must be a vector with two elements s0 and kappa0. If s0 and kappa0 are not given then they are zero by default. If type = 'joint' then priorMu must be a vector of length two with elements m0 and sd0. The user can also specify priorVar, which if given, must be a vector with two elements s0 and kappa0. If s0 and kappa0 are not given then they are zero by default. |
A data frame containing three variables
[1,] | mu | a sample from the posterior distribution of the mean |
[2,] | sig | a sample from the posterior distribution of the standard deviation |
[3,] | mu | a sample from the posterior distribution of the variance = sig^2 |
James M. Curran
## firstly generate some random data mu <- rnorm(1) sigma <- rgamma(1,5,1) y <- rnorm(100, mu, sigma) ## A \eqn{N(10,3^2)} prior for \eqn{\mu} and a 25 times inverse chi-squared ## with one degree of freedom prior for \eqn{\sigma^2} MCMCSampleInd <- normGibbs(y, steps = 5000, priorMu = c(10,3), priorVar = c(25,1)) ## We can also use a joint conjugate prior for \eqn{\mu} and \eqn{\sigma^2}. ## This will be a \emph{normal}\eqn{(m,\sigma^2/n_0)} prior for \eqn{\mu} given ## the variance \eqn{\sigma^2}, and an \eqn{s0} times an \emph{inverse ## chi-squared} prior for \eqn{\sigma^2}. MCMCSampleJoint <- normGibbs(y, steps = 5000, type = 'joint', priorMu = c(10,3), priorVar = c(25,1)) ## Now plot the results oldPar <- par(mfrow=c(2,2)) plot(density(MCMCSampleInd$mu),xlab=expression(mu), main = "Independent") abline(v=mu) plot(density(MCMCSampleInd$sig),xlab=expression(sig), main = "Independent") abline(v=sigma) plot(density(MCMCSampleJoint$mu),xlab=expression(mu), main = "Joint") abline(v=mu) plot(density(MCMCSampleJoint$sig),xlab=expression(sig), main = "Joint") abline(v=sigma)