niw.coeffs: Posterior coefficients of the Normal-Inverse-Wishart...

In mlysy/nicheROVER: Niche Region and Niche Overlap Metrics for Multidimensional Ecological Niches

Posterior coefficients of the Normal-Inverse-Wishart distribution with its conjugate prior.

Description

Given iid d-dimensional niche indicators X = (X_1,\ldots,X_N) with X_i \sim N(\mu, \Sigma), this function calculates the coefficients of the Normal-Inverse-Wishart (NIW) posterior p(\mu, \Sigma | X) for a conjugate NIW prior. Together with niw.mom(), this can be used to rapidly compute the point estimates E[\mu | X] and E[\Sigma | X].

Usage

niw.coeffs(X, lambda, kappa, Psi, nu)

Arguments

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'.

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]=`colMeans(X)` and E[\Sigma | X]=`var(X)`.

Value

Returns a list with elements lambda, kappa, Psi, nu corresponding to the coefficients of the NIW posterior distribution p(\mu, \Sigma | X).

See Also

rniw(), niw.mom(), niw.post().

Examples

# NIW prior coefficients
d <- 3
lambda <- rnorm(d)
kappa <- 5
Psi <- crossprod(matrix(rnorm(d^2), d, d))
nu <- 10

# data
data(fish)
X <- fish[fish$species == "ARCS",2:4]

# NIW posterior coefficients
post.coef <- niw.coeffs(X, lambda, kappa, Psi, nu)

# compare
mu.mean <- niw.mom(post.coef$lambda, post.coef$kappa, post.coef$Psi, post.coef$nu)$mu$mean
mu.est <- rbind(prior = niw.mom(lambda, kappa, Psi, nu)$mu$mean,
                data = colMeans(X),
                post = mu.mean)
round(mu.est, 2)

mlysy/nicheROVER documentation built on Oct. 17, 2023, 2:21 p.m.