Computer scientists invent an efficient new way to count

I found the paper took about as long to read as the blog post and is more informative:

https://arxiv.org/pdf/2301.10191

It is about estimating the cardinality of a set of elements derived from a stream. The algorithm is so simple, you can code it and play with it whilst you read the paper.

The authors are explicit about the target audience and purpose for the algorithm: undergraduates and textbooks.

Python implementation:

  def streaming_algorithm(A, epsilon, delta):
      # Initialize parameters
      p = 1
      X = set()
      thresh = math.ceil((12 / epsilon ** 2) * math.log(8 * len(A) / delta))

      # Process the stream
      for ai in A:
          if ai in X:
              X.remove(ai)
          if random.random() < p:
              X.add(ai)
          if len(X) == thresh:
              X = {x for x in X if random.random() >= 0.5}
              p /= 2
              if len(X) == thresh:
                  return '⊥'

      return len(X) / p


  # Example usage
  A = [1, 2, 3, 1, 2, 3]
  epsilon = 0.1
  delta = 0.01

  output = streaming_algorithm(A, epsilon, delta)
  print(output)

Estimating the amount of unique elements in a set and counting the amount of unique elements in a set are very different things. Cool method, bad headline.

This algorithm seems to resemble HyperLogLog (and all its variants), which is also cited in the research paper. Using the same insight of the estimation value of tracking whether we've hit a "run" of heads or tails, but flipping the idea on its head (heh), it leads to the simpler algorithm described, which is about discarding memorized values on the basis of runs of heads/tails.

This also works especially well (that is, efficiently) in the streaming case, allowing you to keep something resembling a "counter" for the distinct elements, albeit with a error rate.

The benefit of HyperLogLog is that it behaves similarly to a hash set in some respects -- you can add items, count distinct them, and, importantly, merge two HLLs together (union), all the while keeping memory fixed to mere kilobytes even for billion-item sets. In distributed data stores, this is the trick behind Elasticsearch/OpenSearch cardinality agg, as well as behind Redis/Redict with its PFADD/PFMERGE/PFCOUNT.

I am not exactly sure how this CVM algorithm compares to HLL, but they got Knuth to review it, and they claim an undergrad can implement it easily, so it must be pretty good!

When do we stop calling this counting and start calling it estimation?

Does finding the number of unique elements in a set actually require comparison of each element with everything else? Can't you use a hashtable? For every element, add it to the table (ignore if already exists), and finally, take a count of keys.

Is it me or is the description of the algo wrong?

    > Round 1. Keep going through Hamlet, adding new words as you go. If you come to a word that’s already on your list, flip a coin again. If it’s tails, delete the word;

    if hash_set.contains(word) {
        if !keep_a_word(round) {
            hash_set.remove(word);
            continue;
        }
    } else {
        hash_set.insert(word.to_string());
    }

    Total word count: 31955 | limit: 1000
    End Round: 20, word count: 737
    Unique word count: 7233
    Estimated unique word count: 772800512

    hash_set.insert(word.to_string());

    if !keep_a_word(round) {
        hash_set.remove(word);
        continue;
    }

    Total word count: 31955 | 1000
    End Round: 3, word count: 905
    Unique word count: 7233
    Estimated unique word count: 7240

CS guys always wanting to throw away a good system and start from scratch.

The trick, he said, is to rely on randomization.

When the space is full, press pause and flip a coin for each word. Heads, and the word stays on the list; tails, and you delete it.

I wasn't expecting to go that far: randomization. How can you verify if the answer is good? Only approximation, maybe..

I don't know a word or phrase for this, but I really enjoy any examples of "thinking outside the box" like this because it's something I struggle with in my professional career. Learning not only the right ways to solve problems, but figuring out the questions to ask that make solving the problems you have easier or even in some cases possible. In this case, it's hey, we don't need exact numbers if we can define a probabilistic range given defined parameters. Other problems are gonna have other questions. I guess my hope is that if I see enough examples I'll be able to eventually internalize the thought process and apply it correctly.

After skimming Knuth's paper — does the algorithm work if values are hashed, that is, the "uniform deviate" is selected deterministically for each unique value of the stream?

