4
$\begingroup$

Is there an efficient method of sampling an n-choose-k combination at random (with uniform probability, for example)?

I have read this question but it asks for generations of all combinations, not combinations at random.

I general I'm aware of rejection sampling, however it's very inefficient.

I also came across reservoir sampling, but that appears to be primarily geared towards very large or unknown n. I'm more interested in large but finite n (definitely not large enough to not be able to fit in memory. Well. An n-sized array itself will fit in memory, but the state space of all n-choose-k combinations might not).

Is there any survey/review on this topic? Does Knuth cover random n-choose-k sampling in his TAOCP texts?

Thanks in advance.

Edit: To be a bit more specific, a 5-choose-3 space over the string 'ABCDE' would look like this:

['ABC', 'ABD', 'ABE', 'ACD', 'ACE', 'ADE', 'BCD', 'BCE', 'BDE', 'CDE']

(Note: this is combination without replacement). And I want to be able to sample from this space with uniform distribution, using a general algorithm (one that works with arbitrary n and k).

$\endgroup$
7
  • $\begingroup$ oeis.org/wiki/Ranking_and_unranking_functions $\endgroup$
    – D.W.
    Feb 27, 2019 at 19:40
  • $\begingroup$ @D.W. the ranking function is a nice approach. If "the state space of all n-choose-k combinations might not" fit in the memory, however, it might not be an immediate hit. $\endgroup$
    – John L.
    Feb 27, 2019 at 20:01
  • 2
    $\begingroup$ Reservoir sampling solves exactly this problem in $O(n)$ time. Can you elaborate why it won't work for you? $\endgroup$
    – orlp
    Feb 27, 2019 at 20:28
  • $\begingroup$ @Apass.Jack, for what it's worth, ranking/unranking functions don't require storing the set of all such combinations (just an index into that set). In any case, orlp has a reasonable answer, so it's moot here. $\endgroup$
    – D.W.
    Feb 27, 2019 at 20:56
  • $\begingroup$ @D.W. I meant the index might be out of the memory. For example, if $n=1000000$ and $k=100000$. $\endgroup$
    – John L.
    Feb 27, 2019 at 21:13

3 Answers 3

4
$\begingroup$

Here is the simplest algorithm, which is efficient when $k$ is much smaller than $n$ relatively.


Input: two positive integers $n$ and $k$ with $k\le n$
Output: a random permutation of $k$ integers from $1,2,\cdots,n$
Algorithm:

  1. Let arr be an array of size $n$ and a default value that is not True.
  2. Let out be an empty array.
  3. Let i be a random integer in $[0,n)$. If arr[i] is not True, append i+1 to out and set arr[i] to True.
  4. Go back to 3 unless we have selected $k$ elements.
  5. return out.

If $k\le n/2$, then the algorithm above runs in $O(k)$ average time. For example, if $n=1000$ and $k=50$, it will use the random number generator less than 55 times in average. This is much better the reservoir sampling that will use the random number generator about 1000 times.


How to speed up the algorithm if $k$ is not much smaller than $n$?

We will let each element in arr be a pair (i, False). After we have selected $n/4$ elemented, we will compactify arr by removing the elements whose second entries have been changed to True. We will set $n$ to $n-n/4$ and $k$ to $k-n/4$. Repeat the algorithm.

$\endgroup$
4
  • $\begingroup$ The speed-up can be tuned in various ways. If somebody have found a reference along this approach, please edit my answer. $\endgroup$
    – John L.
    Feb 27, 2019 at 21:19
  • $\begingroup$ if $k\ge n/2$, we can run the algorithm with $k$ replaced by $n-k$ and return the elements that are not selected. $\endgroup$
    – John L.
    Feb 27, 2019 at 21:23
  • 1
    $\begingroup$ +1, I also added an answer showing that Knuth basically agrees with you as well. $\endgroup$ Jan 1, 2020 at 9:14
  • $\begingroup$ Why it will run less than 55 times on average? @JohnL. $\endgroup$
    – SteveS
    Jun 17, 2021 at 14:03
2
$\begingroup$

Reservoir sampling solves exactly this problem in $O(n)$ time. When $n$ is really large, but you can calculate the $i$th element just from the index $i$ in $O(1)$ time then you can simulate the sampling by getting the first $k$ elements from a random permutation of $\{1, \dots, n\}$ in $O(k \log n)$ time using Sometimes-Recurse Shuffle.

$\endgroup$
2
$\begingroup$

