ANN Indexes: How to Not Look at a Billion Vectors

Chapter 1: The Recall–Latency–Memory Triangle

Chapter 0 established that you must skip most of the data. This chapter builds the vocabulary for describing how much you skipped and what it cost you — because "it's approximate" is not a specification, and a vendor benchmark that reports only latency is not a benchmark.

Recall, defined precisely enough to argue about

Fix a query q and a corpus. Let G_k(q) be the ground truth: the set of the k genuinely closest vectors, computed once, offline, by brute force. Let R_k(q) be what your index actually returned. Then

Three traps live in that innocent formula, and every one of them has been used, deliberately or not, to make an index look better than it is.

Trap one: k on the left is not always k on the right. The literature writes x-recall@y: did the true top-x appear among the y results returned? "10-recall@10" is the honest, strict version. "1-recall@10" asks only whether the single true nearest neighbour is somewhere in your ten — a far easier bar. On the same index at the same settings these two numbers can differ by fifteen points. Always read which one is plotted.

Trap two: ties. If several corpus vectors sit at exactly the same distance, the ground truth is ambiguous and recall becomes unstable. On real float embeddings exact ties are vanishingly rare; on quantized or binary data they are common enough to matter.

Trap three: the ground truth must use the same metric you deploy. Computing ground truth with L2 and then serving with cosine on un-normalised vectors produces a recall number that measures nothing.

Worked example 3 — recall compounds, and it compounds badly

Here is the calculation that changes how people set their recall target. A retrieval-augmented generation pipeline works like this: retrieve 20 chunks with the vector index, rerank them with a cross-encoder, feed the top 5 to the model.

Suppose answering a particular question genuinely requires three specific chunks — a definition, an example, and a caveat, scattered across the corpus. Your index has recall@20 = 0.90. Treat the three retrievals as roughly independent (they are not exactly, but the direction of the error is what matters):

Twenty-seven percent of multi-hop questions are unanswerable before the model has read a single token, and no amount of prompt engineering will fix it, because the evidence is not in the context. Now push recall to 0.98:

Eight points of recall bought a 4.6× reduction in end-to-end failures. That is why the interesting part of every recall/QPS curve is the far right end, above 0.95, where the curve is steepest and most expensive — and why benchmarks that stop plotting at 0.9 are hiding the region you will actually operate in.

What recall does not capture

Recall is a set metric: it asks which items came back, not in what order. Two failure modes slip straight through it, and both are visible to users.

Ordering inside the returned set. An index that returns the true top 10 in the order 7, 3, 1, 9, 2, 5, 10, 4, 8, 6 scores recall@10 = 1.000 — a perfect score for a result page whose best item is third from last. If your interface shows a ranked list and users read from the top, that is a real regression that recall cannot see.

How badly you missed. Recall counts a miss identically whether the item you returned instead was almost as good or completely unrelated. Two indexes at recall 0.9 can differ enormously in whether the substituted item was a near-duplicate of the right answer or noise.

The other two axes

Latency and throughput are different numbers and conflating them is the most common benchmarking error in this field.

Memory is measured per vector, including everything. The number that matters is total resident bytes divided by N, and it must include the index overhead, not just the payload. A worked comparison at N = 1,000,000 and d = 128 (the classic SIFT1M setup, so these numbers are checkable against published results):

Read the last two rows against the first. Product quantization takes the payload from 512 bytes to 8 — a factor of 64 — and at that point the identifier costs as much as the vector. That is not a rounding detail: FAISS's 64-bit ids are half the memory of an IVF*,PQ8 index, which is why large deployments switch to 32-bit ids or to implicit ids derived from position.

Notice also that HNSW's overhead is not small. At M = 48 the graph costs 384 bytes per vector — three quarters of the vector itself. HNSW is fast partly because it keeps everything in RAM at full precision, and that is exactly why it cannot reach a billion points on one machine.

Why QPS is always plotted on a log axis

A detail of presentation that encodes a fact about the problem. Recall/QPS plots put throughput on a logarithmic axis, and the reason is that throughput varies by three or four orders of magnitude across the parameter sweep of a single index while recall varies over a factor of well under two.

From ef = 10 to ef = 2048, recall moves 0.199 and throughput falls by a factor of 103. On a linear y-axis the entire high-recall region — the only part anyone deploys at — is squashed flat against the bottom and invisible. The log axis is not a stylistic preference; it is the only way to see the region where decisions are made.

The same asymmetry explains a common misreading of vendor claims. "Our index is 3× faster" is a statement about a vertical distance on a log plot, which corresponds to a very small horizontal distance in recall — often smaller than the difference between two reasonable recall definitions. Always ask what recall the 3× was measured at, and whether the comparison held it fixed.

Why it is a triangle and not three dials

The three quantities are not independent; each index family fixes one relationship and lets you trade the other two.

Recall per byte, and recall per millisecond

A useful discipline once you have a measurement harness: stop comparing configurations by recall and start comparing them by what the last point of recall cost. The interesting quantity is a derivative, not a level.

Take a plausible sweep on a 10-million-vector corpus at d = 768 and tabulate the marginal cost of each step up:

The last row is the one that ends arguments. Moving from the reranked PQ configuration to full HNSW buys three tenths of a point of recall for thirty-one gigabytes — roughly a hundred gigabytes of RAM per point of recall. Whether that is a bargain or an absurdity depends on the pipeline, and Worked Example 3's compounding argument is how you decide. But the question is now numeric.

Formalising approximate search — and why the formalism is not used

The theory community defines the problem cleanly. A (1 + ε)-approximate nearest neighbour query returns a point p such that

This is a lovely definition. Locality-sensitive hashing (Indyk & Motwani, 1998) comes with provable bounds in exactly this framework, and for two decades it was the theoretically respectable answer.

Practitioners essentially never use it, and the reason is Chapter 0's concentration result. Work out what ε buys you at d = 768. Suppose the true nearest neighbour is at distance 0.87 (cosine 0.62) and the great mass of random points sits at distance about 1.41 (cosine 0). Set a seemingly tight ε = 0.05:

Any point within 0.914 is an acceptable answer. Now: how many of the billion points lie between 0.87 and 0.914? Because distances concentrate near 1.41 with a standard deviation around 1/√768 ≈ 0.036, the region [0.87, 0.914] is roughly fifteen standard deviations below the mean — and contains, in practice, either the true neighbour and a handful of its genuine siblings, or nothing at all. So ε = 0.05 is either trivially satisfied or unsatisfiable, and it tells you nothing about whether the system is any good. Push to ε = 0.5 and suddenly the acceptable region is [0.87, 1.305], which is most of the corpus: a guarantee that admits garbage.

The index card — four numbers per candidate

When you finish evaluating a configuration, write it down in the same shape every time. Four numbers and one sentence make configurations comparable months later, when nobody remembers what "the fast one" meant.

Two disciplines make the card honest. First, always record the knob — a recall number without its parameter setting is not a measurement, it is a claim. Second, split hot and cold memory: an index using 96 resident bytes plus a 3 KB cold rerank path has a completely different cost structure from one using 3,168 resident bytes, and a single "bytes per vector" column hides exactly that.

Sim — the triangle, explorable

Before the machinery, get a feel for the shape of the trade. The instrument below has two knobs: how many bytes each vector is compressed to, and how many candidate vectors the search is allowed to visit. Everything else — recall, latency, memory, and which named index family lives at your current setting — falls out of those two.

Three things to notice. First, the recall curve saturates: past a few thousand visited vectors, more effort buys almost nothing, because you are already touching everything that could plausibly win. Second, dropping from 32 bytes to 8 does not gently degrade recall — it lowers a hard ceiling that effort cannot climb back over. Third, the region every production system actually lives in is the knee, and the knee moves with your code size.

Worked example 4 — tail latency, and what sharding does to it

Sharding is the standard answer to "the corpus does not fit." It has a cost that no single-machine benchmark can show you, and the arithmetic is short enough to do here.

Suppose one shard has p50 = 4 ms and p99 = 20 ms — a perfectly normal graph index under load, where the variance comes from differing hop counts and cache behaviour. Now fan a query out to S shards and merge the results. The request finishes when the slowest shard finishes, so the probability that the whole request comes in under 20 ms is

At 100 shards, the 20 ms figure that used to be your rare bad case happens on nearly two requests in three. To recover a real p99 for the whole request you now need each shard's p99.99 to be acceptable — a far harder engineering target than its p99.

The measurement harness, which you should build before any index

Everything in this chapter is unfalsifiable until you can run it. The harness is genuinely short:

python — recall + latency sweep, the whole thingimport numpy as np, time

def ground_truth(X, Q, k=10, block=200_000):
    """Exact top-k by brute force, in blocks so it fits in memory."""
    best_s = np.full((len(Q), k), -np.inf, "float32")
    best_i = np.zeros((len(Q), k), "int64")
    for b in range(0, len(X), block):
        S = Q @ X[b:b+block].T                       # (nq, block)
        S = np.concatenate([best_s, S], axis=1)
        I = np.concatenate([best_i, np.arange(b, min(b+block, len(X)))
                                     [None, :].repeat(len(Q), 0)], axis=1)
        top = np.argsort(-S, axis=1)[:, :k]
        best_s = np.take_along_axis(S, top, 1)
        best_i = np.take_along_axis(I, top, 1)
    return best_i

def recall_at_k(pred, truth):
    k = truth.shape[1]
    return np.mean([len(set(p[:k]) & set(t)) / k for p, t in zip(pred, truth)])

def sweep(index, Q, truth, knob, values, k=10):
    for v in values:                                  # THE point: a curve, not a point
        setattr(index, knob, v)
        t = time.perf_counter()
        _, I = index.search(Q, k)
        dt = (time.perf_counter() - t) / len(Q)
        print(f"{knob}={v:5d}  recall={recall_at_k(I, truth):.4f}  {1/dt:8.0f} QPS")

# sweep(ivf,  Q, truth, "nprobe", [1, 2, 4, 8, 16, 32, 64, 128])
# sweep(hnsw, Q, truth, "efSearch", [10, 20, 40, 80, 160, 320, 640])

Ten thousand held-out queries and a few minutes of brute force give you ground truth good enough to make every subsequent decision empirical. Building this first, before choosing an index, is the highest-leverage hour in the whole project.

What an index actually is, as a data structure

It is worth stating the interface, because "index" is used loosely enough to hide real differences.

the contract every index in this lesson implementsclass Index:
    def train(self, sample):
        # learn the partition and/or the codebooks. IVF: k-means over centroids.
        # PQ: m separate k-means runs. HNSW: nothing to train. Cost: minutes to hours.

    def add(self, vectors, ids):
        # assign / encode / link. HNSW inserts one node at a time and is the
        # expensive one here. IVF-PQ's add is embarrassingly parallel.

    def search(self, q, k, **knobs):
        # knobs is the whole point: nprobe for IVF, ef for HNSW, beam width for
        # DiskANN. These move you along the recall/QPS curve WITHOUT rebuilding.

    def remove(self, ids):
        # the operation nobody benchmarks and everybody needs. Often a tombstone.

The line that deserves attention is search(**knobs). An index is not a point on the recall/QPS plane; it is a curve, and the knob is the parameter that traces it. Every serious comparison sweeps that knob. Every marketing claim quotes a single point on it. Chapter 7 is largely about that distinction.

The second line that deserves attention is train. Partition and compression methods have a training stage that sees a sample of the data and learns something about its distribution. That makes them fast and compact — and makes them drift when the data distribution moves. Graph methods have no training stage at all, which is a real and underrated advantage when your corpus is changing under you.