Do they assume any particular distribution of the items?

Otherwise Nassim Taleb would like a word with them.

This Nim program might clarify some things for someone.

    import std/[random, math, sets, times]; type F = float
    
    template uniqCEcvm*(it; xfm; k=100, alpha=0.05): untyped =
      ## Return (unbiased estim. cardinality of `it`, 1-alpha CI ebound*(1 +- e)).
      var x1 = initHashSet[type(xfm)](k); var x  = x1.addr
      var x2 = initHashSet[type(xfm)](k); var y  = x2.addr
      var p  = 1.0          # p always 2^-i means could ditch FP
      var n  = 0            #..RNG & arith for bit shifts/masks.
      let t0 = epochTime()
      for e in it:
        inc n   #; let e = xfm # to compute randState.next here
        x[].excl e          # 'd be nice to reduce to 1 hash by
        if rand(1.0) < p:   #..excl w/prob 1 - p but then cannot
          x[].incl e        #..add new keys.  p -> 0 anyway.
        while x[].len == k: # KnuthLoop not CVMIf guards against
          for e in x[]:     #..unlikely case of freeing nothing.
            if rand(1.0) < 0.5: y[].incl e
          p *= 0.5
          swap x, y; y[].clear # ↓ e.bound conservative by ~15X
      (int(x[].len.F/p), sqrt(12.0*ln(8.0*n.F/alpha)/k.F),
       (epochTime() - t0)*1e9/n.F)
    
    when isMainModule:  # Takes nItem/trial dupMask space nTrial
      when defined danger: randomize()
      import std/[strutils, os, strformat]
      let n   =  if paramCount()>0: parseInt(paramStr(1)) else: 100000
      let msk = (if paramCount()>1: parseInt(paramStr(2)) else: 0xFFFF).uint64
      let k   =  if paramCount()>2: parseInt(paramStr(3)) else: 512 # 32KiB
      let m   =  if paramCount()>3: parseInt(paramStr(4)) else: 35
      var ex  = initHashSet[uint64]()
      var xs  = newSeq[uint64](n)
      for j in 1..m:
        xs.setLen 0; ex.clear
        for i in 1..n: xs.add randState.next and msk
        for e in xs: ex.incl e
        let (apx, err, t) = uniqCEcvm(xs, uint64, k)
        let ap = apx.F; let an = ex.len.F
        echo ex.len," ",apx,&" eb: {err:.4f} |x-a|/x: {abs(ap-an)/an:.4f} {t:.1f}ns"

    51160 49152 eb: 0.6235 |x-a|/x: 0.0392 10.5ns
    51300 51840 eb: 0.6235 |x-a|/x: 0.0105 10.7ns
    51250 55168 eb: 0.6235 |x-a|/x: 0.0764 10.9ns
    51388 52352 eb: 0.6235 |x-a|/x: 0.0188 10.5ns

Note also that the error bound in the CVM paper seems very conservative. My guess is that it could maybe be tightened by 15+X (at the scale of those example numbers), but that might also require some kind of special function evaluation. Anyone have any ideas on that?

Also note that since this is just estimating unique identities, you can always just count unique 64-bit hashes of keys (or maybe even 32-bit hashes depending on your accuracy needs and cardinalities) if you care about the key space/key compares are slow.

At first find this paragraph confusing:

Keep going through Hamlet, adding new words as you go. If you come to a word that’s already on your list, flip a coin again. If it’s tails, delete the word; heads, and the word stays on the list.

Why would you delete the word already on the list by flipping coins? Doesn't this reduce the accuracy by counting less words than expected? And will the word be added to the list later?

After thinking about it for a while and reading the paper, I've finally developed a good mental model for how this algorithm works as below, which should convince even a high schooler why the algorithm works:

1. You are given a streaming set of elements from [n] (a finite set of n distinct elements). Now let's give each element a random uniform real number in [0,1]. This helps us choose the elements we want: if we choose all elements with the number below 0.5, we get about half of the elements.

