Probabilistic Counting with Stochastic Averaging (PCSA)
WIP:
- Need to provide an interactive visualisation somehow?
- Do the analysis of the algorithm.
Introduction
The Flajolet–Martin algorithm (PCSA) estimates how many distinct elements (cardinality) appear in a stream using very little memory by exploiting randomness in hash values and counting trailing, as in 0, 1, 2, … 64 … etc. zero bits in binary hash outputs.
In computer science, the count-distinct problem (also known in applied mathematics as the cardinality estimation problem) is the problem of finding the number of distinct elements in a data stream with repeated elements. This is a well-known problem with numerous applications. The elements might represent IP addresses of packets passing through a router, unique visitors to a web site, elements in a large database, motifs in a DNA sequence, or elements of RFID/sensor networks.
Uses or Use Cases
- Counting unique users/IPs/cookies on high-traffic sites. Used when exact counting is too expensive (billions of events/day).
- For example: “How many unique visitors hit this endpoint today?”. PCSA avoids storing every user ID and instead uses compact bitmaps.
- Estimating: Unique requests, unique trace IDs and unique error signatures.
- Counting: Unique source IPs hitting a server and unique destination ports scanned.
- DDoS detection.
- Port scan detection.
- Traffic cardinality estimation.
- Real-time pipelines: “How many unique users in the last 5 minutes?”.
- Ad Tech & Marketing: Counting unique impressions and counting unique users exposed to an ad
Core Idea
- Hash each element to a (pseudo) random-looking binary string.
- Look at how many trailing zeros are in that hash (e.g., has three trailing zeros).
- Very long runs of trailing zeros are rare; seeing them suggests many distinct elements have been observed.
A rule of thumb: if you see a hash with trailing zeros, that event happens with probability about , so you expect to see it only after on the order of elements.
Basic Single-register Algorithm
Analysis
References
One of the most fascinating algorithms if not the most, that I have come cross. It is a solution to the count-distinct problem.
Given a Multiset $\mathfrak{M}$ of random binary strings each of size $L$. Let $R$, the rank, be the maximum index of the first 1-bit amongst all the binary strings. If $n$ is the number of distinct elements, $2^R \approx n$.
An implementation in Go can be found in https://github.com/banaio/countdistinct/tree/master/pcsa, see the whole repo for other implementations https://github.com/banaio/countdistinct.
tl;dr;
A Set-like data structure with a bounded error in the the answer. See the [Implementation in Go](#Implementation in Go) section for an implementation and runnable code.
Introduction
This family of algorithms—1, 2 and 3—that are based bit-pattern observables are pure magic:
We have seen in the previous section that the result $R$ Of the $COUNT$ procedure has an average close to $\log_2 \varphi n$, with a standard deviation close to $1.12$. Actually the values of
are amazingly close to $n$ as the following instances show:
Meaning the bits needed is $\log_2\log_2 n$; for $n = 2^{32}$ it’s no smaller than 5 since $R$ needs to store an index between 0–31.
I am intentionally using a small sample size for the stream and the size of the hash function as I think it’s easier to see exactly what’s going on, it’s also much easier to show what the functions in paper evaluate to. If you spot an error, drop me an email as I am not aware of how to get a comments-like section in VuePress. The sketching posts by Neustar, e.g., Sketch of the Day: Probabilistic Counting with Stochastic Averaging (PCSA), deserve a mention as well.
Problem and Constraints
The algorithm addresses the count-distinct problem. Let’s recap the points and try to imagine how we’d implement it:
- Items are being pushed to us.
- There’s a lot of them.
- We need to keep track of the items seen, without allowing them to be removed.
- Upon query, reply with an approximate answer of the unique number of items.
- The device that has a limited amount of memory, say, ~1,024 Kilobytes.
So,
- Seems like we’ll be reading from a Stream
S, or $\mathfrak{M}$ as in the papers, in one-pass and data flows in our direction only. Sis large. We’ll see why this is important later on.- We’re designing a Set-like data structure, say,
PCSA, that allows you toPCSA.Add(item). Removals however are not permitted. - The true unique item count, the value returned by
PCSA.count(), should only vary so much from the correct unique count. - We don’t have a lot of memory at hand.
Naive solution
You could use a Set, but there’s really no fun in that. The amount of memory required to store all the possible hash values you could encounter is 512MB memory, ignoring any optimizations such as compressed bitmaps or sorting; $32 * 2^{32}$ each unique hash requires $32$ bits and we can see $2^{32}$ such hash values.
Binary Strings
You could use a Set, but there’s really no fun in that. Instead we’ll see how this family of algorithms solve this problem. Conceptually they’re all somewhat similar, and the main ideas are to:
- Store the maximum consecutive zeros, Least significant bit (LSB), from the binary representation of all the hash values we’ve seen, call it, $R$. That is, hash something, convert hash to a binary string, and count how many successive zeros we see from the beginning till the end, this is the position of the LSB 1-bit.\ If we’re using a 32-bit hash function, the maximum value $R$ can take is $32$ ($\log_2 2^{32}$), and to count to $32$ we need $5$ ($\log_2 32$) bits of memory, i.e., $\log_2 \log_2 2^{32}$. You may have noticed, and briliant just like the insight, the name of the algorithm, LogLog, is derived from the number of bits needed to store the counter $R$—$\log_2 \log_2 2^{32}$.
- Estimate the unique item count as $2^R$. Say, what? If you already have an understanding of why this might work, you might want to look at the section on how the authors reduce the error of the estimate.
Since PCSA observes patterns in the binary strings coming from S, we’ll look at what patterns these binary strings can take. Limiting the hash to 4-bits long, $L = 4; 2^L = 16$, just so that it’s easier to visualise. I might extend it to support dynamic zooming. For now:
- A row for each binary string; $2^L - 1$ in total.
- A column for each bit in the binary string; $L - 1$ in total.
- Darker cells are 0-bits and lighter cells are 1-bits.
- Rows are also clickable if you wish to avoid a random selection.
- The index of the LSB, $\rho(y)$ in the paper, and $R$ will be outlined.
You can randomly select a row via the button.
Analysis
The sections that follow explain how the authors, pg. 186 onwards, define the random variable $R_n$. The mechnical definitions are listed first as $R_n$ is Discrete random variable and the rest of the theorems after. As we go along I’ll try to describe what each part of the formula is doing, I find it easier to reason about the workings when I do this.
Expected value
The Expected value of $R_n$. The weighted sum of the values $R$ can take, $0,1,…31$, against the chance of it occurring:
Variance
The Variance (the 2nd version of the formula is used) and hence the Standard deviation. Like above, though, we square the possible values ($k$) of $R$ against the chance of it occurring and then deduct the squared expectation ($\bar{R}_n$).
Distribution - Theorem 1
The plan is to derive $p_{n,k}$. I cannot see it in the paper but I suspect it’s derived as $p_{n,k} = q_{n,k} - q_{n,k+1}$; the chance that it’s greater than or equal to $k$ then removing the chance that it’s greater than $k+1$, or if you like, get a subset then from this subset get another subset.
The authors provide a proof that given $n$ elements the random variable $R$ will take on a value $k$ or larger as $q_{n,k}$, $v(j)$ is the number of 1’s in the binary form of $j$—the Hamming weight or population count:
Possibly not relevant at this point, but let’s break down that down to make it less intimidating:
- $\sum_{j=0}^{2^k}$; Loop through all the possible binary patterns and summing their chances. The events section below explains why we sum the probabilities instead of multiplying.
- $(-1)^{v(j)}$; Thue–Morse sequence mentioned in the Acknowledgements section, pg. 209, of paper. Produces sequences of $-1$ and $1$. If $v(j) > 0$ this can be rewritten to $\href{http://www.wolframalpha.com/input/?i=e^%28i+n+%CF%80%29}{e^{i n \pi}}$.
- {: .math-red} $\left(\color{red}{1-}\frac{j}{2^k}\right)^n$; The chance of not seeing $j$, the current hash if you like, in the $n$ items. The chance of not seeing $j$ will not effect the chance of us seeing it in the future, i.e., it’s independent so we multiply.
I think, double-check this, the roundabout route of first defining the events to then derive the distribution is required as $R_{n}$ is the maximum of $n$ Geometric random variables—meaning it’s not easy to derive—and the above is sort of calculating the max.
Events
The authors then define an Event for each possible position of LSB 1-bit, and then evaluate the probabilities to attain $q_{n,k}$. You can think of it as a way of segmenting the Stream S that contains $n$ items into the sets below depending on the LSB 1-bit of the hashed value of the item:
- $\class{set-e0}{\boldsymbol E_0} = \color{red}{\texttt{1}}\texttt{….}$: All items having an LSB in the first bit of their hashed value.
- $\class{set-e1}{\boldsymbol E_1} = \texttt{0}\color{red}{\texttt{1}}\texttt{…}$: All items having an LSB in the second bit of their hashed value.
- $\class{set-e2}{\boldsymbol E_2} = \texttt{00}\color{red}{\texttt{1}}\texttt{..}$: All items having an LSB in the third bit of their hashed value.
- $\class{set-e3}{\boldsymbol E_3} = \texttt{000}\color{red}{\texttt{1}}\texttt{.}$: All items having an LSB in the fourth bit of their hashed value.
- $\class{set-k4}{\boldsymbol K_4} = \texttt{0000}\color{red}{\texttt{1}}$: All items having an LSB in position $k$, 4 in this example, or higher. This allows us to capture the pattern $\texttt{0000}$.
Now we can represent all, a re-look at the patterns visualisation above might help, the possible hash values we can see from the Stream using the above disjoint subsets. Specifically, the hash will be contained in any one of the above sets, so adding all the sets together produces a new set with all the possible hash values, the Sample space, e.g., when the number of bits the hash produces is $k=4$, a single draw will come from: $E_0 + E_1 + E_2 + E_3 + K_4$. The addition of each of the subsets is due to the events being disjoint. Repeating this $n$ times as each draw of a hash value is independent, we get a polynomial:
Looking good, yo, we’re making progress. We’ve just defined the entire sample space of the Stream S, in terms of subsets.
Extracting the events
Inclusion–exclusion principle.
Assigning probabilities to events
To get the distribution $q_{n,k}$, we need to assign probabilities to all the subsets of events we just extracted above, see Probability measure.
Asymptotic limits of the distributions - Theorem 2
If you’re not entirely sure what the authors are describing in theorem 2, the explanation in chap. 4 pg. 143 of Mining of Massive Datasets is great. There are three cases to consider:
- $R_n \ll \log_2(n)$: The estimate is too low.
$R_n = \log_2(n)$: The estimate is right.- $R_n \gg \log_2(n)$: The estimate is too high.
Reducing the error
I’ll explain this in another post as enough has been covered in this post.
Summary
…
Implementation in Go
Links
- Flajolet–Martin algorithm
- HyperLogLog
- Analytic Combinatorics:
- Analytic Combinatorics by Philippe Flajolet and Robert Sedgewick, see
APPENDIX C Concepts of Probability Theory. - Robert Sedgewick’s Analytic Combinatorics book site
- Philippe Flajolet’s Analytic Combinatorics book site
- Analytic Combinatorics by Philippe Flajolet and Robert Sedgewick, see
- Philippe Flajolet’s lectures
- Geometric series
- Geometric distribution