That last row is the one to memorise. Partition and graph methods can always be pushed to exact search by turning the knob all the way up; compression cannot. It is a fundamentally different kind of approximation, and it is why every high-recall production system that uses PQ also has a reranking stage that reads the real vectors.

Chapter 2: IVF From Zero

Start from the crudest possible idea and refine it until it is the algorithm that actually ships.

You cannot look at all the vectors. So: sort them into buckets ahead of time, and at query time look in only a few buckets. That is the whole idea. Everything interesting is in three questions — how do you choose the buckets, how do you decide which ones to open, and what does it cost you when the answer is in a bucket you did not open?

Where the name comes from, because it explains the design

The structure is called an inverted file index (IVF), borrowed directly from text search. In a text engine, an inverted index maps each word to the list of documents containing it — the "posting list." A query for "kalman filter" opens two posting lists and intersects them, instead of reading every document.

The vector version replaces "word" with "region of space." You learn a set of representative points — centroids — and each centroid owns the region of space closer to it than to any other centroid. That region is its Voronoi cell. Every database vector is filed under the centroid whose cell it falls in, and the list of vectors filed under centroid i is the posting list, usually called an inverted list. A query computes which cells it is in or near, and reads only those lists.

Worked example 5 — the cell arithmetic, done fully

Take the concrete setting used by every FAISS tutorial so the numbers are checkable: N = 1,000,000 vectors, d = 128, n_list = 4,096, nprobe = 8.

Step one: average list length. k-means makes cells of roughly equal population where the data is dense and roughly equal volume where it is sparse, but on real data the populations are within a factor of a few. To first order:

Step two: the coarse step. Comparing the query against every centroid is itself a brute-force scan — a small one:

Step three: the fine step. nprobe = 8 lists, each about 244 vectors:

Step four: the total, and the speedup.

Now look hard at Step two versus Step three. The coarse step costs more than twice the fine step. Two thirds of the query is spent deciding where to look. That is not a rounding error — it is a structural fact about IVF that drives two real design decisions, and most tutorials skip straight past it.

Deriving the √N rule of thumb

Everyone repeats "set n_list to about √N." Almost nobody says where it comes from. It falls out of the arithmetic above in three lines.

Total cost, holding nprobe fixed, as a function of n_list:

The first term grows with more cells, the second shrinks. Differentiate with respect to n_list and set to zero:

Check it against the numbers we just used: nprobe = 8, N = 10⁶ gives √(8 × 10⁶) = 2,828. The conventional 4,096 sits right next to it. With nprobe = 1 the formula collapses to exactly √N = 1,000. The folk rule is the solution of a one-line optimisation.

Two honest caveats. First, the optimum is shallow — the cost curve is nearly flat between about ½√(nprobe·N) and 2√(nprobe·N), so anything in that band is fine. Second, this optimises cost at fixed nprobe, not recall. More cells at the same nprobe means a smaller fraction of the corpus scanned, which lowers recall. The real tuning loop fixes a recall target and then finds the cheapest (n_list, nprobe) pair that reaches it.

Training cost, and the warning FAISS prints at you

k-means over n_list centroids needs enough training points to place them meaningfully. FAISS warns below 39 points per centroid and recommends 39–256; the widely used practical target is 100–256 per centroid.

This is why billion-scale IVF is trained on a subsample — typically a few tens of millions of vectors — rather than on the full corpus. It is also why changing n_list is a rebuild, not a tuning knob. nprobe is free to change at query time; n_list is baked in at train time. Knowing which parameters live on which side of that line is most of what "operating a vector index" means.

The failure mode, hand-worked: the Voronoi boundary

Here is where IVF loses recall, made completely concrete in two dimensions so you can check every number with a calculator.

Two centroids: c₁ = (0, 0) and c₂ = (10, 0). The Voronoi boundary between them is the vertical line x = 5. Query q = (4.9, 0.2).

Which cell is q in?

q is closer to c₁, so with nprobe = 1 we open list 1 and only list 1.

Now a database point p = (5.05, 0.25), sitting just on the other side of the boundary:

p was filed under c₂. And yet:

p is 0.158 away from q — overwhelmingly the nearest neighbour, thirty times closer than either centroid — and nprobe = 1 will never see it. It is in the list next door. Raise nprobe to 2 and the answer appears immediately, because c₂ is the second-nearest centroid.

Two refinements the literature added, both worth knowing:

Multi-probe by soft assignment at build time. Instead of filing each vector under one centroid, file it under its two or three nearest. Recall at a given nprobe rises; memory rises with it, because vectors are now duplicated. DiskANN's sharded build (Chapter 5) uses exactly this trick, with each point sent to its two nearest shards.

The inverted multi-index (Babenko & Lempitsky, 2012). Rather than one k-means over n_list centroids, run k-means separately on the first and second half of the vector, with K centroids each. The effective number of cells is K², so K = 2¹⁴ gives 2²⁸ ≈ 268 million cells — a partition far finer than you could ever train directly — while the coarse scan costs only 2K = 32,768 comparisons. If that sounds like the trick in the next chapter, it is: it is product quantization applied to the coarse quantizer.

Why k-means, and what else has been tried

The coarse quantizer's job is narrow: given a query, name a small set of regions that probably contain the answer. k-means is the standard choice, but it is not the only one, and knowing the alternatives clarifies what property actually matters.

The property that separates the winners from the losers in that table is simple: does the partition put roughly the same number of points in each region, and do the region boundaries follow the data rather than the coordinate axes? k-means gets both approximately right because that is literally what it optimises. Everything data-independent gets the first one badly wrong on clustered data, and every axis-aligned method gets the second one wrong in high dimensions.

Residuals — the detail that makes IVF-PQ work at all

One design choice inside IVF matters enormously once compression enters, and it is easy to miss. When an IVF index stores compressed codes, it does not encode the vector x. It encodes the residual:

Why? Because the residual is a much smaller, much better-behaved vector. Suppose your corpus has a spread (standard deviation per coordinate) of σ, and your 4,096 cells cut that spread down so that within a cell the spread is σ/4. The quantization error of any fixed-size codebook scales with the variance of what it is quantizing, so encoding residuals instead of raw vectors reduces squared error by roughly a factor of 16 at the same code length.

Concretely, with an 8-byte PQ code you can either spend all 64 bits describing where in the whole corpus the vector is, or spend them describing where within its cell it is — having got the cell for free from the inverted list id you were storing anyway. The second is obviously better, and the effect is large: encoding residuals rather than raw vectors is worth several points of recall at identical memory.

The cost is a small complication at query time. The distance from q to x is now the distance from (q − c) to r, so the PQ lookup tables must be built per probed list, not once per query. With nprobe = 8 you build eight sets of tables instead of one. FAISS calls this by_residual and it is on by default; turning it off is a real speed/recall knob for very large nprobe.

The arithmetic of soft assignment

Filing each vector under its two nearest centroids instead of one is a cheap way to buy recall, and it is worth costing exactly rather than hand-waving.

At N = 10⁶ and n_list = 4,096, single assignment gives lists of 244. Double assignment stores every vector twice, so the total entries double to 2 × 10⁶ and the average list length becomes 488.

Read the third and fourth rows together. Double assignment at nprobe = 4 scans exactly as many vectors as single assignment at nprobe = 8 — identical query cost — but the vectors it scans are a better-chosen set, because every point near the boundary of the four probed cells is present rather than one door away. You have traded memory (which you may have) for boundary coverage (which is where all your recall was going). That is usually a better trade than raising nprobe, and it is the same insight DiskANN uses when it sends every point to its two nearest shards during the merge-based build.

The obvious follow-up: why not assign to three, or five? Because the memory cost is linear and the recall benefit falls off fast — a point sits near one face far more often than near three. Two is the standard because the second copy captures most of the boundary population.

Living with an index that grows

IVF's real operational advantage over graphs is that adding and removing vectors is trivial: assign, append, done. But "trivial" is not "free of consequences", and the failure is slow enough to miss.

The remedy is unglamorous and should be built on day one: a golden query set with brute-force ground truth, recomputed on every rebuild, with recall published as a deploy-time metric. Recall regressions are otherwise completely silent — the system returns confident, plausible, slightly wrong results, and nothing alerts.

Sim — cells, probes, and the boundary problem

Press "Put the query on a boundary" and then step nprobe from 1 to 2. That single step is the two-dimensional version of the difference between 60% and 90% recall on a real corpus. Then press the recall sweep and watch how the measured recall responds to nprobe — the curve rises fast and then flattens, which is exactly the shape you saw in Chapter 1's triangle explorer and exactly the shape you will see on ann-benchmarks plots in Chapter 7.

The continuum, from exact to reckless

One property of IVF is easy to state and worth internalising: nprobe interpolates continuously between a very fast approximation and exact search. Set nprobe = n_list and every list is scanned, so the result is exactly what brute force would return — slower, because you also paid for the centroid scan, but exact. Nothing about the index is lossy.

This is a real operational lever, not a curiosity. It means the same index can serve a cheap default path at nprobe = 8 and an expensive high-recall path at nprobe = 128 for the queries that deserve it — a paying customer, an internal evaluation job, a retry after an empty result. Graph indexes have the same property through ef. Compression does not: there is no code length you can raise at query time.

List imbalance, the problem the average hides

"Average list length 244" is a comforting number and a slightly dishonest one. k-means minimises squared error, not variance in cluster population, so on real data the lists are far from equal. A dense region of the embedding space — boilerplate, near-duplicate documents, one enormous customer — produces cells holding thousands of vectors while sparse regions produce cells holding a handful.

Why it hurts, concretely. Query latency is proportional to the total length of the probed lists, not to nprobe. If your longest list holds 20,000 vectors instead of 244, then any query that lands in it costs eighty times the average. That is not a tail you can average away — it is a systematic, content-dependent latency spike that will correlate exactly with your most common queries, because the biggest cluster is the most popular topic.

Two standard mitigations. Spherical k-means normalises centroids to unit length at every iteration, which is the correct algorithm when your metric is cosine and which produces noticeably more even cells on normalised embeddings. Balanced assignment adds a penalty proportional to current cell population, so a point near two centroids is nudged toward the emptier one — you trade a little quantization error for a lot of latency predictability, which is usually a good trade in a serving system.

IVF-Flat versus IVF-PQ, and what the choice actually decides

IVF is only the routing layer. What sits inside the inverted lists is a separate decision, and it is the one that determines memory and the recall ceiling.

python — the same corpus, three ways, in FAISSimport faiss, numpy as np

d, N, nlist = 128, 1_000_000, 4096
X = np.random.randn(N, d).astype("float32")

# 1. exact -- the baseline you must beat. No training, no tuning.
flat = faiss.IndexFlatL2(d)
flat.add(X)

# 2. IVF-Flat -- routing only. Exact distances inside the probed cells.
ivf = faiss.index_factory(d, f"IVF{nlist},Flat")
ivf.train(X[:400_000])        # ~100 vectors per centroid; a sample suffices
ivf.add(X)
ivf.nprobe = 8                  # the only query-time knob

# 3. OPQ + IVF-PQ + refine -- routing, compression, and an exact rerank.
pq = faiss.index_factory(d, f"OPQ16_128,IVF{nlist},PQ16", faiss.METRIC_L2)
pq.train(X[:400_000])
pq.add(X)
faiss.extract_index_ivf(pq).nprobe = 8

