There is a common pattern you see pop up in codebases everywhere that many of us have written: given a set of objects, pick one at random. Every object should have the same chance of being picked. A rather simple and direct solution is to pick a big random number, then clamp that number to the number of choices by modulo:

Set<T> choices = ten_things();
u64 r = random_u64(); // uniform chance over [0..UINT64_MAX]
T chosen = choices[r % 10];

I fondly remember choosing this option several times when I was much younger, especially when writing C — its standard library does not offer anything beyond rand() — and also my other language de jour, Haskell. After all, it is a fairly immediate and intuitive solution when you have such a problem, and almost any codebase or language, no matter how austere or feeble, will have some mechanism to pick a random number.

But there is a hidden problem with this solution: it doesn’t preserve the underlying uniform distribution of the given random_u64() function. That means that the variability of the chosen object you pick doesn’t respect what we might intuitively think: we might think every item has a 1/10 chance of being chosen. It does not.

(Non-)Preservation of uniform choice

Let’s start with a simpler problem: pick an integer from the interval [0, 9]. There are 10 possible picks. Assuming that our random choice function has a uniform distribution, every integer will have a chance of 1/10 = 10% of being chosen.

Then, let us use that integer to randomly pick a value from a set of 3 objects, using the modulo operator like above. The intuitive (implicit) hope is that this operation picks one of the 3 objects, each with a 1/3 = 33.3% chance; in other words, we we want the resulting probability distribution to “respect” the uniform distribution of the input variable.

Unfortunately, this is not the case; simply tabulating the outcomes makes that obvious:

Input value  -> Choice
-----------------------------
{0, 3, 6, 9} -> choose object #1
{1, 4, 7}    -> choose object #2
{2, 5, 8}    -> choose object #3

In other words: of the given 10 inputs, 4 of them yield choice 1, while choice 2 and 3 only have 3 inputs that do so. This means choice 1 is picked 40% of the time rather than the intended 33%, and likewise 2 and 3 are only picked 30% of the time.

It’s one of the classic blunders, especially when random_u64() is all you have. Now, the correct behavior instead is a random_between(l, h) function which gives each inclusive number between l and h a 1/((h-l)+1) chance of being picked. Such a function is, as we might now expect, a bit more complex!

The issue in my mind is twofold. One, I think it is suboptimal API design to only offer a uniform random function, especially given that one of its most useful variants is easy to implement incorrectly. But it is also a good example of mentally working in the wrong domain in the first place.

Giving an explicit distribution is (hopefully) easier

I think the first thing is the fact that an API like random_u64() is a bit too low level and specialized. We don’t really want to sample a uniform number but a more general operator to express “how often does something happen”.

On one hand, you rarely need to only pick a random number, you normally want to pick some object out of some possible choices; second, it often is useful to weigh certain objects as more or less likely than other objects — the uniform distribution is not the only valuable one! random_u64() does not help with either of these things; a random_between(0, 9) only helps with the first at a glance (more on that shortly).

Most of the time, I think it is useful to instead reach for a more general operator, that makes you confront these issues directly: a random_choice() function that takes a set of discrete options and their probability of being chosen. In other words, you should spell out the probability distribution yourself. With this, you can represent both of the previous behaviors, and many more interesting ones.

Typically when we talk about probabilities we think of them as summing to 1.0, like the following reformulation of the original example; in such a design the erroneous “lopsided” pick is going to stick out like a sore thumb:

T chosen = random_choice([
  (First,  0.4), // 40% chance?!?!?!
  (Second, 0.3), // 30% chance
  (Third,  0.3), // 30% chance
])

I think this is a case (one of many) where a generalized API, a more abstract API, can in fact make much of your code easier to reason about and write correctly too. You should of course still have random_u64() and random_between() as well — but I think random_choice() is the much better default tool because it will help you fall into the pit of success.

Improving the API: relative integer weights

The previous interface is pedagogically nice but floating point instability can muddle with our ability to pick the precise weights we might need, as they must sum exactly to 1.0. Also, the requirement to sum to 1.0 is annoying as you might have to calculate the appropriate scales of each probable choice. And so when confronted with any problem involving floating point, we will use the tried and true method to handle it: get rid of it and do something else.

