BernoulliFactory.jl
This package provides some Bernoulli factory algorithms, and a constrained, unbiased estimation algorithm.
Given a function g such that g() returns a Bernoulli(p) variate, and a function f mapping [0,1] -> [0,1] a Bernoulli factory algorithm should use calls to g produce a Bernoulli(f(p)) variate. Not all functions f admit an algorithm (Keane & O'Brien, 1994), and some algorithms require additional information.
Currently implemented
| Command | Description |
|---|---|
| f(p) = exp(-λp), λ >= 0 | expMinus(g, λ) |
| f(p) = p^a, a >= 0 | power(g, a) |
| f(p) = C*p / (1+C*p), C > 0 | logistic(g, C) |
| f(p) = C*p, ϵ ∈ (0, 1-C*p) | linear(g, C, ϵ) |
| f(p) = C/p, ϵ ∈ (0, p-C) | inverse(g, C, ϵ) |
| f(p1, p2) = C1*p1 / (C1*p1 + C2*p2) | twocoin(g1, g2, C1, C2) |
Constrained, unbiased estimation
Also implemented is a variant of signed estimation, signedEstimate(μ, φ, a, b, δ, c, n):
Let X ~ μ and real-valued φ satisfy
a <= inf_x φ(x) < b < δ <= E[φ(X)] < sup_x φ(x) <= c,
with known (a, b, δ, c). If simulation from μ and evaluation of φ is feasible, one can simulate W satisfying
- E[W] = E[φ(X)]
- Pr(b <= W <= max{2b-a,c}) = 1
The random variable W is the product of two independent random variables, X and Z, and the positive integer parameter n specifies a number of averages to use to define X.
Algorithms used
expMinus is an obvious extension of the approach described in Wastlund (1999) for λ=1.
power is from Mendo (2016).
logistic is from Huber (2017).
linear is from Huber (2016), and for a large range of parameter settings appears to use the least expected flips of all algorithms for this problem.
inverse is described in Lee, Doucet & Łatuszyński (2014).
twocoin is described in Gonçalves, Łatuszyński & Roberts (2017).
signedEstimate is described in Appendix C of Lee, Doucet & Łatuszyński (2014).
Usage
All Bernoulli factory algorithms return a tuple (X, flips) where X is true or false and flips is the number of calls of g() by the algorithms.
signedEstimate returns a tuple (X, flips, calls) where X is true or false, flips is the number of Bernoulli factory flips and calls is the number of calls of μ() by the algorithms. calls is typically much larger thatn flips
References
Gonçalves, F.B., Łatuszyński, K.G. and Roberts, G.O., 2017. Exact Monte Carlo likelihood-based inference for jump-diffusion processes. arXiv:1707.00332
Huber, M., 2016. Nearly optimal Bernoulli factories for linear functions. Combinatorics, Probability and Computing, 25(4), pp.577-591.
Huber, M., 2017. Optimal linear Bernoulli factories for small mean problems. Methodology and Computing in Applied Probability, 19(2), pp.631-645.
Keane, M.S. and O'Brien, G.L., 1994. A Bernoulli factory. ACM Transactions on Modeling and Computer Simulation (TOMACS), 4(2), pp.213-219.
Lee, A., Doucet, A. and Łatuszyński, K., 2014. Perfect simulation using atomic regeneration with application to Sequential Monte Carlo. arXiv:1407.5770
Mendo, L., 2016. An asymptotically optimal Bernoulli factory for certain functions that can be expressed as power series. arXiv:1612.08923
Wästlund, J., 1999. Functions arising by coin flipping. Technical Report, KTH, Stockholm.