# Memory: 512 MB / 520 MB / 24 MB. Recall: 1.00 / ~0.98 / ~0.80 before rerank.
# The 20-point gap in row 3 is the quantization ceiling, and a refine stage
# over the top 500 candidates closes almost all of it.

Read the last comment carefully, because it is the most common misdiagnosis in the field: the third index is not "worse." It is a different point on the triangle, chosen because 24 MB fits somewhere 520 MB does not, and its apparent recall deficit is an artefact of measuring it without the reranking stage it was designed to feed.

IVF in one table

Chapter 3: Product Quantization, Derived

Jégou, Douze and Schmid's 2011 paper is the single highest-leverage idea in this lesson. It is what turns "3 terabytes" into "16 gigabytes," and it does it with an idea you can hold in your head completely. This chapter derives it from nothing and then works a complete numerical example by hand, from raw vectors to the final distance table.

Start with plain vector quantization, and watch it die

A quantizer is a function q that maps any vector to one of a finite set of centroids C = {c₁, …, c_k}, called the codebook. The natural choice is nearest-centroid assignment, and the natural way to learn the codebook is k-means — which is precisely the algorithm that minimises the expected squared error E[‖x − q(x)‖²].

Once you have a codebook, you store only the index of the centroid. With k centroids, an index takes log₂(k) bits. So the compression is spectacular in principle: a 128-dimensional float32 vector costs 512 bytes, and its centroid index with k = 256 costs one byte.

Now push it. We want codes long enough to be useful — say 64 bits, so 2⁶⁴ distinct codes. How big is the codebook?

Nine sextillion bytes. Global installed storage capacity is on the order of 10²² bytes, so the codebook for one 64-bit quantizer would consume a meaningful fraction of all the storage humanity has built. And you would have to train it, which needs more data points than centroids.

The product trick

Here is the move. Split the vector into m contiguous chunks of equal length. For d = 128 and m = 8, each chunk is d* = 16 dimensions:

Now quantize each chunk independently, with its own codebook of k* centroids learned by k-means in the subspace. The code for x is the tuple of m sub-indices:

The reconstruction is the concatenation of the chosen sub-centroids: x̂ = [ c_i1⁽¹⁾ | … | c_im^(m) ]. And now count what you have.

What did you pay? Independence. The product quantizer can only represent points that are a concatenation of per-subspace centroids; it cannot express correlations between subspaces. If dimensions 3 and 97 are strongly correlated and land in different chunks, PQ cannot exploit that. Section "OPQ" below is the fix.

Worked example 6 — the complete hand calculation

Everything below uses d = 8, m = 2, d* = 4, and k* = 4. Real PQ uses k* = 256 so that each sub-index is exactly one byte; we shrink it to 4 (two bits per sub-index) purely so the whole thing fits on a page and you can check every number. Nothing about the structure changes.

The two codebooks, learned by k-means on the first four and last four dimensions respectively:

The database vector to encode:

Encoding, subspace 1. x⁽¹⁾ = (0.9, 0.1, 0.2, 1.1). Squared distance to each centroid, term by term:

The winner is index 0, at squared error 0.07.

Encoding, subspace 2. x⁽²⁾ = (2.8, 0.2, 0.1, 3.2):

The winner is index 1, at squared error 0.13.

The reconstruction and its error. x̂ = (1, 0, 0, 1 | 3, 0, 0, 3), and the total squared error is exactly the sum of the two per-subspace errors — the subspaces are disjoint coordinates, so squared errors add:

Asymmetric distance computation, hand-worked

Now the query arrives:

The obvious thing to do — decode every stored code back to a full vector and compute distances — would throw away the entire speed advantage. Instead, asymmetric distance computation (ADC): leave the query in full precision, and precompute a small table of the distance from each query subvector to each centroid in that subspace. Then the distance to any stored code is a sum of table lookups.

Table for subspace 1, from q⁽¹⁾ = (1.0, 0.2, 0.1, 0.9):

Table for subspace 2, from q⁽²⁾ = (2.9, 0.1, 0.2, 3.1):

Two tables, four entries each: eight numbers, computed once per query. Now the estimated squared distance to any stored vector is two lookups and one addition. For our x with code (0, 1):

No multiplications at all. That is the operational payoff and it is easy to under-appreciate: scanning a million PQ codes is a million table lookups and adds over an 8-byte-per-vector stream, versus a million 128-dimensional dot products over a 512-byte-per-vector stream.

How good is the estimate? Compute the truth. q − x = (0.1, 0.1, −0.1, −0.2 | 0.1, −0.1, 0.1, −0.1):

Estimated 0.13, true 0.11. The estimate is 18% high. That direction is not an accident, and the exact decomposition explains it.

Why ADC is biased upward, exactly

Write the true distance in terms of the reconstruction, by adding and subtracting x̂:

The first term is exactly what ADC computes, 0.13. The last is the quantization error, 0.20. Check the middle term by hand. q − x̂ = (0, 0.2, 0.1, −0.1 | −0.1, 0.1, 0.2, 0.1) and x − x̂ = (−0.1, 0.1, 0.2, 0.1 | −0.2, 0.2, 0.1, 0.2):

The identity closes exactly. Now the statistical reading: over many database points, the residual (x − x̂) points in an essentially arbitrary direction relative to (q − x̂), so the cross term averages to zero and

ADC over-estimates distances by, on average, the mean squared quantization error. The paper notes you can subtract that constant to get an unbiased estimator — and also that it does not matter for ranking, because subtracting the same constant from every candidate changes nothing. What does matter is the variance: the noise floor it adds is what caps recall. Shorter codes, larger quantization error, higher noise, lower ceiling. That is Chapter 1's ceiling, now with a mechanism.

The alternative that loses: symmetric distance computation

SDC quantizes the query too, then looks up precomputed centroid-to-centroid distances. Encode our query: q⁽¹⁾ is nearest to c⁽¹⁾₀ (squared distance 0.06, the smallest entry in T₁) and q⁽²⁾ is nearest to c⁽²⁾₁ (0.07). So code(q) = (0, 1) — the same code as x. Therefore:

Estimated distance zero, true distance 0.332. Worse, every database vector whose code is (0, 1) now gets distance exactly 0, so SDC cannot rank within a cell at all — it has thrown away all the resolution that the full-precision query still had. The paper's conclusion follows immediately: ADC has strictly lower distance-estimation error, and SDC is worth using only when the query itself must be stored compressed (for instance when you are indexing the queries too).

Cost accounting for a real scan

Back to d = 128, m = 8, k* = 256, scanning one IVF list of 244 vectors, or all 1,953 vectors of an nprobe = 8 search.

Notice the last two rows are the real story. The arithmetic is 16× cheaper, which is nice; the memory traffic is 64× cheaper, which is decisive on a bandwidth-bound workload. And the tables themselves are 8 × 256 = 2,048 floats = 8 KB — they live in L1 for the entire scan.

Modern implementations push this further with 4-bit PQ and SIMD lookups. With k* = 16 instead of 256, a whole sub-table is 16 bytes and fits in a single SIMD register, so the lookup becomes a vector shuffle instruction (vpshufb on x86, tbl on ARM) that scores 16 or 32 codes per instruction rather than one at a time. FAISS calls these "fast scan" indexes and writes them as PQ64x4fs; they are several times faster than the 8-bit version at slightly lower per-code accuracy, which is usually recovered by using twice as many subquantizers at the same total bytes.

OPQ — fixing the independence assumption for free

PQ assumes the chunks are independent and, implicitly, that the variance is spread evenly across them. Real data violates both. If your embedding has most of its energy in the first 30 coordinates — which is exactly what happens after any PCA-like transform, and approximately what happens in many learned embeddings — then chunk 1 has huge variance and chunk 8 has almost none. You have spent equal bits on both, which is a bad allocation.

Optimized product quantization (Ge, He, Ke & Sun, 2013; independently "Cartesian k-means") learns an orthogonal rotation matrix R and encodes R·x instead of x, choosing R to balance variance across the chunks and decorrelate them. Because R is orthogonal it preserves all distances, so nothing is lost — it merely presents the data to PQ in a friendlier basis.

For d = 128 that is 16,384 multiply–adds once per query — half the cost of building the PQ tables, and utterly negligible against the scan. In exchange you typically gain two to five points of recall. It is one of the highest-return-per-line changes available, which is why FAISS index strings so often start with OPQ.

The ceiling, and the rerank that breaks it

PQ's recall ceiling is real and unavoidable: two different vectors with the same code are indistinguishable, forever. The standard escape is a two-stage search:

The recall of the final answer is now bounded by the recall of the shortlist at r, not at k — and shortlist recall at r = 500 is dramatically higher than top-10 recall from codes alone, because PQ's estimates are noisy enough to reorder within the top 10 but rarely bad enough to eject a true neighbour from the top 500. This is why "IVF-PQ has bad recall" is almost always a statement about a missing rerank stage rather than about PQ.

Worked example 7 — a billion vectors on one machine

Chapter 0's impossible corpus, re-costed with PQ at d = 128 (the classic SIFT1B / BIGANN setting):

A 32× factor between the first and last row, and a 32× factor again between the last row and the second. The choice of m is the single biggest lever you have on the cost of a large vector service, and it moves recall smoothly: reported BIGANN results put m = 8 in the rough neighbourhood of 0.2–0.4 for 1-recall@1 without reranking and considerably higher for recall@100, with each doubling of m buying a substantial chunk of what remains. Add reranking and the numbers move again. Measure on your own data — the shape holds, the exact values do not transfer.

Sim — build the distance table yourself

Query A is the one worked by hand above; its ADC estimate for vector x₁ should read 0.06 + 0.07 = 0.13 against a true 0.11. Query B is chosen to sit between two centroids in the second subspace, where the quantization error is large and the ranking genuinely breaks — that is the failure PQ actually exhibits in production, and seeing it on six vectors makes the million-vector version unsurprising.

Choosing m, as a bit-allocation problem

m is the one parameter that trades memory against the recall ceiling, and it is worth reasoning about rather than copying. Two constraints bound it from either side.

From below: each subquantizer must have enough dimensions to be meaningful. With k* = 256 centroids covering a d*-dimensional subspace, the centroids per dimension is 256^1/d*. At d* = 4 that is 4 centroids' worth of resolution per axis — already coarse. At d* = 2 it is 16, which sounds better but wastes bits describing a nearly-empty plane. The community converged on d* between 4 and 16 for a reason.

From above: memory. Each subquantizer costs one byte (at k* = 256), so m is your byte budget.

Product quantization for inner product, which is not the same thing

Everything above computed squared Euclidean distance. If your objective is maximum inner product, the table is built differently — and, pleasingly, more simply. Because the inner product decomposes additively over the coordinates just as squared distance does:

So the lookup table holds inner products instead of squared distances, and the per-candidate operation is still m lookups and m−1 adds. The only change is the sign of what you sort by. FAISS exposes this as METRIC_INNER_PRODUCT, and it is the substrate on which Chapter 6's anisotropic weighting is applied — ScaNN keeps this scoring machinery and changes only the objective that trains the codebooks.

Write it yourself — product quantization in twenty lines

python — a complete, working PQ. Reproduces the hand calculation exactly.import numpy as np
from scipy.cluster.vq import kmeans2

