A log conconcave function is majorized with a piecewise envelop, which on the original scale is piecewise exponential. As the resulting extremely precise envelop adapts, the rejection rate dramatically decreases.

Envelop(lines::Vector{Line}, support::Tuple{Float64, Float64})

A piecewise linear function with k segments defined by the lines L_1, ..., L_k and cutpoints c_1, ..., c_k+1 with c1 = support[1] and c2 = support[2]. A line Lk is active in the segment [ck, ck+1], and it's assigned a weight wk based on exp_integral. The weighted integral over c1 to ck+1 is one, so that the envelop is interpreted as a density.

Line(slope::Float64, intercept::Float64)

Basic ensamble-unit for an envelop.

Objective(logf::Function, support:)
Objective(logf::Function, grad::Function)

Convenient structure to store the objective function to be sampled. It must receive the logarithm of f and not f directly. It uses automatic differentiation by default, but the user can provide the derivative optionally.

RejectionSampler(f::Function, support::Tuple{Float64, Float64}, init::Tuple{Float64, Float64}) RejectionSampler(f::Function, support::Tuple{Float64, Float64}[ ,δ::Float64]) An adaptive rejection sampler to obtain iid samples from a logconcave function supported in support = (support[1], support[2]). f can either be the either probability density of the function to be sampled, or its logarithm. For the latter, use the keyword argument log=true. The functions can be unnormalized in the sense that the probability density can be specified up to a constant. The adaptive rejection samplings algorithm requires two initial points init = init[1], init[2] such that (d/dx)logp(init[1]) > 0 and (d/dx)logp(init[2]) < 0. These points can be provided directly (typically, any point left and right of the mode will do). It is also possibe to specify and search range and delta for a greedy search of the initial points. The support must be of the form (-Inf, Inf), (-Inf, a), (b, Inf), (a,b), and it represent the interval in which f has positive value, and zero elsewhere.

The alternative constructor uses a search_range, a value δ for the distance between points in the search, and min/max slope values in absolute terms.

Keyword arguments

• max_segments::Int = 10 : max size of envelop, the rejection-rate is usually slow with a small number of segments
• max_failed_factor::Float64 = 0.001: level at which throw an error if one single sample has a rejection rate exceeding this value
• logdensity::Bool = false: indicator fo whether f is the log of the probability density, up to a normalization constant.
• search_range::Tuple{Float64,Float64} = (-10.0, 10.0): range in which to search for initial points
• min_slope::Float64 = 1e-6: minimum slope in absolute value of logf at the initial/found points
• `max_slope::Float64 = Inf: maximum slope in absolute value of logf at the initial/found points

add_segment!(e::Envelop, l::Line)

Adds a new line segment to an envelop based on the value of its slope (slopes must be decreasing always in the envelop). The cutpoints are automatically determined by intersecting the line with the adjacent lines.

eval_envelop(e::Envelop, x::Float64)

Eval point a point x in the piecewise linear function defined by e. Necessary for evaluating the density assigned to the point x.

exp_integral(l::Line, x1::Float64, x2::Float64)

Computes the integral $LaTeX \int_{x_1} ^ {x_2} \exp\{ax + b\} dx.$ The resulting value is the weight assigned to the segment [x1, x2] in the envelop

intersection(l1::Line, l2::Line)

Finds the horizontal coordinate of the intersection between lines

run_sampler!(rng::AbstractRNG, sampler::RejectionSampler, n::Int)
run_sampler!(sampler::RejectionSampler, n::Int)

It draws n iid samples of the objective function of sampler, and at each iteration it adapts the envelop of sampler by adding new segments to its envelop.

sample_envelop(rng::AbstractRNG, e::Envelop)
sample_envelop(e::Envelop)

Samples an element from the density defined by the envelop e with it's exponential weights. See Envelop for details.