integrate.gh {ecoreg}

Univariate Gauss-Hermite integration

Description

Computes the integral of a univariate function, or several univariate functions simultaneously, using Gauss-Hermite quadrature.

Usage

integrate.gh(h, n=1, points = 10, mu = 0, scale = 1, ...)

Arguments

h

The function to be integrated. May either have a scalar first argument and return a scalar result, or have a first argument of length n and return a vector of n results, corresponding to n independent functions.

n

The dimension of the result returned by h.

points

Number of Gauss-Hermite quadrature points.

mu

Mode of the function, to centre the quadrature points around.

scale

Scale of the quadrature points.

...

Other arguments to be passed to h.

Details

The integral is more accurate if the standard quadrature points are shifted and scaled to match the mode and scale of g(x), that is the objective function divided by the standard normal density. The scale is estimated by 1/\sqrt{-H}, where H is the Hessian at the maximum of g(x).

Value

The integral of h(x) between -Inf and Inf, of length n. In the usual application of Gauss-Hermite quadrature, h(x) is equivalent to a function g(x)\phi(x), where \phi(x) is the standard normal density function.

Author(s)

C. H. Jackson [e-mail]

The Gauss-Hermite polynomial values and weights are calculated using the gauss.hermite function copied from the rmutil package by J. K. Lindsey.

References

Liu, Q. and Pierce, D. A. (1994) A note on Gauss-Hermite quadrature. Biometrika, 81 (624-629)

See Also

gauss.hermite

Examples

## Want the integral of h over the real line
g <- function(x) 4 * exp( - ((1 - x)^2 + 1))
h <- function(x) g(x) * dnorm(x)
integrate(h, -Inf, Inf)
integrate.gh(h)
## Not very accurate with default 10 points. Either use more quadrature points, 
integrate.gh(h, points=30)
## or shift and scale the points.
opt <- nlm(function(x) -g(x), 0, hessian=TRUE)
integrate.gh(h, mu=opt$estimate, scale=1/sqrt(opt$hessian))