class PQ:
    def __init__(self, m, ks=256):
        self.m, self.ks = m, ks

    def fit(self, X):
        n, d = X.shape
        self.ds = d // self.m                          # d* -- must divide exactly
        self.cb = np.zeros((self.m, self.ks, self.ds), "float32")
        for j in range(self.m):                       # m INDEPENDENT k-means runs
            sub = X[:, j*self.ds:(j+1)*self.ds]
            self.cb[j], _ = kmeans2(sub, self.ks, iter=20, minit="points")
        return self

    def encode(self, X):                              # -> (n, m) uint8 codes
        codes = np.empty((len(X), self.m), "uint8")
        for j in range(self.m):
            sub = X[:, j*self.ds:(j+1)*self.ds]
            # (n, ks) squared distances to this subspace's centroids
            D = ((sub[:, None, :] - self.cb[j][None])**2).sum(2)
            codes[:, j] = D.argmin(1)
        return codes

    def tables(self, q):                             # -> (m, ks), built ONCE per query
        return np.stack([((q[j*self.ds:(j+1)*self.ds] - self.cb[j])**2).sum(1)
                         for j in range(self.m)])

    def adc(self, T, codes):                         # -> (n,) estimated squared distances
        # m gathers and a sum. No multiplications touch the database at all.
        return sum(T[j][codes[:, j]] for j in range(self.m))

# Sanity check against Worked Example 6: set m=2, ks=4, plug in the two
# codebooks by hand, and adc() must return 0.06 + 0.07 = 0.13 for x1.

That is the entire technique. The adc method is four tokens of numpy and it is the thing that lets a laptop hold a billion vectors. Type it out once and the rest of this lesson stops being abstract.

Product quantization in one table

Chapter 4: HNSW — Navigation Instead of Partition

IVF asks "which region?" and then scans. Graph methods ask a different question: "which neighbour should I step to?" — and answer it repeatedly until stepping stops helping. This chapter builds HNSW from two much older ideas, both of which you can hold entirely in your head, and then traces a complete search by hand.

Idea one: the skip list

Pugh's skip list (1990) solves a problem that looks unrelated: searching a sorted linked list. A linked list holding 1, 4, 7, 9, 13, 18, 22, 30 costs O(n) to search — you must walk every node, because a linked list has no random access.

The fix: build express lanes above it. Each element is promoted to the next level up with probability p (typically ½), independently. Level 0 has everything; level 1 has about half; level 2 about a quarter; and so on.

Searching descends: start at the top-left, walk right while the next element is still ≤ your target, then drop down a level and repeat. Each level you skip past roughly half the remaining elements, so the search is O(log n) expected.

Trace it. Sixteen elements, p = ½. Expected occupancy: level 0 = 16, level 1 = 8, level 2 = 4, level 3 = 2, level 4 = 1. Looking for 13:

Four comparisons for sixteen elements, against an expected eight for a plain scan. At a million elements the gap is 20 versus 500,000.

Idea two: navigable small worlds

In a metric space you cannot sort, so "walk right" is replaced by greedy routing: from your current node, look at its neighbours, move to whichever is closest to the query, and stop when no neighbour is closer than where you stand.

Whether that works depends entirely on the graph's edges. Connect each node only to its ten nearest neighbours and greedy routing takes O(N^1/d) hops — you crawl across the space one small step at a time. Connect nodes at random and greedy routing has no gradient to follow at all.

Kleinberg's 1999 analysis of small-world networks gives the sharp answer. On a d-dimensional lattice with local links plus one long-range link per node drawn with probability proportional to r^−d (r being the distance), greedy routing finds the target in O(log² N) steps. Change the exponent in either direction — more long links or fewer — and the bound degrades to a polynomial. Navigability is not a property of "having some long edges"; it is a property of having them at every scale in the right proportion.

Malkov's NSW (2011–2014) got that link distribution for free from an incremental construction: insert points one at a time, connecting each new point to its M nearest among those already inserted. Early insertions happen when the graph is sparse, so their "nearest neighbours" are far away — they become the long-range links. Later insertions land in a dense graph and produce short links. The scale distribution emerges from the insertion order.

What was still wrong. Flat NSW routes in a polylogarithmic number of hops, but the first phase — getting from a random entry point into the right neighbourhood — is a long, low-information walk over high-degree nodes, and its cost grows with N. The long links and the short links are tangled in one graph, so every hop has to consider both.

HNSW: put the long links on their own floor

Malkov & Yashunin's insight (arXiv:1603.09320) is to separate the scales explicitly, exactly as a skip list does. Build a hierarchy of graphs:

The upper layers are a coarse-to-fine zoom; layer 0 is where the accuracy is bought. The number of layers is O(log N) in expectation, the work per upper layer is O(1) in expectation, and the layer-0 beam search dominates the cost.

The knobs, and exactly what each one does

The memory arithmetic, done properly. N = 1,000,000, d = 128, M = 16, 4-byte neighbour ids.

Layer 0 gives every node up to 2M = 32 links:

Upper layers: the expected number of nodes summed over all layers above zero is

each with up to M = 16 links:

Plus the vectors themselves at 512 B each = 512 MB. Total 644 MB, or 644 bytes per vector, of which 132 bytes (20.5%) is graph. Push M to 48 and the graph alone becomes 384 bytes per vector, three quarters of the payload.

The neighbour-selection heuristic — the part everyone skips

When inserting element e, you run a search to collect efConstruction candidates and must choose M of them as neighbours. The obvious choice — the M nearest — produces a graph that fails in a specific, instructive way.

Picture two dense clusters, A and B, with a gap between them. Insert a point in cluster A. Its M nearest candidates are all in cluster A, because A is dense. Every point in A links only within A. The same happens in B. The result is a graph with no edges crossing the gap: greedy routing that enters A can never reach B, and recall for queries near B that should return points in A collapses. The graph is locally excellent and globally disconnected.

HNSW's Algorithm 4 fixes this with a diversity rule. Walk the candidates in increasing distance from e. Accept a candidate c only if

In words: keep c only if it is closer to me than to anything I have already agreed to link. If c is closer to an existing neighbour r than to me, then routing from me toward c can go via r, so the direct edge is redundant — and spending a precious slot on it costs me an edge in a direction I do not yet cover.

The effect is that the M chosen neighbours point in M different directions rather than all clustering in the densest one. In the two-cluster picture, once the first few in-cluster neighbours are accepted, further in-cluster candidates get pruned, and a point across the gap survives — because nothing already accepted is closer to it than e is. The bridge gets built.

Worked example 8 — a complete search trace, by hand

Twelve points in the plane. Coordinates:

Layers: layer 2 = {A, H}; layer 1 = {A, C, F, H, J}; layer 0 = all twelve. Entry point = A. Query q = (8.2, 8.4), k = 3, ef = 3.

The true answer, computed by brute force: L at 1.000, H at 1.414, J at 2.236.

Layer 2 (edges: A—H). Start at A.

H improves enormously, so move to H. H's only layer-2 neighbour is A, already worse. Greedy stops. One hop covered 89% of the distance to the query. That is what a sparse top layer is for.

Layer 1 (edges: A—C, C—F, C—H, F—J, F—H, H—J). Enter at H, d = 1.414. Neighbours of H at this layer are C, F, J:

None beats 1.414. Greedy stops immediately; H remains the entry point for layer 0. Cost: three distance computations.

Layer 0 (edges: A—B, A—C, B—C, B—G, C—D, C—F, D—E, D—I, E—F, E—H, E—I, F—G, F—J, G—K, G—J, H—I, H—J, H—L, J—K, J—L). Beam search with ef = 3.

Initialise: candidate queue = {H: 1.414}, result set W = {H: 1.414}.

Round 1. Pop the closest candidate, H. Expand its neighbours E, I, J, L:

W has room for 3. Insert L (1.000) and J (2.236); W = {L: 1.000, H: 1.414, J: 2.236}. E at 4.561 and I at 6.450 are worse than W's worst (2.236) and W is full, so they are discarded and not queued. Candidate queue = {L: 1.000, J: 2.236}.

Round 2. Pop L (1.000). Its neighbours are H and J — both already visited. Nothing changes.

Round 3. Pop J (2.236). The stopping rule is "halt when the closest unexpanded candidate is farther than W's worst"; 2.236 equals W's worst, so we proceed. J's neighbours are F, G, H, K, L:

All worse than 2.236; discarded. H and L already visited. Queue is empty. Halt.

The cost. Distances computed: A, H at layer 2; C, F, J at layer 1; E, I, L (J already known) at layer 0 from H; then F, G, K from J. That is ten distinct nodes out of twelve. B and D were never touched.

Setting M and efConstruction without guessing

The two build parameters are usually copied from a blog post. They deserve twenty minutes of thought, because they are the ones you cannot change later without a rebuild.

Memory is a closed-form function of M, so start there and work backwards from your budget. With 4-byte neighbour ids:

The first term is the vector, the second is layer 0's 2M links, and the third is the expected upper-layer links (recall that N/(M−1) nodes exist above layer 0, each with up to M links). Evaluate at d = 768:

efConstruction has a cleaner rule, because it costs nothing at query time: raise it until the recall/QPS curve stops moving. In practice 100 is thin, 200 is a good default, 500 is thorough, and beyond that the curve is usually indistinguishable while the build takes twice as long. The one signal that you set it too low: a curve that flattens below your target recall no matter how high ef goes. That is a graph whose edges are simply wrong, and no amount of query-time beam width repairs it — the same "ceiling versus effort" diagnosis as Chapter 1, one level up.

Sim — watch the layers do their jobs

Why greedy routing works at all — the ideal graph, and the one we can afford

Greedy routing has an obvious failure mode: you can walk into a local minimum, a node closer to the query than all of its neighbours but not the global answer, and stop. What property of the graph rules that out?

The exact answer is known. On the Delaunay graph — the graph connecting points whose Voronoi cells share a face — greedy routing is guaranteed to reach the true nearest neighbour from any starting point, with no local minima. That is a beautiful theorem and a useless one, because in d dimensions the Delaunay graph's average degree grows exponentially with d. At d = 128 nearly every point is a Delaunay neighbour of nearly every other point, and the graph is the complete graph — storing it costs more than the data and traversing it is a full scan.

Because it is an approximation, local minima do exist, and this is exactly what ef fixes. A beam of width ef does not stop at the first node with no improving neighbour; it keeps ef candidates alive and continues expanding until none of them can improve the result set. A local minimum traps a beam of width 1 and is stepped over by a beam of width 16. That is the mechanism behind the recall/QPS curve: ef is buying you escapes from local minima, and each additional escape is rarer and more expensive than the last — which is exactly why the curve saturates.

Insertion, which is just a search with a write at the end

HNSW has no training stage, and its build is a loop of searches. Worth walking once, because every build-time parameter lives here.

HNSW insertion — Algorithms 1 and 4 of the paper, in wordsdef insert(e, M, efConstruction, mL):
    l = floor(-ln(uniform(0,1)) * mL)     # this element's top layer
    ep = entry_point

    # PHASE 1: zoom down through layers ABOVE l with a width-1 greedy walk.
    # Cheap: these layers are sparse, and we only need a good entry point.
    for lc in range(top_layer, l, -1):
        ep = greedy_search(ep, e, layer=lc, ef=1)

    # PHASE 2: at each layer from l down to 0, do a real beam search and link.
    for lc in range(min(top_layer, l), -1, -1):
        W  = search_layer(ep, e, layer=lc, ef=efConstruction)   # candidates
        nb = select_neighbours_heuristic(e, W, M)               # the DIVERSITY rule
        link(e, nb, layer=lc)                                  # bidirectional
        for n in nb:                                           # repair over-full nodes
            if degree(n, lc) > (2*M if lc == 0 else M):
                relink(n, select_neighbours_heuristic(n, neighbours(n, lc),
                                                     2*M if lc == 0 else M), lc)
        ep = W
    if l > top_layer: entry_point, top_layer = e, l