Is there any survey/review on this topic? Does Knuth cover random n-choose-k sampling in his TAOCP texts?

I first looked at The Art of Computer Programming Volume 4A ("Combinatorial Algorithms Part 1"), specifically Section 7.2.1.3 "Generating all combinations" (which comes under Section 7.2.1 "Generating Basic Combinatorial Patterns", itself under 7.2 "Generating All Possibilities"). This is focused on ways to generate all the n-choose-k combinations (so it turns out not to be in that section, though you can, for instance, look up how to generate the $m$th of the $\binom{n}{k}$ combinations), but the very first paragraph says "and we learned in Section 3.4.2 how to choose [combinations] at random".

This Section 3.4.2 "Random Sampling and Shuffling" is in Chapter 3 "Random Numbers", which opens Volume 2 "Seminumerical Algorithms" of TAOCP. It's a short section, and the very first page has an "Algorithm S" (so outside this section it will be referred to as "Algorithm 3.4.2S"), which is the following (rewriting to the notation here):

  • Set $m$ to $0$ (this denotes the number of items we have selected so far).

  • For $t$ from $0$ to $n-1$,

    • With probability $p = (k - m)/(n - t)$ (which is the number of items that still need to be selected, divided by the number of items still available to consider), select the current item $t$ (and thus increment $m$).

    • If $m$ is now equal to $k$, we are done: terminate.

This may call the random-number generator up to $n$ times, and in fact on average (see Exercise 5 of that section) will call it $(n + 1) k/(k+1)$ times. If we don't know $n$ in advance, we can use the algorithm on the next page (Algorithm 3.4.2R = reservoir sampling), as you mentioned.

But the section also says "Significant improvements to both Algorithms S and R are possible, when $k$/$n$ [in the notation of the question here on StackExchange] is small, if we generate a single random variable to tell us how many records should be skipped instead of deciding whether or not to skip each record. (See exercise 8.)", and if you look at exercise 8, it gives the algorithm due to Jeffrey Vitter (incidentally, a PhD student of Knuth), which is more efficient.

Rewriting in current terminology, the algorithm is basically:

  • Compute the (random) index $X$ of the first item that will be picked.

  • Pick the $X$th item, then repeat for the remaining subproblem (picking $k-1$ items from $n - X - 1$ items).

The details of how to compute $X$ in constant time (if it's worth it) are in that exercise and its solution, which also makes external references to:

  1. "Precise details are worked out carefully in CACM 27 (1984), 703–718," (Vitter, "Faster methods for random sampling")

  2. "and a practical implementation appears in ACM Trans. Math. Software 13 (1987), 58–67." (Vitter, "An efficient algorithm for sequential random sampling")

  3. "A similar approach speeds up the reservoir method; see ACM Trans. Math. Software 11 (1985), 37–57." (Vitter, "Random sampling with a reservoir")


Anyway, regardless of all this, Knuth himself thinks that the best solution is simply rejection sampling (when $k$ is "reasonably large yet small compared to" $n$). Namely:

  • We can assume that $k \le n/2$ (otherwise the problem of sampling $k$ items out of $n$ items reduces to that of selecting the $n - k \le n/2$ items that are not in the sample).

  • Maintain a set $S$ of items we have picked.

  • While we have generated fewer than $k$ items (i.e. size of $S$ is less than $k$),

    • Pick one of the $n$ items uniformly at random.

    • If we had already picked it earlier then do nothing, else add it to the set $S$ of picked items.

That's it. (Note that if we don't ensure $k \le n/2$ we might end up needing to try about $n$ times towards the end; see coupon collector's problem).

The remaining question is how to maintain the set $S$. The answer here by John L. suggests using a boolean array of size $n$, which is great if you can afford it. Else you may still be able to use your language library's standard hash set or equivalent. In Exercise 16 of that section and its solution, Knuth describes a way to do it in (on average) $O(k)$ time (specifically, fewer than $2k$ random numbers generated) using $O(k)$ units of memory (in fact, an array of $2k$ integers), using a data structure he calls an ordered hash table. The solution is described tersely there, but it refers to CACM 29 (1986), 366–367, which is the Programming Pearls column in which Jon Bentley introduced Knuth's "Literate Programming" and has Knuth's 2-page Pascal program for the problem.

$\endgroup$
1
  • $\begingroup$ Forgot to link the last sentence: it (program using ordered hashing) is here. $\endgroup$ Jan 2, 2020 at 0:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.