Efficient and Accurate Parallel Inversion of the Gamma Distribution

Abstract

A method for parallel inversion of the gamma distribution is described. This is very desirable for random number generation in Monte Carlo simulations where gamma variates are required. Let α be a fixed but arbitrary positive real number. Explicitly, given a list of uniformly distributed random numbers our algorithm applies the quantile function (inverse CDF) of the gamma distribution with shape parameter α to each element. The result is, therefore, a list of random numbers distributed according to the said distribution. The output of our algorithm has accuracy close to a choice of single-or double-precision machine epsilon. Inversion of the gamma distribution is traditionally accomplished using some form of root finding. This is known to be computationally expensive. Our algorithm departs from this paradigm by using an initialization phase to construct, on the fly, a piecewise Chebyshev polynomial approximation to a transformation function, which can be evaluated very quickly during variate generation. The Chebyshev polynomials are high order, for good accuracy, and generated via recurrence relations derived from nonlinear second order ODEs. A novelty of our approach is that the same change of variable is applied to each uniform random number before evaluating the transformation function. This is particularly amenable to implementation on SIMD architectures, whose performance is sensitive to frequently diverging execution flows due to conditional statements (branch divergence). We show the performance of a CUDA GPU implementation of our algorithm (called Quantus) is within an order of magnitude of the time to compute the normal quantile function.

Quantile mechanics: Issues arising from critical review

“…The general quality of results, convergence to the exact values and high precision of the QM approach in quantile approximation presents the method as a major improvement over the others. This was the outcome in the approximation of the student's t distribution (Shaw and Brickman, 2009), normal distribution (Luu, 2016) and gamma distribution (Luu, 2015;.…”

“…The QM approach has been used to develop numerical algorithms for quantile approximation, which was found to be efficient, robust, and fast and save computing time. For example the algorithm obtained for the normal QF is found to be better than those available in literature, see Shaw and Brickman (2009) and Luu (2015). The details on this aspect were discussed extensively in Luu (2016).…”

“…As noted by Luu (2016) that quantile approximations of distributions with shape parameters are cumbersome and require two processes. The QM approach performed better than the Cornish Fisher expansion method in quantile approximation of the student's t distribution (Shaw and Brickman, 2009;Shaw and McCabe, 2009), gamma distribution (Steinbrecher and Shaw, 2008;Luu, 2015; and beta distribution (Steinbrecher and Shaw, 2008).…”

Efficient procedure to generate generalized Gaussian noise using linear spline tools

“…This effect is useful for sensitivity analysis. Contrast this to rejection method, see [17]. Mention must be made that some researchers avoid using this method for simulation, because for many distribution functions we do not have an explicit expression for the inverse of cumulative distribution function.…”

Approximations for the inverse cumulative distribution function of the gamma distribution used in wireless communication

“…The solutions are given per value of the shape parameter because distributions with shape parameters are often difficult in obtaining their approximation [ 42 ]. The existing approximations are often slow [ 43 , 44 ], plagued with slow convergence [ 45 , 46 , 47 ] and cumbersome in dealing with the extreme tails of complex distributions [ 48 , 49 , 50 , 51 , 52 ]. On the other hand, the second-order nonlinear ODE generated using QM is often complex.…”