Three things this makes visible. One: efConstruction only appears in Phase 2, so it costs build time and nothing at query time — there is never a reason to be stingy with it beyond your patience. Two: the relink step means inserting one element can rewrite its neighbours' edge lists, which is why parallel insertion needs locking and why HNSW builds do not scale linearly with cores. Three: the same diversity heuristic is applied both when choosing e's neighbours and when repairing an over-full node — so the property is maintained, not merely established.

The build cost follows: N insertions, each a search of roughly the same cost as a query at ef = efConstruction, plus the repair work. At efConstruction = 200 that is roughly 200/ef_query times the cost of a query, per inserted vector — which is why a million-vector build is minutes and a billion-vector build is not something you do on a laptop.

The ef ≥ k rule, and the mistake it causes

One small constraint trips up nearly everyone once. The beam is the result set, so an index cannot return more results than its beam held. Setting ef = 10 and asking for k = 100 does not return 100 items; it returns 10, or throws, depending on the library.

The consequence for a reranking pipeline is easy to miss. If you fetch 500 candidates for a reranker, ef must be at least 500 — and the recall you should be measuring is 10-recall@500, not 10-recall@10. Teams routinely set ef = 64 out of habit, request 500, silently receive 64, and then conclude the reranker is not helping. Check the length of what actually came back, on every path, at least once.

The implementations you will actually meet

HNSW is one algorithm and half a dozen widely used implementations, and they differ in ways that matter more than the parameter values you copy between them.

Where HNSW hurts

Against all of that: when the data fits in RAM, HNSW is very hard to beat. On the standard million-scale benchmarks it sits at or near the top of the recall/QPS frontier, reaching recall above 0.95 while evaluating a fraction of a percent of the corpus, with one intuitive runtime knob and no training stage to go stale. If your corpus is under ten million vectors and mostly static, the correct first move is HNSW, and the correct second move is to measure before doing anything cleverer.

Chapter 5: DiskANN — When the Vectors Live on SSD

Subramanya, Devvrit, Kadekodi, Krishnaswamy and Simhadri's NeurIPS 2019 paper makes a claim that sounds impossible after Chapter 4: one billion points, a single machine with 64 GB of RAM, more than 5,000 queries per second, 95%+ 1-recall@1, mean latency under 5 milliseconds. This chapter derives how, and the derivation is almost entirely arithmetic about a storage device.

Start with why the obvious thing fails

The naive plan: build an HNSW graph over a billion points and keep the vectors on SSD, paging them in as the walk needs them. Cost it out with d = 128.

The memory that would be needed if you kept it all resident — Chapter 4's arithmetic at a billion:

What happens if you page instead. A random 4 KB read from a good NVMe SSD takes about 100 microseconds at queue depth 1. A graph search visits, say, 120 nodes to reach 95% recall. And — this is the crucial part — the visits are sequentially dependent: you cannot know which node to fetch next until you have read the current one and computed distances to its neighbours.

Twelve milliseconds is already over a tight latency budget, and it is a floor: no amount of CPU helps, and the drive is doing 120 IOPS worth of work per query at queue depth 1, so throughput is dreadful too. Paging a graph index is not a small inefficiency; it is a different performance regime.

The three-part answer

Worked example 9 — the memory budget, line by line

A billion 128-dimensional points on a 64 GB machine.

Note how neatly the node record fits the hardware. 772 bytes is under a 4 KB sector, so a single read never straddles a boundary; DiskANN packs several nodes per sector and reads whole sectors, which is what the drive wants anyway. If d were 768, the node record would be 3,072 + 256 = 3,328 bytes — still one sector, but only one node per sector, which is exactly why high-dimensional corpora often get dimension-reduced (or Matryoshka-truncated) before a disk index is built.

Vamana, and the α that makes it different

Vamana builds a directed graph with maximum out-degree R (typically 64–128). The loop is simple:

Vamana, in words and then in numbers# 1. Start from a random R-regular graph, or an empty one.
# 2. Compute the MEDOID s = the point minimising the sum of distances to all
#    others (approximated on a sample). This is the fixed entry point for
#    every query forever -- a deliberate choice, unlike HNSW's top layer.
# 3. For each point p, in random order:
#      V = GreedySearch(s -> p)         # the set of nodes VISITED, not just the top-k
#      out(p) = RobustPrune(p, V, alpha, R)
#      for each j in out(p):            # keep the graph roughly undirected
#          out(j) = out(j) + {p}
#          if |out(j)| > R: out(j) = RobustPrune(j, out(j), alpha, R)
# 4. Run the whole pass TWICE: first with alpha = 1, then with alpha = 1.2.

RobustPrune is where α lives. Given a candidate set V, sort by distance from p and repeatedly take the closest surviving candidate p*, add it to the output, then delete every remaining candidate v satisfying

Read that as a detour test. The edge p → v is redundant if you could get from p to v by hopping through p* and the detour is not much longer than the direct edge. With α = 1 the rule is exactly HNSW's diversity heuristic (and the classical relative-neighbourhood graph): prune v whenever the detour is no longer than the direct edge. With α > 1 the rule becomes more permissive: keep the direct edge unless the detour is shorter by a factor of α.

The hand-worked pruning decision. Take p = (0, 0), the closest accepted neighbour a = (1, 0) so d(p, a) = 1, and a candidate v = (1.5, 2.4).

Now apply the rule at two settings:

The ratio d(p,v)/d(a,v) = 2.830/2.452 = 1.154. Any α above 1.154 keeps this edge; anything below prunes it. Geometrically, v sits far from p but sideways relative to a, so routing through a saves only 13% of the distance. α = 1.2 says: unless the detour saves you 20%, keep the direct road.

The other deliberate difference from HNSW: Vamana is flat. No layers. The long edges that the hierarchy was providing are now inside the single graph, put there on purpose by α rather than emerging from insertion order. A flat graph is much easier to lay out on disk — there is exactly one record per point — and the fixed medoid entry point removes the need to store or traverse upper layers at all.

The medoid, and why one fixed entry point is enough

HNSW spends a hierarchy on the problem of "where should a search start?" Vamana answers it with a single fixed point and no hierarchy at all, which deserves an explanation rather than an assertion.

The medoid is the data point minimising the total distance to all other points — the most central actual point, as distinct from the centroid, which is an average and generally not in the dataset. Computing it exactly is quadratic, so it is estimated on a sample:

estimating the medoid — enough precision for an entry pointS = sample(X, 100_000)              # a sample is plenty; this need not be exact
c = S.mean(axis=0)                    # the centroid -- not necessarily a real point
medoid = argmin(dist(x, c) for x in S)   # the real point nearest to it

Why is one entry point sufficient? Because α-pruning has already put long-range edges into the graph deliberately. In flat NSW the entry point mattered enormously, since reaching a distant region required many short hops. In a Vamana graph the first few hops from the medoid are long by construction — that is precisely what α = 1.2 retained — so a search reaches any region of the space in a handful of steps regardless of where it started.

Three consequences, all of them practical:

The cost is that the medoid must be recomputed when the data distribution shifts substantially — another item for the rebuild checklist, and a cheap one.

Beam search: turning dependent reads into parallel ones

Greedy search is a chain: read a node, learn its neighbours, pick one, read it. Each link in the chain is a full SSD latency, and the drive is idle in between. The fix is to widen the frontier.

Beam width W means: at each round, take the W closest unvisited candidates and issue all W reads simultaneously. The drive processes them concurrently — NVMe is built for exactly this, with deep hardware queues — so W reads cost about the same wall-clock time as one.

Worked example 10 — latency and throughput versus beam width. Assume a search needs about 120 node visits for 95% recall, and a 4 KB random read costs 100 µs.

Two things fall straight out. First, beam width is a latency-versus-throughput dial, not a recall dial: it makes each query faster while making the fleet slower, because wider beams re-read nodes the narrow search would have skipped. Second — and this is the satisfying part — the paper's headline of "over 5,000 QPS" is not a mystery number. A single high-end NVMe drive sustains roughly 500,000 to 1,000,000 random 4 KB IOPS. Divide:

The headline throughput is the drive's IOPS rating divided by the graph's hop count. That is the whole system in one division, and it tells you exactly which two numbers to improve: buy more IOPS, or build a graph with fewer hops. α-pruning is the second option.

Two-stage scoring: approximate for routing, exact for ranking

During the walk, when a node's neighbours are considered, their distances are estimated from the in-memory PQ codes — free, no I/O. Only nodes actually selected for expansion have their sectors read, and those reads deliver full-precision vectors, which are used to compute exact distances for the result set.

So the search is guided by cheap approximations and answered with exact numbers. The quantization error affects which nodes get visited (a routing decision, where being slightly wrong usually still lands you in the right neighbourhood) but never the final ranking (where being wrong is fatal). This is the deepest structural idea in the paper, and it is the same asymmetry Chapter 3 described for ADC, applied one level up.

Worked example 11 — how much RAM the codes really need, and why m is chosen there

DiskANN's RAM budget is dominated by one line: the PQ codes. So m is not a compression decision here, it is a machine size decision, and the two constraints pull in opposite directions.

The tension is unusual and worth naming: shorter codes save RAM but cost disk I/O, because a badly routed walk visits more nodes and every extra node is a 100 µs read. Compression here is not a pure win against a recall ceiling as it is in IVF-PQ — the final ranking uses exact vectors from disk, so recall is not capped by m at all. What m controls is how many sectors you have to read to find the answer. That is a genuinely different role for the same technique, and it is why DiskANN's m (32) is larger than a typical IVF-PQ m (8–16) despite both being "PQ."

Building a billion-point graph on a machine that cannot hold it

Vamana's construction is itself memory-hungry — it needs random access to vectors. The paper's answer is a merge-based build:

The reported build for a billion SIFT descriptors took on the order of two days on a single 64 GB machine. That is the honest headline cost of DiskANN: queries are cheap, construction is a batch job you schedule. If your corpus turns over daily, this is the wrong index — which is precisely why FreshDiskANN followed, adding a small in-memory index for recent writes that is periodically merged into the on-disk graph, the same read-optimised-plus-write-buffer split that LSM-tree storage engines use.

The on-disk record, byte by byte

The layout is the design, so it is worth seeing as bytes. For d = 128 and R = 64:

The co-location is not a micro-optimisation; it halves the number of round trips, which is the only currency here. Separating vectors and adjacency into two files would mean two dependent reads per hop, doubling latency and halving throughput at a stroke.

Note also what is not in the record: nothing about layers, because Vamana is flat, and no PQ code, because those live in RAM in one contiguous array indexed by node id. The two representations of the same point are physically separate and used at different moments.

Caching the top of the graph — a free 20–30% of the I/O

Every single query starts at the medoid and takes its first few hops through the same small set of nodes. Those nodes are read on every query, which is exactly the access pattern a cache is for.

The arithmetic. Cache every node within 4 hops of the medoid. With out-degree R = 64 the reachable set grows fast, but the graph is not a tree and overlap is heavy; in practice a few million nodes covers it. At 772 bytes per node:

— against a 64 GB budget already holding 32 GB of PQ codes, so it fits comfortably. And it removes the first 3–4 hops of every query from the I/O path entirely. If a query takes 120 hops, saving 4 sounds trivial — but those 4 are sequential and at the head of the dependency chain, and with beam width 4 they represent a full round of latency each. The paper reports this cache as a real contributor, and it is the cheapest optimisation in the system.

The search loop, with the two precisions made explicit