2. Assume for a moment we have unbounded memory. Now we maintain a table of already seen elements, and for each element we keep the last real number attached with that element: so if an element A appears three times with the number 0.3, 0.7, 0.6, we keep A with the number 0.6.

3. Our memory is bounded! So we keep only the elements below a threshold 2^-r, that is, 1, 0.5, 0.25, etc. So at first we will keep all elements, but when we don't have enough memory, we will need to filter existing elements according to the threshold. It's easy to see that when we decrease threshold, only elements already in memory will meet the criteria and continue to exist in memory. No other elements in the stream will be below threshold. Also note that whether an element is in memory depends only on its last occurence. It can exist in memory for a while, get dropped because a later element does not meet the threshold, and get back in.

4. We don't actually have a real number attached to every element! But we can pretend as if there is one. For each new element X from the stream, we replace the number in memory with its attached number, and we only care if its attached number is below 2^-r or not. If it is, it should be in our memory, and if it's not, it should be out. Once the number is in memory, it's a random number in [0,2^-r] and we care no further.

When increasing r, we only care about whether the number is in [0,2^{-r-1}]. This has exactly 1/2 probability. So each number has 1/2 probability of getting out, and 1/2 probability of being kept in memory.

5. Now it's easy to see that whether an element is in the memory depends solely on the real number attached to its last occurence. That is, each distinct element has exactly 2^-r probability of being in memory. 2^r multiplied by the number of elements in memory gives a good estimate of number of distinct elements.

They criticized previous approaches as relying on hash functions. A simplified version of the approach is as follows:

1. Make a hash function that maps distinct elements to different uniform real numbers. So all As compute to 0.3, all Bs compute to 0.7, etc.

2. Use that hash function to transform all elements. The same element maps to the same real number, and different elements map to different real numbers.

3. Keep a set of N smallest distinct real numbers. This works by inserting numbers from the stream to the set. It's easy to see that adding the same element multiple times has exactly no effect; it's either too big to be in the set or the same as an existing number in the set.

4. This gives the Nth smallest real number in the set as K. The number of distinct elements is approximated as N/K.

This algorithm is remarkable, but not in that it's efficient or it's new. Both algorithm using a hash function and algorithms without hash functions has been around and is optimal. I assume this algorithm is the first optimal algorithm without hash functions that is very easy to understand even by undergraduates.

Can some native speaker tell me: is "count" the new "guess"?

I wish we had a theory connecting randomness with time/space complexity.

Intuitively I think the key to many computational limitations is making use of randomness.

Ah, so Thanos was just conducting a census.

"In a recent paper" - really? The paper first appeared in ESA 2022. The revised version (with some errors fixed) is from May 2023. To be fair, compared to the rate of progress in this area between Flajolet & Martin (1985) and the previous SOTA (Blasiok, 2018), a couple of years of sitting on this news is not a lot.

On the topic of counting things, I would like to mention this efficient and easily-implemented algorithm for finding the top-k items in a stream, which I think is perhaps not as well known as it should be:

A Simple Algorithm for Finding Frequent Elements in Streams and Bags

Karp, Shenker & Papadimitriou

https://www.cs.umd.edu/~samir/498/karp.pdf

Computer scientists invent a memory-efficient way to estimate the size of a subset

Am I correct in thinking the accuracy of this method has more to do with the distribution in the sample than the algo itself. Meaning, given a 100% unique list of items, how far off would this estimate be?

I was involved with implementing the DNF volume counting version of this with the authors. You can see my blog post of it here:

https://www.msoos.org/2023/09/pepin-our-probabilistic-approx...

And the code here: https://github.com/meelgroup/pepin

Often, 30% of the time is spent in IO of reading the file, that's how incredibly fast this algorithm is. Crazy stuff.

BTW, Knuth contributed to the algo, Knuths' notes: https://cs.stanford.edu/~knuth/papers/cvm-note.pdf

He actually took time off (a whole month) from TAOCP to do this. Also, he is exactly as crazy good as you'd imagine. Just mind-blowing.

I have questions for Quanta Magazine -- this is probably a shot in the dark, but if anyone can ask them a few questions, here would be three of mine:

