# Algorithm for efficient weighted sampling from a collection that can efficiently be updated

I'm writing a Monte Carlo simulation in which I have to maintain a large collection of items. This collection contains a great many duplicates, and it will most likely be best to store some or all of these in the form of a number that records the number of duplicates, rather than storing each one individually. (The items are essentially strings, so they have a non-negligible memory cost.)

On each iteration I remove a random item from the collection, and possibly add one or two new items. Each item has a weight, and when sampling from the container it's important that the probability of drawing a given item is equal to its weight multiplied by the number of items of that type.

The items that get added on each iteration may or may not be duplicates of ones that are already in the container. The items will probably be roughly Pareto distributed, with a few items having loads of duplicates and many having none, but it's hard to tell in advance.

This seems like a fairly common thing to want to do, and it also seems not entirely trivial to do it efficiently. The obvious way is to store the items in some kind of hash table, pairing each with an integer representing its frequency. The problem is that then it's not very efficient to sample from the container, since you'd have to iterate over the items in an essentially random order.

On the other hand, if they were stored in something like a balanced binary tree then sampling would be very efficient. However, since the container's contents are changing on every update, this would involve rebalancing the tree all the time, which (at least if done in a naïve way) would be very inefficient.

It seems to me that I need a data structure has a trade-off between speed of sampling and speed of updating. Perhaps some kind of partially-sorted heap-like thing that can be re-sorted with a low average cost when only a few items have changed. But I don't know of such a structure.

If there is a standard algorithm / data structure for this purpose then what's it called? If there is there a readily available implementation in Python and/or C++ then that would be a huge bonus.

(Note that user thus spake a.k. gave an elegant answer to the question in the non-weighted case, but I've realised that for my application, the weights are really necessary.)

You can do this with any balanced binary search tree data structure by additionally maintaining the total weight below each node. To randomly sample, compute a uniform random number between $0$ and the weight of the root node, and traverse down through the tree until you find the leaf whose range contains the random number.

Unfortunately, while maintaining the total weight below a given node is straightforward, it does mean you'll likely have to implement your own balanced binary tree structure. Copying and modifying an existing data structure might be fastest in terms of programmer time; you could even use libstdc++'s red black tree implementation.

Assuming C++, if you're willing to store duplicate items, a time efficient approach for the unweighted case is to store them unsorted in a std::vector or std::deque, v.
You can then efficiently add new items by inserting them at the end of v.
To draw a sample, simply pick a random index i into v and use

sample = v[i];
v[i] = v.back();
v.pop_back();


For the additional storage cost you'll get (amortized in the case of std::vector) constant time iterations.

• This might be quite a sensible way to do it actually. I'm not sure whether my model will generate enough duplicate entries that the extra memory use will become an issue, but it's easy enough to try, and I can optimise later if it does. I'll leave the question open for now, though, in case someone knows of a solution that doesn't take the extra memory. – Nathaniel Apr 14 '14 at 2:25
• One issue is that the items I'm storing have a non-negligible memory cost of their own. Of course I can avoid this by storing pointers in the array, but then I need to have a separate structure (hash table?) so that I can figure out if I've seen that particular item before and get the pointer to the "canonical" instance of it. (Not an un-resolvable issue, just one that's worth mentioning.) – Nathaniel Apr 14 '14 at 7:19
• That would be a sensible way to cut down on storage expense. std::unordered_set would be appropriate, if you're using C++11. – thus spake a.k. Apr 14 '14 at 8:04

The "standard" data structures for this in C++ are the multiset

http://www.cplusplus.com/reference/set/multiset/

based on a balanced binary tree and the unordered_multiset

http://www.cplusplus.com/reference/unordered_set/unordered_multiset/

based on a hash table.

I would try those first and see if their performance is adequate for your requirements.

• How does one randomly sample from these structures? Unless I'm missing something, in either case there seems to be no option but to iterate over all of the items in an essentially arbitrary order, which is $O(n)$. – Nathaniel Apr 13 '14 at 14:30
• Not at all. You can use the equal_range() method for these classes to find just the members in the set with a particular key value. – Bill Greene Apr 13 '14 at 17:42
• Generating a random sample requires more than accessing by key. Note that if a bag contains a monkey, a marmoset, two anteaters, a platypus, an alligator, a vole, a mantaee, a tortoise, a weasel and 90 cats then the probability of drawing a cat should be 0.9. You don't want to iterate through all of the other items on the 90% of sampling runs that will result in a cat, so you really need a way of accessing by frequency, not just by key. – Nathaniel Apr 14 '14 at 2:20

You say that you expect to have a lot of items only once, and a few many times.

While in general a self-balancing search tree like the one Geoffrey suggested would be fastest, there is an easier to implement alternative.

Simply keep the $k$ items that only occur once in a vector, like thus spake a.k. suggested, and the other $m$ unique items in an unordered map. You then keep track of the total number of items in each structure, and when you have to sample the map, you accept the $O(m)$ time. However as you expect $k >> m$, this is likely not a problem.

Suppose you have I unique items, each one with multiplicity m_i. You can create a hash table that maps i -> (item, m_i). When you want to sample, construct a partition of [0,1] (i.e. a sorted array of numbers from 0 to 1) where each interval corresponds to an item, and its width is proportional to the multiplicity E.g. if you have 1 cat, 1 dog, 2 frogs, and 3 rats, the array will be [0, 1/7, 2/7, 4/7, 1]. Then, draw from the uniform distribution, and use binary search to find the interval it falls in; this is the sampled item.

• The problem with this is that I'd have to rebuild those data-structures on every iteration, whenever I add or remove something. I can see a couple of ways to improve on that (keep the array unnormalised, store a separate array of indices into the hash to avoid rebuilding it), but deleting an item would still be $O(n)$. – Nathaniel Apr 23 '14 at 1:40