DiskANN beam search — where each precision is useddef search(q, k, L, W):            # L = candidate list size, W = beam width
    cand = PriorityQueue([(pq_dist(q, medoid), medoid)])   # RAM, approximate
    seen, out = {medoid}, []

    while cand.has_unexpanded():
        # pick the W closest UNEXPANDED candidates -- all by PQ, zero I/O
        batch = cand.closest_unexpanded(W)

        # ONE round trip for all W. This is the whole point of the beam.
        sectors = ssd.read_many([node_offset(n) for n in batch])

        for n, sec in zip(batch, sectors):
            vec, nbrs = parse(sec)
            out.append((exact_dist(q, vec), n))     # EXACT -- for the answer
            for m in nbrs:
                if m not in seen:
                    seen.add(m)
                    cand.push((pq_dist(q, m), m))     # APPROXIMATE -- for routing
        cand.truncate(L)

    return sorted(out)[:k]                      # ranked by EXACT distances only

Every call to pq_dist is a table lookup over data already in RAM and decides only where to go. Every call to exact_dist uses bytes that arrived from the SSD and decides only what to return. Trace those two functions through the loop and you have the paper.

The IOPS budget under concurrency

Worked example 10 costed a single query. Serving many at once changes which resource binds, and the arithmetic is worth finishing.

Assume 120 reads per query, beam width 4, a drive rated at 600,000 random 4 KB IOPS and 3 GB/s of sequential bandwidth.

Notice that IOPS runs out before bandwidth does, by a comfortable margin — 4 KB reads at 600k IOPS is only 2.46 GB/s against a 3 GB/s rating. That is the signature of a random-access workload, and it tells you which drive specification to shop on. It also explains why widening the beam is not free at scale: beam width 8 needs roughly 60% more reads per query, so the same drive supports proportionally fewer queries per second even though each one finishes faster.

When DiskANN is the wrong answer

Updates, and why the paper needed a sequel

Everything above assumes the graph is built once and then read forever. Real corpora are not like that, and the gap between "billion-point search on one node" and "billion-point search on one node with a live write path" took a second paper to close.

The three hard operations, in increasing difficulty:

The practical planning rule that falls out: your rebuild cadence must be shorter than the time it takes the delta to become a significant fraction of the base. If 1% of a billion points changes per day, the in-memory delta is 10 million points — a perfectly reasonable in-memory Vamana index — and a weekly compaction keeps it that way. If 20% changes daily, no compaction schedule saves you and a disk-resident graph is the wrong architecture.

DiskANN against HNSW, honestly

DiskANN is slower per query than in-memory HNSW, by roughly an order of magnitude, and that is the correct trade: it is 20× cheaper in hardware for the same corpus. Nobody runs DiskANN on a million vectors. Everybody should consider it before provisioning a 40-machine HNSW fleet.

Chapter 6: ScaNN — Quantize for the Decision, Not the Vector

Every quantizer so far has minimised the same thing: reconstruction error, E[‖x − x̂‖²]. It is the objective k-means minimises, it is what PQ inherits, and it feels obviously right — a good compression is one that loses little.

Guo, Sun, Lindgren, Geng, Simcha, Chern and Kumar's ICML 2020 paper (arXiv:1908.10396) points out that it is not obviously right, because reconstruction error is not the thing you care about. You care about whether the ranking survives. And it turns out that some reconstruction errors damage the ranking severely while others, of exactly the same magnitude, do not damage it at all.

First, the objective this is really about

ScaNN targets maximum inner product search (MIPS): given q, find argmax_x ⟨q, x⟩. This is the native objective of two-tower recommenders and retrieval models, where the score is a dot product and the vectors are deliberately not normalised — norm carries information (popularity, confidence, document length).

Quantization enters as an approximation to the score. You store x̂ instead of x, so you compute ⟨q, x̂⟩ instead of ⟨q, x⟩. Write the residual r = x − x̂. Then the error in the score is exactly

So the quantity we actually want small is ⟨q, r⟩ — the projection of the residual onto the query — not ‖r‖. Those are very different demands, and the difference is a direction.

The decomposition, built slowly

We do not know q at training time. But we know something about which q matters. A database point x only ever affects a result list when it scores highly, and it scores highly only for queries pointing roughly along its own direction. Queries pointing elsewhere give x a low score, x is nowhere near the top-k, and an error in its score changes nothing.

Make that precise. Split the residual into its component along x and the rest:

Now take the extreme case where the query is exactly aligned with x, q = x / ‖x‖. Then ⟨q, r_⊥⟩ = 0 by construction, and the entire score error comes from the parallel part:

Move q away from x's direction and the perpendicular part starts to contribute — but in high dimensions, a query that scores x highly cannot be far from x's direction, because inner products with random directions concentrate near zero (Chapter 0's result, arriving for the third time). So conditioned on "x is a plausible top-k answer for q," the parallel component dominates the error, and it dominates more as the dimension grows.

This is why the paper's derived weights depend on both the dimension and a score threshold T: the higher the bar for "this item matters," the tighter q is pinned to x's direction, and the more the parallel weight should exceed the perpendicular one. In practice the whole thing collapses to a single tunable ratio.

Set h_∥ = h_⊥ and you recover ordinary k-means exactly, since ‖r_∥‖² + ‖r_⊥‖² = ‖r‖². Anisotropic quantization is a strict generalisation with one extra knob, and the knob has a clear meaning: how much more do I fear a length error than a sideways error?

Worked example 12 — two errors of identical size, one fatal

This is the whole paper in eight lines of arithmetic. Take a two-dimensional database point:

Consider two candidate reconstructions, chosen so that their reconstruction errors are exactly equal.

Candidate A — a pure length error. x̂_A = (3.6, 4.8):

Candidate B — a pure sideways error. x̂_B = (3.8, 3.4):

Ordinary k-means is completely indifferent between A and B. Both have squared error 1.0. Now score them with a query aligned to x — the query for which x is a genuine top result. Take the unit vector q = (0.6, 0.8):

Identical reconstruction error; one is perfect and the other is off by a full unit on a score of 5. And now the other half of the argument — check a query perpendicular to x, where B's error should show up. Take q′ = (0.8, −0.6):

B is wrong by 1.000 here. But look at where: it turned a score of 0 into a score of 1, on an item that was never going to be in the top-k of that query. Compare with A's failure, which turned a score of 5 into a 6 — on a query where x is a contender, and where inflating the score by 20% can push x past a genuinely better item. The same magnitude of error costs everything in one place and nothing in the other.

The loss, evaluated. With h_⊥ = 1 and h_∥ = η:

At η = 4, candidate A is charged four times as much as B, and any k-means-style optimisation over this loss now prefers B decisively. Change the objective, and the centroids move to different places — places that preserve rankings rather than places that preserve coordinates.

Why the weighting grows with dimension

The paper derives h_∥ and h_⊥ in closed form from an assumption that queries are uniform on the sphere, conditioned on the score exceeding a threshold T. The intuition is one you already have from Chapter 0.

A uniformly random unit vector in d dimensions has inner product with any fixed direction distributed with standard deviation 1/√d. At d = 768 that is 0.036. So conditioning on "this query scores x above a meaningful threshold" is an extremely tight constraint: it selects queries lying in a narrow cap around x's direction. Inside that cap, the perpendicular components of q are small, so ⟨q, r_⊥⟩ is small, while ⟨q, r_∥⟩ is nearly the full ‖r_∥‖.

Which gives a clean, checkable prediction: anisotropic quantization should help more in high dimensions and more when you only care about a short result list — and less when d is small or you need a long, well-ordered tail. That is exactly the reported behaviour.

Worked example 12b — the threshold, and why it sets the weight

The claim that "queries aligned with x are the ones that matter" deserves numbers rather than assertion. Stay in two dimensions with x = (3, 4), ‖x‖ = 5, and parameterise the query by the angle φ it makes with x. Take a unit query q(φ), so the true score is

Candidate A's residual is parallel with ‖r_∥‖ = 1, so its score error is ⟨q, r_∥⟩ = 1 · cos φ. Candidate B's is perpendicular with ‖r_⊥‖ = 1, so its score error is 1 · sin φ. Tabulate:

Read the table as a weighting problem. B's error is largest exactly where it does not matter (φ near 90°, score near zero) and vanishes exactly where it would (φ near 0°). A's error does the opposite: it is at its maximum precisely on the queries for which x is the answer. Weight each row by the probability that x is actually a contender — which is concentrated at small φ — and A's expected damage dwarfs B's.

Now bring in dimension. In two dimensions, a uniformly random query has a decent chance of landing at small φ. In 768 dimensions the angle between a random unit vector and any fixed direction concentrates sharply around 90° — the cosine has standard deviation 1/√768 = 0.036, so a random query sits at φ ≈ 88–92° essentially always. Conditioning on "the score exceeds threshold T" therefore selects an astonishingly narrow cap around x's direction, and inside that cap the table above is dominated entirely by its top row.

That is the paper's closed-form result, arrived at by counting. And it produces a falsifiable prediction, which is the mark of a real derivation: anisotropic weighting should help most in high dimensions with a short result list, and least in low dimensions with a long one. Which is what the ablations show.

What it does to the codebook

It is worth being concrete about how the learned centroids differ, because "a weighted loss" sounds like it should change little.

Plain k-means places a centroid at the mean of its assigned points. Averaging vectors that point in slightly different directions produces a result that is shorter than the vectors themselves — the classic shrinkage of averaging. Shortening is a residual along the data direction, which is exactly the parallel component. So ordinary k-means has a systematic bias toward the single error type that damages inner products most, and it has it structurally, not by accident.

The anisotropic loss penalises that shrinkage heavily, and the fitted centroids push outward relative to their k-means positions — better preserving the norms and radial positions of the points they represent, at the cost of tolerating more sideways error. The reconstruction error goes up. The ranking quality goes up too. If you ever measure a quantizer by mean squared error and find that the better-ranking one scores worse, this is why, and it means the metric is wrong rather than the quantizer.

Worked example 13 — a ranking that actually flips

The previous example showed an error of 1.0 on a score of 5.0. That only matters if it changes an ordering, so let us make it change one.

Two database items and one query, in two dimensions:

The true scores put y ahead:

Now quantize each of them with a unit-length residual, but in different directions. For x, take the parallel-error reconstruction from Worked example 12:

For y, build a unit residual that is exactly perpendicular. A vector orthogonal to y = (3.4, 4.2) is (−4.2, 3.4); normalise it by ‖y‖ = 5.404 to get the unit perpendicular u = (−0.777, 0.629). Set ŷ = y − u:

Both reconstructions carry exactly the same reconstruction error. Score them:

The estimated ranking now puts x above y, reversing the truth. Both reconstructions have exactly unit error, so a reconstruction-error objective sees them as equally good approximations — and yet one of them silently swapped the top two results. Repeat that across a corpus and you have an index whose top-1 accuracy is far worse than its mean squared error suggests.

The anisotropic loss charges x̂ η times more than ŷ and pushes the codebook away from exactly this failure. That is the entire mechanism, and it is why the paper's improvement shows up most strongly in recall@1 and least in recall@100 — deep in the list the ordering errors are between items nobody looks at.

Why rescoring is not optional here

A fair objection: if the final stage rescores candidates with exact inner products, does the quantizer's ranking error matter at all? It matters for exactly one thing, and it is the thing that determines recall — which candidates reach the rescoring stage.

If quantization pushes the true best item out of the top-r shortlist, no rescoring can bring it back. Chapter 1's shortlist-recall framing applies verbatim: the quantizer's job is not to rank correctly, it is to not eject the winner from the candidate pool. Anisotropic weighting improves precisely that, because the errors it suppresses are the ones capable of demoting a genuine top item by a large amount.

Worked example 13b — the reduction people try instead, and why it fails

There is a well-known trick for turning MIPS into ordinary nearest-neighbour search, and it is worth working through because its failure is instructive and because you will meet people proposing it.

Let M = max_x ‖x‖. Augment every database vector with one extra coordinate, and the query with a zero:

Every augmented database vector now has norm exactly M, because ‖x‖² + (M² − ‖x‖²) = M². And so:

The first two terms are constants, so minimising the augmented Euclidean distance is exactly maximising the original inner product. A clean reduction: now use any L2 index you like. Why does nobody ship this?

The arithmetic. Take M = 10 and two database vectors with ‖x₁‖ = 2 and ‖x₂‖ = 3. Their extra coordinates are

Let q be a unit vector with ⟨q, x₁⟩ = 1.8 and ⟨q, x₂⟩ = 1.5. The augmented squared distances are

Correct ordering — and look at the margin. The entire signal is 0.6 out of roughly 97.7, a relative difference of 0.6%. Everything an index would exploit — cell membership, graph structure, quantized codes — is now dominated by a constant term contributing 101 units, with the answer hiding in the sixth part of a percent.

What ScaNN is, as a shipped system

The anisotropic loss is the paper's contribution; the library around it is a full three-stage pipeline, and the stages should look familiar by now.

When it was published, this combination took the top of the ann-benchmarks glove-100-angular leaderboard, roughly doubling the throughput of the previous best at comparable recall. The partitioning and the SIMD scoring were engineering; the anisotropic loss was the idea, and it is the part that transfers.

The rest of the system, briefly

The anisotropic loss is applied twice, at two granularities, which is easy to miss:

And the reason the library is fast rather than merely accurate is the 4-bit LUT16 scoring from Chapter 3: with k* = 16, each subquantizer's lookup table is 16 bytes, which fits in a single SIMD register, so a shuffle instruction scores many candidates at once. Accuracy from the loss, throughput from the register width, correctness from the rescoring pass. Three separate ideas, and only the first one is the paper.

The generalisable lesson

"But my vectors are normalised" — does any of this apply?

A fair question, since most modern text encoders emit unit vectors and cosine similarity. If every ‖x‖ = 1, is there a parallel error to worry about?

Yes, and the reason is that the reconstruction is not normalised. The database vector x lies on the unit sphere; its quantized reconstruction x̂ is a sum of codebook centroids and lands wherever it lands — generally inside the sphere, because quantization is a form of averaging and averaging shortens. So the residual r = x − x̂ has a radial component even when x itself carries no norm information.

If the shrinkage were identical for every x it would cancel in the ranking and nothing would matter. It is not identical: vectors in dense regions are reconstructed accurately and barely shrink, while vectors in sparse regions shrink more. So the radial error is a density-dependent bias that systematically demotes items in sparse parts of the space — which, unhelpfully, is exactly where the distinctive, specific documents live.

Setting the weight in practice

The paper gives h_∥ and h_⊥ in closed form from a score threshold T, but the practical recipe is shorter and you should know both.

Step 3 is where people abandon the technique by accident. If your quantizer evaluation dashboard reports mean squared reconstruction error — and most do, because it is the easy number — then anisotropic quantization looks like a strict regression and gets reverted. Delete that metric or clearly mark it as not the objective.

The same idea, in four other places

"Optimise the compression for the decision, not for the reconstruction" is a template, and recognising it elsewhere is most of this chapter's long-term value.

Where ScaNN fits — and where it does not

Chapter 7: Benchmark Literacy

There is no theorem you can cite instead of measuring, so the professional skill in this field is reading measurements correctly. This chapter is about the plots, the traps in them, and the single most important question nobody asks: should I be using an index at all?

The plot, and how to read it

ann-benchmarks (Aumüller, Bernhardsson & Faithfull) is the community's shared yardstick. Its protocol is deliberately austere: standard datasets with precomputed brute-force ground truth, each algorithm in its own container, single-threaded, one query at a time, and every algorithm swept across its full parameter grid.

The canonical chart is recall on the x-axis, QPS on a logarithmic y-axis. Up and to the right is better. Every point on a curve is one parameter setting; the curve is generated by sweeping the runtime knob.

How to actually read two curves:

Three deliberate choices in the standard protocol

ann-benchmarks makes three decisions that surprise people, and each one is defensible and consequential.

The first row is the one to hold onto. A single-threaded QPS number is a fine comparative measure and a poor absolute one. Your production question is queries per second per dollar under concurrency, and the ordering of algorithms on that axis can differ from the ordering on the published plot — particularly for anything whose advantage is parallelism.

The billion-scale sibling, big-ann-benchmarks, changes the rules deliberately: it fixes a hardware budget, allows batch queries, and separates tracks for in-memory, out-of-core, and custom-hardware systems. That is a different and equally valid framing — "what is the best system under this budget" rather than "what is the best algorithm per core" — and the winners differ accordingly. When someone quotes a leaderboard, ask which question it was answering.

The axes the plot does not have

Four quantities decide real deployments and appear on no standard chart.

Build time. HNSW on a million vectors is minutes; a billion-point Vamana graph is days. If your corpus is re-embedded every time you ship a new model, build time is the system's cost. It can differ by three orders of magnitude between two curves that look identical.

Memory. A curve that reaches 0.99 recall at 20,000 QPS while using 900 bytes per vector loses to one at 8,000 QPS and 24 bytes per vector, if the corpus is large enough that the first one needs a fleet.

Mutability. Insert cost, delete semantics, and how gracefully the index degrades between rebuilds. Benchmarks build once and query forever. Almost no production system does.

Tail latency. QPS is one over the mean. Graph searches have variable hop counts, so p99 can be several times p50 — and if you shard across S machines and merge, the query waits for the slowest of S, so the fan-out p99 approaches the p99 of the max, which is much worse than the p99 of one machine.

Recall definitions, one more time, because this is where papers cheat

If a paper reports 1-recall@10 and you need 10-recall@10, the number does not transfer. If you are building a pipeline with a reranker, measure recall at the shortlist size, not at the final k — measuring the wrong one is the most common way teams conclude their index is bad when their metric is.

Worked example 14 — when brute force wins, calculated

The most valuable arithmetic in this chapter, because it saves teams from building infrastructure they do not need. Exact search is memory-bandwidth-bound (Chapter 0), so the latency is one division:

Take a single modern core streaming at roughly 20 GB/s, and a whole server at 200 GB/s:

And exact search scales further than people expect on a GPU, where bandwidth is an order of magnitude higher. An accelerator with 1.5 TB/s of memory bandwidth scans a million 768-dimensional vectors in

— exactly, with no index, no build, no staleness, and trivially correct metadata filtering. FAISS's GPU flat index is a genuinely competitive production choice into the tens of millions.

The trap benchmarks hide completely: filtered search

Real queries almost never say "nearest neighbours." They say "nearest neighbours where tenant = 42 and updated_at > last month and language = 'de'." No standard benchmark measures this, and it breaks graph indexes badly.

Let s be the selectivity — the fraction of the corpus passing the filter. Three strategies:

The post-filter arithmetic. With selectivity s = 0.001 and a search returning 100 candidates, the expected number that survive the filter is

Nine times out of ten you return an empty page. To expect ten survivors you would need to retrieve 10,000 candidates — which for most indexes costs more than scanning the 0.1% subset exactly. Notice this is not a tuning problem; it is arithmetic, and it tells you the answer: at low selectivity, pre-filter and brute force.

There is a subtler failure for graph indexes even with generous post-filtering. Greedy routing follows edges, and if most nodes fail the predicate, the walk keeps landing on ineligible nodes. The graph's navigability was built for the whole corpus, and the filtered subgraph may not be connected at all — so recall over the filtered set can collapse far below what the unfiltered recall would suggest. Partitioning the index by the high-cardinality filter dimension (a separate index per tenant) is unglamorous and usually correct.

Measuring throughput without fooling yourself

Latency and QPS numbers are easy to produce and easy to produce wrongly. Five mistakes account for nearly all of the bad numbers in circulation.

The second row deserves emphasis, because it is the most common serious error and it produces beautiful numbers. If your query set is a random sample of the indexed vectors, the nearest neighbour of each query is that query, sitting at distance exactly zero, and any index that touches it at all will find it. Recall comes out near 1.000 for every configuration, the curves look flat and wonderful, and the whole experiment measured nothing.

Recall is a proxy. Sometimes measure the thing itself

Recall against brute-force ground truth answers "does my index agree with exact search?" That is the right question for tuning an index and the wrong question for deciding whether the product got better, because exact search is not the ground truth of usefulness — it is the ground truth of your embedding model.

Ground truth for filtered queries, which nobody computes and everybody needs

If your production queries carry predicates, then unfiltered ground truth measures a system you do not run. The fix is mechanical and takes an afternoon:

That last plot is the useful artefact, and it almost always looks the same: recall is fine for selectivity above a few percent, and falls off a cliff below it. The cliff's location is where you should switch from post-filtering to pre-filtering, and it is a number you can measure rather than argue about. Wiring the switch is then a one-line policy: if the estimated matching-set size is under some threshold, scan it exactly; otherwise search and post-filter.

Sim — the tuning playground

Press "Target recall 0.99" and compare the QPS the two families give up to get there. Then drag N down to 10⁴ and watch the brute-force line rise above both curves entirely — at which point the correct index is a for-loop.

Why public datasets flatter every algorithm

SIFT and GloVe are a decade old and were not chosen because they resemble your data. They were chosen because they are public, million-scale, and everyone already had them. Four differences matter enough to change your conclusions.

The last row is the one that ambushes people. Real corpora contain thousands of near-identical chunks — licence headers, navigation text, templated product copy. They collapse into a single dense knot that dominates one IVF cell, fills every graph neighbourhood in that region with clones, and makes the ground-truth top-10 a set of ten essentially identical items whose ordering is numerical noise. Recall computed against that ground truth measures nothing useful. Deduplicate before you benchmark, or your numbers will describe your boilerplate.

The three plots to make, and what each one answers

Plot 2 is the one that changes decisions in practice, because it converts recall into money. Draw it once with your real corpus size and the argument about whether to use PQ ends immediately in either direction.

Reading four curves, in order

Suppose you have swept four candidates on your own data and produced one plot. Here is how the reasoning actually goes, in the order it should go.

Three habits follow. Fix the recall target before looking at the plot, so the crossing points resolve themselves. Never compare a compressed index without its rerank stage to an uncompressed one — you are comparing a component to a system. And always annotate each curve with bytes per vector, because the y-axis you can see is not the axis that will decide.

A recall regression, debugged

The most common real-world task is not choosing an index; it is explaining why the one you have got worse. A playbook, in the order that finds the cause fastest:

One more trap: the k inside your ground truth

Ground truth is computed for some fixed k — usually 100, because that covers most downstream uses. Two things go wrong with that constant.

You cannot measure recall above the k you computed. If ground truth holds the true top 100 and you want 10-recall@500, the metric is still well defined (are the true top 10 inside your 500?) but 100-recall@500 is not, because you do not know items 101–500 of the truth. Compute ground truth at the largest k any consumer will ever need, once, and slice down as required.

Ties at the boundary make recall jittery. If items 100 and 101 sit at nearly identical distances — overwhelmingly likely with near-duplicate content — then which one is "the truth" is decided by floating-point noise, and recall@100 fluctuates by a point between runs for no reason. The clean fix is to record the ground-truth distances alongside the ids and count a returned item as correct if its distance is within a tiny epsilon of the k-th true distance. It costs nothing and removes an entire category of phantom regressions.

A benchmarking checklist you can actually run

Chapter 8: The Composed Stack

Nobody ships "an index." They ship a pipeline in which each stage narrows the candidate set and the next stage scores what survives more expensively than the last. Every technique in this lesson is a stage in that pipeline, and the composition is the actual engineering artefact.

The pattern, stated once

The canonical five-stage stack

Read the funnel widths. Each stage costs roughly the same total time as the others, because each one is (many more items) × (much cheaper per item). That balance is not an accident either — it is what you get when you tune a pipeline properly, and a stage that costs 10× the others is a stage that is misconfigured.

The funnel, with real numbers on every arrow

Abstract funnels are unconvincing. Here is the same five-stage pipeline on 100 million vectors at d = 768, with the width and the cost of each stage.

Read the "items in" column top to bottom: 65,536 → 24,400 → 500 → 50. Each stage is roughly two orders of magnitude narrower than the one before and roughly two orders of magnitude more expensive per item. That is not a coincidence — it is what a tuned funnel looks like, and a stage that violates it is misconfigured.

Reading a FAISS index string, which is the real literacy test

FAISS encodes an entire pipeline in one string. Being able to decode one on sight is the most compact demonstration that you understand this lesson. Take:

That one line is: rotate, partition with a graph-accelerated coarse quantizer, scan 32-byte codes with SIMD, rerank exactly. Four of this lesson's chapters, composed, in 27 characters.

Three more strings, for practice

Decoding index strings is the fastest way to check whether the composition really landed. Cover the right column and work through these.

The fifth row is worth dwelling on. Scalar quantization — one byte per dimension, learned only as a min and max per dimension — gets you 4× compression with almost no accuracy loss and none of PQ's machinery. It is not glamorous and it is frequently the correct first compression step. Reach for PQ when 4× is not enough.

Worked example 15 — the same corpus, costed three ways

100 million vectors from a modern text encoder, d = 768, cosine, target recall@10 ≥ 0.95, 500 QPS.

Option A — pure HNSW in RAM.

Latency around 1–2 ms, recall easily above 0.95, one machine, no reranking needed. At roughly $4/hour for a memory-optimised 384 GB instance that is about $35,000 per year, plus a replica for availability.

Option B — OPQ + IVF-PQ + rerank. Truncate to 256 dimensions with a Matryoshka-trained encoder, then PQ to 64 bytes:

Plus full-precision vectors on SSD or object storage for the rerank stage — 307 GB of cold bytes, of which a query reads perhaps 500 × 3 KB = 1.5 MB. Latency around 3–6 ms. At roughly $0.25/hour that is about $2,200 per year, plus cheap storage.

Option C — brute force on a GPU. 307 GB does not fit in accelerator memory, so this needs sharding across several devices or streaming from host memory. Sharding across 4 accelerators with 80 GB each: exact search, zero recall loss, per-query scan of 307 GB across 4 devices at 1.5 TB/s each

Far too slow for 500 QPS, and vastly more expensive. GPU brute force is a great answer at 10⁶–10⁷; at 10⁸ with d = 768 it has run out of road.

A 16× cost difference between A and B for two points of recall and four milliseconds. Which is correct depends entirely on what those two points are worth — Chapter 1's compounding argument is how you find out. The point of this table is that "which index?" is a budget question with a numeric answer, not a matter of taste.

Cost per query, in actual money

The last conversion that makes architecture arguments settle: put every stage in dollars per million queries. Using the 100M-vector Option B configuration and round cloud prices.

The corollary for capacity planning: the memory-optimised choice matters most when the corpus is large enough that the first row is not fourteen cents. Redo the table at your own N — if the search row is dominated by everything else, take the simplest index that meets your recall target and spend the attention elsewhere.

GPU graph indexes: CAGRA

One more branch worth knowing. NVIDIA's CAGRA (in the cuVS library) is a graph index designed for GPUs rather than ported to them. The key difference is in the search: a CPU beam search expands one node at a time with a narrow frontier, which suits a few fast cores; CAGRA instead evaluates a large batch of candidates per iteration so that hundreds of threads have work, trading a less efficient search order for massive parallelism.

The reported results are large speedups in build time over CPU HNSW — graph construction is embarrassingly parallel and had been the painful part — and substantial throughput gains at high recall for batched queries. There is also a practical bridge: a CAGRA graph can be converted into an HNSW-compatible structure, so you can build on a GPU in minutes and serve on CPUs. If your bottleneck is rebuild cadence rather than query latency, that conversion is the interesting part.

Its limitation is the obvious one: accelerator memory is small and expensive, so CAGRA lives happily up to tens of millions of vectors, and beyond that you are back to partitioning and compression.

The write path, which is half the system and gets a tenth of the attention

Everything so far described the read path. The write path decides your rebuild cadence, your staleness, and most of your operational pain, and it has a shape of its own.

The staleness budget follows directly. If step 2 embeds 2,000 chunks per second and a document produces 20 chunks, you ingest 100 documents per second. A corpus of 100 million documents therefore takes 10⁶ seconds — 11.6 days — for a full re-embed on one worker. That single number decides whether a model upgrade is an afternoon or a quarter, and it is almost never in the capacity plan.

Worked example 16 — a capacity plan on one page

Target: 100M documents, 20 chunks each = 2 × 10⁹ vectors at d = 768. 2,000 QPS, p99 under 150 ms end to end, full corpus refresh monthly.

Choosing, as a decision procedure

The millisecond budget of a real query

Vector search is never the whole request, and looking at where the time actually goes reorders most optimisation priorities. A representative retrieval-augmented request, 100 million documents, OPQ+IVF-PQ with rerank:

Hybrid retrieval, because dense search alone is not the system either

One more composition, because it is nearly universal in production and it is a different axis from everything so far. Dense vector search fails in a specific, predictable way: exact tokens. Product codes, error identifiers, rare proper nouns, version numbers — an embedding smooths precisely the surface detail that makes those queries answerable, so "error TS2571" retrieves documents about TypeScript errors in general and not the one page that names it.

Sparse lexical retrieval (BM25 over an inverted index — the original meaning of the term this chapter has been borrowing) fails in the complementary way: it cannot match paraphrase. The standard answer is to run both and fuse the rankings, most often with reciprocal rank fusion, which needs no score calibration between two systems whose scores are not comparable:

Worked micro-example. A document ranked 1st by BM25 and 30th by the dense index, against one ranked 8th by both, with K = 60:

Doc B wins. That is the intended behaviour: RRF rewards agreement across independent systems over a single spectacular result, which is exactly the right prior when one system is prone to lexical over-matching and the other to semantic drift. The additive form also makes the fusion trivially extensible — a third retriever joins by adding a term.

The systems consequence for this lesson: your vector index does not need to be the only recall path, which means its recall target can be lower than you feared. A dense index at 0.93 recall inside a hybrid system can produce better end-to-end results than one at 0.99 alone, at a fraction of the memory. Chapter 1's warning that the right recall target is a property of the pipeline, not the index, cashes out here.

What breaks first, as you scale

A rough ladder. Each rung is where a real system typically stops working, and the fix is usually architectural rather than a parameter.

Operational realities the papers do not cover

Chapter 9: Connections and Cheat Sheet

Four papers, fourteen years, one problem. Here is everything on one page, followed by where each thread continues.

The four papers, in one table each

The lesson in ten paragraphs, one per chapter

The symbol table

The equations that carry the weight

The numbers worth remembering

Glossary — every bolded term, in one line each

What to memorise, and what to look up

The left column is a page you could write from memory in ten minutes and it decides nearly every architectural question in this space. The right column changes every release and should never be memorised.

The failure modes, in one place

The knob cheat sheet

Nine things that are widely believed and wrong

Where each thread continues

Build it yourself — the weekend recipe

How to read the four papers

All four are short. If you read them in this order, each one answers a question the previous one raised.

Two supporting reads that repay the time: Johnson, Douze & Jégou on FAISS, which is the engineering context for all four, and the ann-benchmarks paper, which is the methodology you will be judged by.

The first three weeks, if you are building this now

If week 1 ends with "brute force is fine," that is a successful outcome and you should write it down so the question does not get reopened every quarter. Include the corpus size at which you expect it to stop being fine, and the arithmetic behind it.

References

1. Jégou, H., Douze, M., Schmid, C. "Product Quantization for Nearest Neighbor Search." IEEE TPAMI 33(1):117–128, 2011. The origin of everything in Chapter 3.
2. Malkov, Y. A., Yashunin, D. A. "Efficient and Robust Approximate Nearest Neighbor Search Using Hierarchical Navigable Small World Graphs," 2016 — arXiv:1603.09320; IEEE TPAMI 42(4), 2020. Chapter 4.
3. Subramanya, S. J., Devvrit, Kadekodi, R., Krishnaswamy, R., Simhadri, H. V. "DiskANN: Fast Accurate Billion-point Nearest Neighbor Search on a Single Node." NeurIPS, 2019. Chapter 5.
4. Guo, R., Sun, P., Lindgren, E., Geng, Q., Simcha, D., Chern, F., Kumar, S. "Accelerating Large-Scale Inference with Anisotropic Vector Quantization." ICML, 2020 — arXiv:1908.10396. Chapter 6.
5. Ge, T., He, K., Ke, Q., Sun, J. "Optimized Product Quantization." IEEE TPAMI, 2014. The rotation that fixes PQ's independence assumption.
6. Babenko, A., Lempitsky, V. "The Inverted Multi-Index." CVPR, 2012. Product quantization applied to the coarse quantizer.
7. Johnson, J., Douze, M., Jégou, H. "Billion-scale Similarity Search with GPUs," 2017 — arXiv:1702.08734. FAISS, and the reference implementation of most of this lesson.
8. Malkov, Y., Ponomarenko, A., Logvinov, A., Krylov, V. "Approximate Nearest Neighbor Algorithm Based on Navigable Small World Graphs." Information Systems, 2014. The flat NSW that HNSW improved on.
9. Fu, C., Xiang, C., Wang, C., Cai, D. "Fast Approximate Nearest Neighbor Search With The Navigating Spreading-out Graph" (NSG). VLDB, 2019. The monotonic-graph line of work that Vamana's pruning descends from.
10. Indyk, P., Motwani, R. "Approximate Nearest Neighbors: Towards Removing the Curse of Dimensionality." STOC, 1998. LSH, and the (1+ε) framework Chapter 1 argues against using in practice.
11. Aumüller, M., Bernhardsson, E., Faithfull, A. "ANN-Benchmarks: A Benchmarking Tool for Approximate Nearest Neighbor Algorithms." Information Systems, 2020. Chapter 7's plots.
12. Simhadri, H. V. et al. "Results of the NeurIPS'21 Challenge on Billion-Scale Approximate Nearest Neighbor Search," 2022 — arXiv:2205.03763. The billion-scale tracks.
13. Ootomo, H. et al. "CAGRA: Highly Parallel Graph Construction and Approximate Nearest Neighbor Search for GPUs," 2023 — arXiv:2308.15136. Chapter 8's GPU branch.
14. Singh, A., Subramanya, S. J., Krishnaswamy, R., Simhadri, H. V. "FreshDiskANN: A Fast and Accurate Graph-Based ANN Index for Streaming Similarity Search," 2021 — arXiv:2105.09613. Updates, which Chapter 5 flagged as the missing piece.