1. Have independent researchers evaluated the effectiveness of Quanta Magazine at their mission? [1]

2. In the case of this article, did the editors consider including visualization? If not, why not?

3. Does Quanta have interactive designers and/or data scientists on staff? How many? Why not more?

## Background & Biases

I'm trying to overcome my disappointment in Quanta by being curious about their mission, organization, and constraints. I am glad they exist, but sometimes your "closest allies" can be one's harshest critics.

I'm greatly troubled by the level of scientific, mathematical, and rational understanding among the denizens of the world. But for some reason, the state of science writing bothers me more. It would seem that I hold out hope that science writers could do better. This may be unfair, I admit.

Anyhow, rather than just bash Quanta for being, say, not as good at the best math blogs or YouTube channels (such as 3Blue1Brown), I really want to figure out (a) if I'm missing something; or (2) if they are actively trying to improve; and (3) what we can all learn from their experience.

[1] From https://www.quantamagazine.org/about/ : "Quanta Magazine is an editorially independent online publication launched by the Simons Foundation in 2012 to enhance public understanding of science. Why Quanta? Albert Einstein called photons “quanta of light.” Our goal is to “illuminate science.”"

From the paper [0]:

We state the following well-known concentration bound, Chernoff bound, for completeness.

Which variant of the Chernoff bound is this? This is almost the (looser variant of the) multiplicative form, but it's not quite right (per the use of 1+delta instead of a single parameter beta). In particular, that bound is only guaranteed to hold for delta >= 0 (not beta = 1 + delta > 0 as asserted in the paper)

[0] https://arxiv.org/pdf/2301.10191

[1] https://en.wikipedia.org/wiki/Chernoff_bound#Multiplicative_...

