BayesMultiMode (version 0.7.4)

bayes_mode: Bayesian mode inference

Description

Bayesian inference on the modes in a univariate mixture estimated with MCMC methods, see Cross2024;textualBayesMultiMode. Provides posterior probabilities of the number of modes and their locations. Under the hood it calls the function mix_mode() to find the modes in each MCMC draw.

Usage

bayes_mode(
  BayesMix,
  rd = 1,
  tol_mixp = 0,
  tol_x = sd(BayesMix$data)/10,
  tol_conv = 1e-08,
  inside_range = TRUE,
  range = c(min(BayesMix$data), max(BayesMix$data)),
  conditional_nb_modes = NULL
)

Value

A list of class bayes_mode containing:

data

From BayesMix.

dist

From BayesMix.

dist_type

From BayesMix.

pars_names

From BayesMix.

modes

Matrix with a row for each draw and columns showing modes.

p1

Posterior probability of unimodality.

p_nb_modes

Matrix showing posterior probabilities for the number of modes.

p_mode_loc

Matrix showing posterior probabilities for mode locations.

mix_density

Mixture density at all mode locations in each draw.

algo

Algorithm used for mode estimation.

range

Range outside which modes are discarded if inside_range is TRUE.

conditional_nb_modes

From BayesMix.

BayesMix

BayesMix.

Arguments

BayesMix

An object of class bayes_mixture generated with either bayes_fit() or bayes_mixture().

rd

(for continuous mixtures) Integer indicating the number of decimal places when rounding the distribution's support. It is necessary to compute posterior probabilities of mode locations.

tol_mixp

Components with a mixture proportion below tol_mixp are discarded when estimating modes; note that this does not apply to the biggest component so that it is not possible to discard all components; should be between 0 and 1; default is 0.

tol_x

(for continuous mixtures) Tolerance parameter for distance in-between modes; default is sd(data)/10 where data is the vector of observations from BayesMix. If two modes are closer than tol_x, only the first estimated mode is kept.

tol_conv

(for continuous mixtures) Tolerance parameter for convergence of the algorithm; default is 1e-8.

inside_range

Should modes outside of range be discarded? Default is TRUE.

range

limits of the support where modes are saved (if inside_range is TRUE);

conditional_nb_modes

Mcmc draws are filtered to include those with only conditional_nb_modes number of modes; default is c(min(BayesMix$data), max(BayesMix$data)). This sometimes occurs with very small components when K is large.

Details

Each draw from the MCMC output after burnin, \(\theta^{(d)}, \quad d = 1,...,D\), leads to a posterior predictive probability density/mass function: $$p(y | \theta^{(d)}) =\sum_{k=1}^{K} \pi_k^{(d)} p(y | \theta_k^{(d)}).$$ Using this function, the mode in draw \(d\) \(y_{m}^{(d)}\), \(m = 1,..., M^{(d)}\), where \(M^{(d)}\) is the number of modes, are estimated using the algorithm mentioned in the description above.

After running this procedure across all retained posterior draws, we compute the posterior probability for the number of modes being \(M\) as: $$P(\#\text{modes}=M)=\frac{1}{D}\sum_{d=1}^{D}1(M^{(d)} = M).$$ Similarly, posterior probabilities for locations of the modes are given by: $$P(y=\text{mode})=\frac{1}{D}\sum_{d=1}^{D} \sum_{m=1}^{M^{(d)}} 1(y = y_m^{(d)}),$$ for each location \(y\) in the range \([\min(y),\max(y)]\). Obviously, continuous data are not defined on a discrete support; it is therefore necessary to choose a rounding decimal to discretize their support (with the rd argument).

References

Examples

# Example with galaxy data ================================================
set.seed(123)

# retrieve galaxy data
y <- galaxy

# estimation
bayesmix <- bayes_fit(
  data = y,
  K = 5, # not many to run the example rapidly
  dist = "normal",
  nb_iter = 500, # not many to run the example rapidly
  burnin = 100
)

# mode estimation
BayesMode <- bayes_mode(bayesmix)

# plot
# plot(BayesMode, max_size = 200)

# summary
# summary(BayesMode)

# Example with DNA data ================================================
set.seed(123)

# retrieve DNA data
y <- d4z4

# estimation
bayesmix <- bayes_fit(
  data = y,
  K = 5, # not many to run the example rapidly
  dist = "shifted_poisson",
  nb_iter = 500, # not many to run the example rapidly
  burnin = 100
)

# mode estimation
BayesMode <- bayes_mode(bayesmix)

# plot
# plot(BayesMode, max_size = 200)

# summary
# summary(BayesMode)

# Example with a Student t ================================================
mu <- c(0.5, 6)
sigma <- c(1, 2)
nu <- c(5, 5)
p <- c(0.8, 0.2) #'
data <- c(
  sn::rst(p[1] * 1000, mu[1], sigma[1], nu = nu[1]),
  sn::rst(p[2] * 1000, mu[2], sigma[2], nu = nu[2])
)

fit <- c(eta = p, mu = mu, sigma = sigma, nu = nu, xi = c(0, 0))
fit <- rbind(fit, fit)

pdf_func <- function(x, pars) {
  sn::dst(x, pars["mu"], pars["sigma"], pars["xi"], pars["nu"])
}

dist_type <- "continuous"

bayesmix <- bayes_mixture(fit, data,
  burnin = 1,
  pdf_func = pdf_func, dist_type = dist_type, loc = "mu"
)

BayesMode <- bayes_mode(bayesmix)

# plot
# plot(BayesMode, max_size = 200)

# summary
# summary(BayesMode)