REVIEW 4 major objections 5 minor 37 references
The acceptance-complement method can generate exact random variates from any log-concave density in one pass, with a fixed bound on run time.
Reviewed by Pith at T0; open to challenge. T0 means a machine referee read the full paper against a public rubric. the ladder, T0–T4 →
T0 review · deepseek-v4-flash
2026-08-01 07:56 UTC pith:7Q7WQSFZ
load-bearing objection A genuinely new universal acceptance-complement sampler for log-concave densities with bounded deterministic time, but several load-bearing inequalities are asserted rather than proved. the 4 major comments →
The acceptance-complement method revisited
The pith
A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.
Core claim
The central claim is that the acceptance-complement method, usually a special-purpose tool, can be made universal. For any log-concave density f with known mode m and a black-box evaluator of f, the paper builds a symmetric 16-interval grid around m, a lower step function r, and an upper envelope Q with 0 ≤ r ≤ f ≤ Q and ∫(Q−r) ≤ 1. One uniform draw then chooses between sampling from q = Q−r and sampling from r; both branches return an exact variate in one pass, no rejection loop. The bound on ∫q comes from a tail inequality for log-concave densities. The same machinery gives a loopless gamma generator for a ≥ 5, a gamma-function-free version for a ≥ 6.2829, and smaller shapes via G_a = U^{1
What carries the argument
The central object is a squeeze pair (r, Q) built from point evaluations of f on a grid centered at the mode: r uses inner grid values as a lower step function, Q uses outer grid values plus exponential tails as an upper envelope. The method samples from q = Q−r with probability A = ∫q and from r otherwise. For log-concave densities, the tail inequality f(x) ≤ f(m) exp(−(x−m)f(m)(1−exp(1−(x−m)f(m)))) keeps A below 1 for n = 7, δ = 2/5, giving a deterministic one-pass guarantee. For gamma, a four-piece envelope (two flats plus two exponential tails) plays the same role, with A < 1 for a ≥ 5; a representation A = E[φ(|N|/√(a−1))] removes the gamma function from the decision step for a ≥ 6.2829
Load-bearing premise
The bounded-time guarantee rests on the tail inequality f(x) ≤ f(m) exp(−(x−m)f(m)(1−exp(1−(x−m)f(m)))) for log-concave densities, whose proof is a terse extremal argument; if it fails, the 16-interval envelope's area can exceed 1, and separately the gamma algorithm for a ≥ 5 depends on an asserted inequality A < 1 that is stated without proof.
What would settle it
Evaluate both sides of the tail inequality numerically over a dense grid of (x−m)f(m) values for several log-concave densities (e.g., Laplace, logistic, gamma with shape just above 1) and look for a violation; equivalently, compute ∫q for the paper's grid with δ = 2/5, n = 7 on those densities and check whether it stays below 1. For the gamma claim, compute A in Eq. (2) for a = 5 and nearby values to verify A < 1, since the paper gives no derivation.
If this is right
- For any log-concave density, each generated variate costs one density evaluation plus bounded overhead, with 15 grid values precomputed once per density.
- Gamma variates of shape a ≥ 5 are generated by a loopless algorithm when Γ(a) is available; for a ≥ 6.2829, no special-function evaluation is needed, and smaller shapes are exact via G_a = U^{1/a}G_{a+1}.
- Beta, beta-prime, and variance-gamma variates become one-liners as deterministic functions of gamma and normal one-liners.
- Batch simulation from log-concave densities becomes competitive with rejection samplers, because the 16-point table is built once and per-sample work is constant.
- The same finite-partition construction extends to any unimodal density with known mode and explicit tail bounds, and to multimodal densities with known extrema.
Where Pith is reading between the lines
- The paper leaves open the smallest universal partition size; testing grids with fewer than 16 intervals on parametric log-concave families could reveal a tighter constant.
- The gamma-function-free trick—expressing a normalizing constant as an expectation over a simple random variable—may generalize to other densities with hard-to-compute normalizing constants, provided a bounded φ with E[φ(Z)] = A can be found.
- If the tail inequality admits a sharper form, the grid could be coarsened, making the universal sampler faster; conversely, a counterexample to the inequality would force a redesign of the 16-interval envelope.
- Because the algorithms are loopless, they are natural candidates for hardware or machine-code implementations, a consequence the paper notes but leaves undeveloped.
Editorial analysis
A structured set of objections, weighed in public.
Referee Report
Summary. The paper revisits the acceptance-complement method for random variate generation, showing that to sample from a density f it suffices to construct r ≤ f and q ≥ f − r with ∫q ≤ 1, yielding a one-iteration (loopless) sampler. It gives a four-part acceptance-complement generator for gamma(a) for a ≥ 5 (Algorithm 5), a gamma-function-free variant for a ≥ 6.2829 (Section 4), and a universal sampler for all log-concave densities with known mode and black-box density evaluations (Algorithm 7, Theorem 1) using a 16-interval grid. The central claim is that deterministic, uniformly bounded execution time is achievable for these families with a finite number of density evaluations per sample.
Significance. If the claims are correct, this is a useful conceptual advance: it replaces geometrically distributed rejection loops with constant-time acceptance-complement samplers, and it provides the first universal bounded-time method for log-concave densities. The paper is commendably explicit: Algorithms 5–7 are concrete, the constructions rely on a parameter-free tail bound from Devroye [5] rather than fitted constants, and no empirical tuning is involved. The main ideas are elegant and the spot checks in the reader's report support the numerical claims. However, several load-bearing inequalities are asserted without proof, and the proof of Theorem 1 is a sketch at exactly the point where correctness depends on a sharp tail bound and a numerical verification. These gaps are local and fillable, but they must be addressed before the paper is publishable.
major comments (4)
- [§5, Inequality (8)] The proof of Theorem 1 rests on the tail inequality (8), but its derivation is only a sketch. The extremal density 'when log f is linear ... and zero outside [m,x]' is discontinuous at the mode, so attainment requires a limiting discussion, and no rigorous argument is given. A recoverable route is to use concavity of log f to obtain the chord lower bound ∫_m^x f ≥ (x−m)(f(m)−f(x))/log(f(m)/f(x)), then use ∫ f = 1 and inequality (7) to derive (8). Since (9) and the bound on ∫q are direct consequences of (8), this point is load-bearing for Theorem 1 and Algorithm 7.
- [§5, after Eq. (9)] The proof asserts, without derivation, that 'The first n for which this bound drops below 1 for some x is n=7. Any choice 0.3467≤δ≤0.4643 will do. The minimal value is approximately 0.9540 and occurs for δ≈0.3968.' This numerical verification is exactly what guarantees ∫q ≤ 1 and hence the validity of Algorithm 7. Please provide a complete argument—for example, a monotonicity and calculus proof for n=7 and δ=2/5, or a certified numerical bound—rather than an unstated computer check.
- [§3, after Eq. (2)] Algorithm 5 is only valid if A < 1 on the stated range a ≥ 5. The text says 'Combining this, one can show that A < 1 for a ≥ 5' with no proof. If A ≥ 1, the proposed q does not satisfy ∫q ≤ 1 and the acceptance-complement construction fails. Because this inequality is load-bearing for the gamma algorithm, the proof must be supplied using the stated Stirling bound.
- [§4, after Eq. (6)] The gamma-function-free modification requires 0 ≤ φ ≤ 1. The text asserts the upper bound is '< 1 for a ≥ 6.2829' without proof. This inequality is load-bearing for the modified Algorithm 5 when Γ(a) is not available. Provide a derivation of the threshold and the bound, including the behavior of ψ(√(a−1)).
minor comments (5)
- [§5, displayed integral in proof of (8)] The lower limit in the displayed integral is written as 0; it should be m (or the mode should be normalized to 0 before the calculation).
- [§4, Eq. (4)] The text in (4) and the surrounding discussion state the condition as 'if W ≤ 1−A', while Algorithm 5 uses 'if W ≥ 1−A'. The replacement and the equivalence calculation should use W ≥ 1−φ(Z), giving probability E[φ(Z)] = A. As written, the directions are inconsistent and confusing.
- [Algorithm 6] The line 'generate I such that P{J=j}' mixes I and J; it should read 'generate J such that P{J=j}'.
- [§5, after Eq. (9)] The phrase 'for some x' in the numerical bound should be 'for some δ'; the independent variable in that bound is δ.
- [Algorithm 7] In the probability statement for J, define r_{n+1}=r_{−(n+1)}=0 explicitly in the algorithm itself, as is done in the surrounding text, to make the weights (Q_j−r_j)/A well defined.
Circularity Check
No circularity: constants are analytic bounds, and the cited log-concave inequality is independent support.
full rationale
No load-bearing step in the paper reduces to its input by construction or to a fitted value. The universal sampler is built from inequality (7), cited from Devroye [5]; this is a parameter-free theorem with stated assumptions (log-concavity, known mode) and does not presuppose the acceptance-complement result, so under the review rules it is independent support rather than a self-citation loop. The sharper tail inequality (8) is derived by an extremal integral argument, and the choices n=7 and delta=2/5 are fixed by an analytic upper-bound calculation (displayed bound ~0.954<1); they are not fitted to the target distribution or to any data. The gamma algorithms similarly use analytic inequalities A<1 and phi<1 to choose thresholds and are validated by the self-contained algebra of Section 7. The unproved/numerically asserted nature of some bounds, flagged by a skeptical reader, is a correctness or verification risk, not a circularity: no equation is defined in terms of the target output and no 'prediction' is forced by a fitted parameter.
Axiom & Free-Parameter Ledger
free parameters (2)
- δ (universal grid width multiplier) =
0.4 (2/5)
- n (number of grid intervals each side of mode) =
7
axioms (6)
- domain assumption Infinite-precision i.i.d. uniform[0,1] random variates are available.
- standard math The gamma density with a≥1 is log-concave and unimodal at a−1.
- domain assumption Inequality (7): f(x)≤f(m)min(1,e^{1−|x−m|f(m)}) for log-concave f.
- domain assumption Tail inequality (8): f(x)≤f(m)exp(−(x−m)f(m)(1−e^{1−(x−m)f(m)})) for (x−m)f(m)>1.
- standard math Reciprocal gamma function Hankel contour representation used to replace A by Eφ(Z).
- standard math A<1 for a≥5 and φ(Z)<1 for a≥6.2829.
read the original abstract
We revisit the acceptance-complement method in random variate generation and show how it can replace the rejection method in many examples. While the rejection method has geometrically distributed execution times, the acceptance-complement method has a constant (deterministic) run time and qualifies as a ``one-liner''. We show how this method can be used to efficiently generate random variates from several distributions, such as the gamma and beta. In addition, we show that there is an acceptance-complement method that is valid for all log-concave densities with known location of the mode and black-box type access to the density.
Reference graph
Works this paper leans on
-
[1]
J. H. Ahrens and U. Dieter. Computer generation of Poisson deviates from modified normal distributions.ACM Transactions on Mathematical Software, 8(2):163–170, 1982. doi: 10.1145/355993.355997. page19
arXiv 1982
-
[2]
R. W . Bailey. Polar generation of random variates with thetdistribution.Mathematics of Computation, 62:779–781, 1994
1994
-
[3]
J. M. Chambers, C. L. Mallows, and B. W . Stuck. A method for simulating stable random variables.Journal of the American Statistical Association, 71:340–344, 1976
1976
-
[4]
I. Deak. An economical method for random number generation and a normal generator. Computing, 27(2):113–121, 1981. doi: 10.1007/BF02243545
-
[5]
L. Devroye. A simple algorithm for generating random variates with a log-concave density. Computing, 33:247–257, 1984
1984
-
[6]
Devroye.Non-Uniform Random V ariate Generation
L. Devroye.Non-Uniform Random V ariate Generation. Springer-V erlag, New Y ork, 1986
1986
-
[7]
L. Devroye. Random variate generation in one line of code. In J.M. Charnes, D.J. Morrice, D.T. Brunner, and J.J. Swain, editors,1996 Winter Simulation Conference Proceedings, pages 265–272, San Diego, CA, 1996. ACM
1996
-
[8]
L. Devroye. Random variate generation for the generalized inverse Gaussian distribution. Statistics and Computing, 24:239–246, 2014
2014
-
[9]
W .R. Gilks. Derivative-free adaptive rejection sampling for Gibbs sampling. In J. Bernardo, J. Berger, A.P. Dawid, and A.F.M. Smith, editors,Bayesian Statistics 4. Oxford University Press, 1992
1992
-
[10]
Gilks and P
W .R. Gilks and P. Wild. Adaptive rejection sampling for Gibbs sampling.Applied Statistics, 41:337–148, 1992
1992
-
[11]
Gilks and P
W .R. Gilks and P. Wild. Algorithm as 287: Adaptive rejection sampling from log-concave density function.Applied Statistics, 41:701–709, 1993
1993
-
[12]
Gilks, N.G
W .R. Gilks, N.G. Best, and K.K.C. T an. Adaptive rejection Metropolis sampling.Applied Statistics, 44:455–472, 1995
1995
-
[13]
Extended one-liners for the gamma, poisson, and binomial distributions, 2026
Dylan Greaves. Extended one-liners for the gamma, poisson, and binomial distributions, 2026
2026
-
[14]
Dylan Greaves. Extended one-liners for the beta, gamma, and Dirichlet distributions with shape parameters below one.arXiv, 2604.11199, 2026
Pith/arXiv arXiv 2026
-
[15]
Hörmann, J
W . Hörmann, J. Leydold, and G. Derflinger.Automatic Nonuniform Random V ariate Generation. Springer-V erlag, Berlin, 2004
2004
-
[16]
M.C. Jones. Student’s simplest distribution.Journal of the Royal Statistical Society Series D, 51:41–49, 2002. page20
2002
-
[17]
M. Kanter. Stable densities under change of scale and total variation inequalities.Annals of Probability, 3:697–707, 1975
1975
-
[18]
Ekatherina A. Karatsuba. On the asymptotic representation of the Euler gamma function by Ramanujan.Journal of Computational and Applied Mathematics, 135(2):225–240, 2001. doi: 10.1016/S0377-0427(00)00586-0. MR 1850542
-
[19]
Kronmal and Jr
Richard A. Kronmal and Jr. Arthur V . Peterson. A variant of the acceptance-rejection method for the computer generatlon of random variables.Journal of the American Statistical Association, 76:446–451, 1981
1981
-
[20]
Kronmal and Jr
Richard A. Kronmal and Jr. Arthur V . Peterson. An acceptance-complement analogue of the mixture-plus-acceptance-rejection method for generating random variables.ACM Transactions on Mathematical Software, 10(3):271–281, 1984
1984
-
[21]
Leydold and W
J. Leydold and W . Hörmann. Black box algorithms for generating non-uniform continuous random variates. In W . Jansen and J.G. Bethlehem, editors,COMPSTAT 2000, pages 53–54, 2000
2000
-
[22]
Leydold and W
J. Leydold and W . Hörmann. Universal algorithms as an alternative for generating non- uniform continuous random variates. In G.I. Schuler and P.D. Spanos, editors,Monte Carlo Simulation, pages 177–183. 2001
2001
-
[23]
E.A. Luengo. Gamma pseudo-random number generators.ACM Computing Surveys, 55(4): 85, 2022
2022
-
[24]
A simple method for generating gamma variables
George Marsaglia and W ai W an T sang. A simple method for generating gamma variables. ACM Transactions on Mathematical Software, 26(3):363–372, 2000
2000
-
[25]
Cristinel Mortici. Improved asymptotic formulas for the gamma function.Computers & Mathematics with Applications, 61(11):3364–3369, 2011. doi: 10.1016/j.camwa.2011.04.036
-
[26]
Cristinel Mortici. Ramanujan’s estimate for the gamma function via monotonicity arguments.The Ramanujan Journal, 25(2):149–154, 2011. doi: 10.1007/s11139-010-9265-y
-
[27]
Cristinel Mortici. On Ramanujan’s large argument formula for the gamma function.The Ramanujan Journal, 26(2):185–192, 2011. doi: 10.1007/s11139-010-9281-y
-
[28]
A remark on Stirling’s formula.The American Mathematical Monthly, 62(1):26–29, 1955
Herbert Robbins. A remark on Stirling’s formula.The American Mathematical Monthly, 62(1):26–29, 1955. doi: 10.2307/2308012
doi:10.2307/2308012 1955
-
[29]
Schmeiser and R
B. Schmeiser and R. Lal. Squeeze methods for generating gamma variates.Journal of the American Statistical Association, 75:679–682, 1980. page21
1980
-
[30]
D. T eichroew. The mixture of normal distributions with different variances.The Annals of Mathematical Statistics, 28(2):510–512, 1957. doi: 10.1214/aoms/1177706981. URLhttps: //doi.org
arXiv 1957
-
[31]
T emme.Special Functions: An Introduction to the Classical Functions of Mathematical Physics
Nico M. T emme.Special Functions: An Introduction to the Classical Functions of Mathematical Physics. Wiley, New Y ork, 1996
1996
-
[32]
V arious techniques used in connection with random digits.Collected W orks, 5:768–770, 1963
John von Neumann. V arious techniques used in connection with random digits.Collected W orks, 5:768–770, 1963. Also in Monte Carlo Method, National Bureau of Standards Series, vol. 12, pp. 36-38, 1951
1963
-
[33]
A. J. W alker. New fast method for generating discrete random numbers with arbitrary frequency distributions.Electronics Letters, 10(8):127–128, 1974. doi: 10.1049/el:19740097. URLhttps://doi.org
-
[34]
A. J. W alker. An efficient method for generating discrete random variables with general distributions.ACM Transactions on Mathematical Software, 3(3):253–256, 1977. doi: 10. 1145/355744.355749. URLhttps://doi.org/10.1145/355744.355749
arXiv 1977
-
[35]
B. Xi, K.M. T an, and C. Liu. Logarithmic transformation-based gamma random number generators.Journal of Statistical Software, 55(4):1–17, 2013
2013
-
[36]
V . M. Zolotarev. On the representation of stable laws by integrals.Selected Translations in Mathematical Statistics and Probability, 6:84–88, 1966
1966
-
[37]
V . M. Zolotarev.One-Dimensional Stable Distributions. American Mathematical Society, Providence, R.I., 1986
1986
discussion (0)
Sign in with ORCID, Apple, or X to comment. Anyone can read and Pith papers without signing in.