findMode {esaddle}

Finding the mode of the empirical saddlepoint density

Description

Given a sample from a d-dimensional distribution, the routine finds the mode of the corresponding Extended Empirical Saddlepoint (EES) density.

Usage

findMode(
  X,
  decay,
  init = NULL,
  method = "BFGS",
  hess = FALSE,
  sadControl = list(),
  ...
)

Arguments

X

an n by d matrix containing the data.

decay

rate at which the SPA falls back on a normal density. Should be a positive number. See Fasiolo et al. (2016) for details.

init

d-dimensional vector containing the starting point for the optimization. By default it is equal to colMeans(X).

method

optimization method used by stats::optim(), see ?optim for details. By default it is "BFGS".

hess

if TRUE also an estimate of the Hessian at the mode will be returned.

sadControl

list corresponding to the control argument in the dsaddle function.

...

Extra arguments to be passed to the optimization routine stats::optim.

Value

A list where mode is the location of mode of the empirical saddlepoint density, logDens is the log-density at the mode and hess (present only if argument hess==TRUE) is the approximate Hessian at convergence. The other entries are the same as for stats::optim.

Author(s)

Matteo Fasiolo <[e-mail]>.

References

Fasiolo, M., Wood, S. N., Hartig, F. and Bravington, M. V. (2016). An Extended Empirical Saddlepoint Approximation for Intractable Likelihoods. ArXiv http://arxiv.org/abs/1601.01849.

Examples

# library(esaddle)
set.seed(4141)
x <- rgamma(1000, 2, 1)

# Fixing tuning parameter of EES
decay <-  0.05

# Evaluating EES at several point
xSeq <- seq(-2, 8, length.out = 200)
tmp <- dsaddle(y = xSeq, X = x, decay = decay, log = TRUE)  # Un-normalized EES

# Plotting true density, EES and normal approximation
plot(xSeq, exp(tmp$llk), type = 'l', ylab = "Density", xlab = "x")
lines(xSeq, dgamma(xSeq, 2, 1), col = 3)
suppressWarnings( rug(x) )
legend("topright", c("EES", "Truth"), col = c(1, 3), lty = 1)

# Find mode and plot it
res <- findMode(x, init = mean(x), decay = decay)$mode
abline(v = res, lty = 2, lwd = 1.5)