# Random draws from the posterior distribution with Normal-Inverse-Wishart (NIW) prior. Given iid $d$-dimensional niche indicators $X = (X_1,\ldots,X_N)$ with $X_i \sim N(\mu, \Sigma)$, this function generates random draws from $p(\mu,\Sigma | X)$ for the Normal-Inverse-Wishart (NIW) prior. ```r niw.post(nsamples, X, lambda, kappa, Psi, nu) ``` ## Arguments - `nsamples`: The number of posterior draws. - `X`: A data matrix with observations along the rows. - `lambda`: Location parameter. See 'Details'. - `kappa`: Scale parameter. Defaults to `kappa = 0`. See 'Details'. - `Psi`: Scale matrix. Defaults to `Psi = 0`. See Details. - `nu`: Degrees of freedom. Defaults to `nu = ncol(X)+1`. See Details. ## Returns Returns a list with elements `mu` and `Sigma` of sizes `c(nsamples, length(lambda))` and `c(dim(Psi), nsamples)`. ## Details The NIW distribution $p(\mu, \Sigma | \lambda, \kappa, \Psi, \nu)$ is defined as $$ \Sigma \sim W^{-1}(\Psi, \nu), \quad \mu | \Sigma \sim N(\lambda, \Sigma/\kappa). $$ The default value `kappa = 0` uses the Lebesque prior on $\mu$: $p(\mu) \propto 1$. The default value `Psi = 0` uses the scale-invariant prior on $\Sigma$: $p(\Sigma) \propto |\Sigma|^{-(\nu+d+1)/2}$. The default value `nu = ncol(X)+1` for `kappa = 0` and `Psi = 0` makes $E[\mu|X]=\code{colMeans(X)}$ and $E[\Sigma | X]=\code{var(X)}$. ## Examples ```r # compare the default non-informative prior to an arbitrary informative prior # for simulated data # simulate normal data with mean and variance (mu0, Sigma0) d <- 4 mu0 <- rnorm(d) Sigma0 <- matrix(rnorm(d^2), d, d) Sigma0 <- Sigma0 %*% t(Sigma0) N <- 1e2 X <- matrix(rnorm(N*d), N, d) # iid N(0,1) X <- t(t(X %*% chol(Sigma0)) + mu0) # each row is N(mu0, Sigma) # informative prior parameters lambda <- rnorm(d) kappa <- 20 Psi <- crossprod(matrix(rnorm(d^2), d, d)) nu <- 5 # iid draws from informative prior pi(mu, Sigma) nsamples <- 2e3 siw0 <- rniw(nsamples, lambda, kappa, Psi, nu) # iid draws from posterior p(mu, Sigma | X) with informative prior siw1 <- niw.post(nsamples, X, lambda, kappa, Psi, nu) # iid draws from posterior p(mu, Sigma | X) with default noninformative prior siw2 <- niw.post(nsamples, X) # compare # prior and posterior densities of mu clrs <- c("orange", "red", "blue", "black") ii <- 1 par(mar = c(4.2, 4.2, 2, 1)+.1) niche.par.plot(list(siw0, siw1, siw2), col = clrs[1:3], plot.index = ii, ylab = "Density") abline(v = mu0[ii], col = clrs[4]) # true value of mu legend(x = "topright", legend = c(parse(text = paste0("pi(mu[", ii, "])")), parse(text = paste0("p(mu[", ii, "]*\" | \"*X)*\", Informative Prior\"")), parse(text = paste0("p(mu[", ii, "]*\" | \"*X)*\", Noninformative Prior\"")), parse(text = paste0("\"True value of \"*mu[", ii, "]"))), fill = clrs) # prior and posterior densities of Sigma ii <- 1 jj <- 2 par(mar = c(4.2, 4.2, 2, 1)+.1) niche.par.plot(list(siw0, siw1, siw2), col = clrs[1:3], plot.index = c(ii,jj), ylab = "Density") abline(v = Sigma0[ii,jj], col = clrs[4]) legend(x = "topright", legend = c(parse(text = paste0("pi(Sigma[", ii, "*", jj, "])")), parse(text = paste0("p(Sigma[", ii, "*", jj, "]*\" | \"*X)*\", Informative Prior\"")), parse(text = paste0("p(Sigma[", ii, "*", jj, "]*\" | \"*X)*\", Noninformative Prior\"")), parse(text = paste0("\"True value of \"*Sigma[", ii, "*", jj, "]"))), fill = clrs) ``` ## See Also `rniw()`, `niiw.post()`.

niw.post function