emix: Compute expectations via weighted mixtures

In jmhewitt/bisque: Approximate Bayesian Inference via Sparse Grid Quadrature Evaluation (BISQuE) for Hierarchical Models

Description

Approximates expectations of the form

E[h(θ)] = \int h(θ) f(θ) dθ

using a weighted mixture

E[h(θ)] \approx ∑_{j=1}^k h(θ^{(k)}) w_k

Usage

emix(h, params, wts, ncores = 1, errorNodesWts = NULL, ...)

Arguments

h	Function for which the expectation should be taken. The function should be defined so it is can be called via f(params, ...). Additional parameters may be passed to h via ....
params	Matrix in which each row contains parameters at which h should be evaluated. The number of rows in params should match the number of mixture components k.
wts	vector of weights for each mixture component
ncores	number of cores over which to evaluate mixture. this function assumes a parallel backend is already registered.
errorNodesWts	list with elements inds and weights that point out which params get used to compute an approximation of the quadrature error.
...	additional arguments to be passed to h

Examples

# density will be a mixture of betas
params = matrix(exp(2*runif(10)), ncol=2)

# mixture components are equally weighted
wts = rep(1/nrow(params), nrow(params))

# compute mean of distribution by cycling over each mixture component
h = function(p) { p[1] / sum(p) }

# compute mixture mean
mean.mix = emix(h, params, wts)

# (comparison) Monte Carlo estimate of mixture mean
nsamples = 1e4
component = sample(x = 1:length(wts), size = nsamples, prob = wts, 
                   replace = TRUE)
x = sapply(component, function(cmp) {
  rbeta(n = 1, shape1 = params[cmp, 1], shape2 = params[cmp, 2])
})
mean.mix.mc = mean(x)

# compare estimates
c(emix = mean.mix, MC = mean.mix.mc)

jmhewitt/bisque documentation built on Feb. 9, 2020, 2:36 a.m.