# Complexity of generating non-uniform random variates

What can we say about the complexity of generating (negative) binomial and (negative) hypergeometric random variates? In particular, it is possible to generate (negative) binomial and (negative) hypergeometric variates in (expected) constant time (i.e. independent of the distribution parameters)?

There is quite a bunch of literature; however, it’s hard to understand a lot of the papers. Moreover, I found some statements that seem contradictory to me (probably due to wrong understanding).

For example, Stadlober (or similar here) mentions a "generalization [of the ratio-of-uniforms approach] to any unimodal discrete distribution". The ratio-of-uniforms approach has been called uniformly fast, which is a synonym for constant time, I suppose (?). That would mean that we can generate random variates for each discrete distribution in (expected) constant time.

However, in another paper, I found the following theorem: On a RAM with word size $$w$$, any algorithm sampling a geometric random variate $$Geo(p)$$ with parameter $$p \in (0,1)$$ needs at least expected runtime $$\Omega(1 + \log(1/p)/w)$$.

Apparently, this means that it is not possible to generate negative binomial variates in time independent of the distribution parameters.

• It's very unlikely that "uniformly fast" means "constant time". "Constant time" is such a widely understood and widely known phrase that nobody's going to invent a new name for it (either because they don't know that there's already a standard name or because they think there's something wrong with the standard name). Another indication that it doesn't mean "constant time" is that you've found some papers that claim that it can be done "uniformly fast" and some that say it can't be done in constant time! Jan 18 '19 at 14:48
• @DavidRicherby From the first paper/link: "The acceptance/rejection approach of Von Neumann leads to uniformly fast algorithms, i.e., algorithms with bounded execution time over the defined parameter range." I'm not quite sure how to interpret this sentence (that's why I wrote " I suppose"), but to me it seems like the running time is bounded (i.e. does not grow infinitely, i.e. is constant in O-notation) no matter what distribution parameters we choose. Jan 18 '19 at 15:04
• @DavidRicherby By now, I'm relatively sure that uniformly fast means constant time, since Lemma 4.11 in nrbook.com/devroye/Devroye_files/chapter_ten.pdf states that an algorithm is uniformly fast, and in the proof they show O(1). Jan 18 '19 at 17:17
• @DavidRicherby And actually, the second paper mentioned in the question does allow constant time variate generation, as long as the word size is logarithmic in the mean of the distribution. I guess it's a common assumption in algorithm analysis that all the numbers "we deal with" fit into a single word, isn't it? Jan 18 '19 at 17:24
• OK -- seems very strange to make up a new term for "constant time" but you've read the papers and I haven't! And, yes, numbers fitting in a word is a common assumption, e.g., in the word RAM model as you mention. Jan 18 '19 at 17:31

The binomial, negative binomial, and hypergeometric distributions are discrete distributions. Knuth and Yao (1976) gave complexity results for discrete distributions in general. Given a stream of i.i.d. unbiased random bits, any algorithm that samples from a discrete distribution will require at least—$$b = -\sum_{x\in\Omega} \mathbb{P}(X=x)*\log_2(\mathbb{P}(X=x))$$random bits on average to sample a variate $$X$$ from that distribution (where $$\Omega$$ is the distribution's support), and they gave an algorithm that comes within 2 bits of this lower bound (and showed that coming within 2 bits is the best possible for any algorithm).
I also want to note one more thing. If the distribution of $$X$$ has a continuous and bounded quantile function (inverse cumulative distribution function) $$F^{-1}$$, I suspect that the bit complexity of sampling $$X$$ also depends on the modulus of continuity of $$F^{-1}$$, call it $$\omega(h)$$. In this sense, Lipschitz continuous quantile functions (where $$\omega(h) = O(h)$$) are the simplest cases of such functions and in that case it's simple to describe an algorithm to sample $$X$$ to a desired accuracy. For such functions, the algorithm uses—$$ceil(\ln(\max(1, L))-\ln(\epsilon))/\ln(\beta))$$random digits, where $$L$$ is the Lipschitz constant, $$\epsilon$$ is the desired accuracy, and $$\beta$$ is the digit base (such as 2 for binary). I suspect that with other moduli of continuity, the complexity will be as follows: $$ceil(-\ln(\omega^{-1}(\epsilon))/\ln(\beta)).$$