You will find that most standard libraries, like Python, do exactly this, with relative integer weights. The “relative” part means the numbers don’t sum to 100; instead you simply sum them and the weight of any choice is the ratio of that choice’s chance to the sum. The previous weights would be represented with (4, 3, 3) instead. 4 + 3 + 3 = 10 and so 3 / 10 = 0.3 as we expect.

One of the nice things about relative integer weights is they behave more intuitively, so you’re more likely to use it correctly. If you count the number of things in a box and see 10 blue things, and 5 red things, you can simply write out your choices at [(Blue, 10), (Red, 5)] and the resulting choice you will be distributed as you expect.

As another small note: implementing this API is much easier to do intuitively if you have random_between(l, h), not just random_u64(), which I think is another sign that random_u64() is too low level. Let’s say we have the integer weights 15 + 12 + 3 = 30. Then simply draw a random number choice = random_between(0, 29). If the chosen value is between [0, 14] it is the first option, between 15, 26 is the second option, and [27, 29] the third. In other words, you can simply map your integer weights onto to subsets of the interval [l, h] and then choosing uniformly from that interval. This is so useful and simple that it should, in my opinion, also be included with the other functions.

Working in the wrong mental domain

The other more subtle problem to me is recognizing what domain you’re operating in. In the previous examples, the things we want to talk about are not numbers, but random variables. We normally call these random variables by capital letters, X, Y, Z, etc.

Much like ordinary numbers, random variables do have an algebra to them; you can add random variables and subtract them and they feature exponents. But these operators also impact the underlying distributions and their expected outcomes, like they did here. But the most important thing is that non-linear operators do not respect the expected value of random variables. In particular, a function f over the expected value E[f(X)] is not always the same as f(E[X]).

That is what happened here; if we aren’t careful it’s easy to think we’re trying to retain some property of the value x = random_u64(), but we are actually trying to retain the properties of the underlying random_u64() function.

Our case

We are customers of Antithesis. In their platform, we effectively find bugs by exercising important code paths at random — it’s a big deterministic fuzzer, for a whole operating system. And like a fuzzer, it is useful to think about code coverage as a guiding light: “did this code get tested” is equivalent to asking “did the fuzzer find this path”, and if the fuzzer did not find that path, that code was not tested.

Increasing coverage is a matter of expressing positive and negative code paths, and then finding them under test. A negative path is a bad case that, if found, is a bug, such as an assert!() getting tripped — for instance, we might read some data from a cache, keyed by the hash of the data itself, and then assert!() that the key matches a a recomputed hash of the data.

A positive path might be ensuring that an upload function which does something different for “small” and “large” uploads actually tests both cases. In such a case, expressing your distribution directly is pretty convenient:

enum BlobSize {
  Small,  // 1KiB: fast path
  Medium, // 1MiB: streaming write path
  Large,  // 2MiB: more
  XLarge, // 8MiB: more more more
}

let choices = vec![
  (Small,  88),
  (Medium, 04),
  (Large,  04),
  (XLarge, 04),
];

if let Some(size) = weighted_choice(choices) {
   // ... upload blob of chosen size ...
};

In this case, under test, ~7/8 uploads are small blobs of data, and ~1/8 uploads are large blobs (of three various sizes). Being able to play with these numbers is extremely valuable to help explore your state space! Extending this example to express e.g. a Poisson distribution is quite easy too.

Another useful case here I’ve found is extending the set of choices under some set of conditions; for example by using choice.push() to add something to the vector above, when some “good” or “bad” thing happens, to steer the search — this would require re-balancing the probabilities in the above example, but you can imagine structuring this example differently to avoid that.

Which leads us to this post, and something that happened last week: upon review of our Antithesis-related code recently, we ultimately were able to remove all our uses of random_u64() that had creeped in with (weighted) random_choice(), as we only ever wanted discrete weighted choices in practice anyway; we’ve put a moratorium on introducing new calls to it until further notice. I figured this lesson was important enough to not leave buried in a commit message somewhere since I keep forgetting it.

The next time you see some code like the original sample, ask yourself: can I just represent the desired distribution directly? If so you will hopefully find a much richer tool to play with!