ars: Adaptive Rejection Sampling

Description Usage Arguments Details Value References Examples

Description

Generate samples from a distribution via the adapative rejection sampling algorithm. Adaptive rejection sampling dynamically builds a “rejection region,” so that the user does not need to explicitly supply one.

Usage

1
2
ars(n, g, dg = NULL, initialPoints = NULL, leftbound = -Inf,
  rightbound = Inf)

Arguments

n

Number of samples.

g

(Unnormalized) density function to sample from. g must be (non-strictly) log-concave.

dg

Derivative of the density function. If not supplied, numeric differentiation is attempted.

initialPoints

A vector of the initial abscissae to generate the envelope and squeezing functions. If not supplied, an optimization routine will attempt to find initial points.

leftbound

The lower bound of the domain of g.

rightbound

The upper bound of the domain of g.

Details

The function g must be an log-concave (unnormalized) density function. Accurate leftbound rightbound values must also be supplied.

The initialPoints are used to construct the envelope and squeezing functions described in the paper (Gilks & Wild, 1992) below. If not supplied, an algorithm is run which attempts to find the mode of g and generates initial points near the mode. If the user wishes to supply initial points (recommended), at least one initial point must be given where the derivative of g is less than zero, and another where the derivative of g is greater than zero, unless the function is monotonic.

Value

Returns a vector of n samples from g.

References

Gilks, W. R., & Wild, P. (1992). Adaptive rejection sampling for Gibbs sampling. Applied Statistics, 337-348.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
# Sample 10 points from Normal(0,1)
ars(10, dnorm, initialPoints=c(-1,1))

# Sample 15 points from Uniform[0,1]
ars(15, dunif, initialPoints=c(0.2, 0.3, 0.8), leftbound=0, rightbound=1)

# Define a quadratic distribution
f <- function (x) {
 return(x^2)
}

df <- function (x) {
 return(2 * x)
}

# Sample 5 points from a quadratic distribution
# Note that ars() takes care of normalization
ars(5, f, df, initialPoints = c(1,2), leftbound=0, rightbound=10)

dylandaniels/ars documentation built on May 15, 2019, 7:23 p.m.