edit: to be clear: I'm not sure at all whether this breaks the proof of correctness, although I'm having a bit of difficulty following the actual details (I think I'd need to work through the intermediate steps on paper).

"Computer scientists invent an efficient new way to approximate counting large unique entries" - fixed the title

I wonder how the error changes if you start a new round just before finishing the stream, compared to if the last word in the stream just fills the buffer and ends a round.

HyperLogLog uses additions, it keeps sums. Thus, you can subtract one HLL sums from other. This is useful if stream supports deletion. Streams with deletions can be found in log-structured merge trees, for one example, so one can estimate count of distinct elements in all of the LSM tree hierarchy.

The algorithm in the paper does not allow for deletions.

Also, if one counts statistics of the stream of large elements (say, SHA-512 hashes, 64 bytes per hash), this algorithm requires some storage for elements from this stream, so memory requirement is O(table size * element size).

Imagine that you’re sent to a pristine rainforest to carry out a wildlife census.

alt-f4

Huh, that's a clever twist on reservoir sampling. Neat.

What if you’re Facebook, and you want to count the number of distinct users who log in each day, even if some of them log in from multiple devices and at multiple times?

Seems like a bad example of when this algorithm could be useful in practice.

If you already know you will want this info when designing the login process, it's simple: keep track of the date of last login for each account, and increment the unique user counter when the stored value is different from the current one.

And even if not, it should still be possible to "replay" a stream of login events from the database later to do this analysis. Unless maybe you already had years worth of data?

While the algorithm is interesting and I commend the authors for the algorithm, it should be noted that this is a guesstimate only. So, it's not really an efficient way to count distinct values but an efficient way to get an approximate count of distinct values.

It is quite amazing that after we had so many talented researchers working for decades on this problem, the authors were able to come up with such an elegant and simple solution. Another proof that research is very necessary and useful.

Feels like reservoir sampling to me

There is a big practicality problem I see with this algorithm. The thresh defined in the paper relies on the length of the stream. It seems to me that in a scenario where you have a big enough data set to desire a solution that doesn't just store every unique value you don't know the length. I did not make it through all the proofs but I believe they use the fact that the defined threshold has the length in it to prove error bounds. If I were to use this in a scenario where I need to know error bounds i would probably ballpark the length of my stream to estimate error bounds and then use the algorithm with a ballpark threshold depending on my systems memory.

Another practical thing is the "exception" if nothing is removed on line 6 in the original algorithm. This also seems needed for the proof but you would not want in production, though the chance of hitting it should be vanishingly small so maybe worth the gamble?

Here is my faithful interpretation of the algorithm. And then a re-interpretation with some "practical" improvements that almost certainly make the provability of the correctness impossible.

    func CountUnique(scanner *bufio.Scanner, epsilon float64, delta float64, m int) int {

    X := make(map[string]bool)
    p := 1.0
    thresh := int(math.Ceil((12 / (epsilon * epsilon)) \* math.Log(8*float64(m)/delta)))

    for scanner.Scan() {
        a := scanner.Text()
        delete(X, a)
        if rand.Float64() < p {
            X[a] = true
        }

        if len(X) == thresh {
            for key := range X {
                if rand.Float64() < 0.5 {
                    delete(X, key)
                }
            }
            p /= 2
            if len(X) == thresh {
                panic("Error")
            }
        }
    }

    return int(float64(len(X)) / p)

  func CountUnique2(scanner *bufio.Scanner, thresh int) int {

     //threshold passed in, based on system memory / estimates
    X := make(map[string]bool)
    p := 1.0

    for scanner.Scan() {
        a := scanner.Text()
        delete(X, a)
        if rand.Float64() < p {
            X[a] = true
        }

        if len(X) >= thresh {  // >= instead of == and remove the panic below
            for key := range X {
                if rand.Float64() < 0.5 {
                    delete(X, key)
                }
            }
            p /= 2
        }
    }

    return int(float64(len(X)) / p)

I tested it with Shakespeare's work. The actual unique word count is 71,595. With the second algorithm it is interesting to play with the threshold. Here are some examples.

threshold 1000 Mean Absolute Error: 2150.44 Root Mean Squared Error: 2758.33 Standard Deviation: 2732.61

threshold 2000 Mean Absolute Error: 1723.72 Root Mean Squared Error: 2212.74 Standard Deviation: 2199.39

threshold 10000 Mean Absolute Error: 442.76 Root Mean Squared Error: 556.74 Standard Deviation: 555.53

threshold 50000 Mean Absolute Error: 217.28 Root Mean Squared Error: 267.39 Standard Deviation: 262.84

This actually comes at a good time. I'm currently refactoring a system that counts visits using inserts into a table, with a tuple of date and IP. I was planning to replace it with a HLL approach, but this is really interesting.

Given the topic of the paper[0], the footnote is especially charming:

the authors decided to forgo the old convention of alphabetical ordering of authors in favor of a randomized ordering, denoted by r⃝. The publicly verifiable record of the randomization is available at https://www.aeaweb.org/journals/policies/random-author-order...

[0]: https://arxiv.org/pdf/2301.10191

edit: formatting

The counting/estimstion technique is rather like a floating point number. An integer exponent k and a mantissa of a population.

very interesting solution. A perfect example of how even some very complicated problems can sometimes have simple solutions. Will definitely try to write an implementation of this.

The only minor downside to this is that it's obviously not deterministic since it depends on chance. But for many applications where the dataset is so big it doesn't fit in memory, that's probably just a tiny rounding error anyway.

With other count-distinct algorithms, you can do unions and intersections on the "sets". (theta sketches, even bloom filiters)

Can you do this with CVM?

The algorithm uses less memory, but more CPU time because of rather frequent deletions, so it's a tradeoff, not just generally better algorithm, as article may suggest.

Wondering if this approach could I found myself be applied to CFD (computational flow dynamics) methods to reduce the total volume of data points while still getting approximately close to the correct final answer?

Does anyone else notice that quanta consistently has good illustrations and diagrams for their articles? Good visualizations and graphics can do extraordinary things for the understandability of a topic -- people like 3Blue1Brown understand this, and I'm glad quanta puts effort into it.

Knuth talks about this paper here - "The CVM Algorithm for Estimating Distinct Elements in Streams" - https://cs.stanford.edu/~knuth/papers/cvm-note.pdf

New interview question: Tell me how you would estimate the number of unique words in Hamlet using only 100 elements.

nice read!