What it actually does
Do not build a model of the world. Keep every example you have seen, and when a new point arrives, find the handful that look most like it and copy their answer — a majority vote for a class, an average for a number. Everything hard about the method is hidden in "look most like it": the distance function is the model, and choosing it badly is a modelling error that no amount of tuning k will repair. The trade is stark and worth stating up front: training costs nothing, and inference costs everything.
- Similar inputs have similar outputs — the smoothness assumption, and the only one the method actually makes about the world
- The distance metric reflects the similarity that matters for the task, which is a modelling claim and not a technical detail
- Features are scaled comparably, since distance is measured in feature units
- The intrinsic dimension is low enough that distances still discriminate, whatever the nominal feature count says
- The training set covers the region queries arrive in; there is no extrapolation mechanism whatsoever
Loss function
There is no objective minimized at training time — this is the one model in the category that fits nothing. What the formula states is the prediction rule: over the k nearest points, return the answer with the lowest weighted loss. Substituting squared error makes the minimizer the weighted mean, and substituting 0-1 loss makes it the weighted majority vote, which is why the same algorithm serves regression and classification without modification. It is local empirical risk minimization over a neighbourhood, computed fresh for every query.
- the indices of the k training points closest to x under the chosen metric
- the loss for answering c when the truth is y_i — squared error for regression, 0-1 for classification
- optional distance weight; uniform weights give the plain vote or mean
- a small constant so a query landing exactly on a training point does not divide by zero
How it is minimized
None — the model is the training set. The engineering problem is the neighbour search, not a fit
- a candidate set of k training indices; the rule picks the set with the smallest total distance
- Mahalanobis-form distance; M = I is Euclidean, M diagonal is per-feature scaling, general M is a learned metric
- the matrix that encodes what "similar" means — the closest thing kNN has to a parameter vector
- the stored training set, which is the model in its entirety
There is nothing to minimize by gradient or by factorization, so the usual optimization vocabulary does not apply — which is exactly why kNN is worth studying next to a linear model. What replaces it is a set of decisions with real consequences: the metric M (a diagonal M is what feature scaling secretly is, and a learned M is metric learning), the value of k, and the data structure used to find neighbours. A KD-tree gives O(log n) queries in low dimensions and collapses to a full scan somewhere around twenty features; a ball tree tolerates general metrics; above that, exact search stops being the goal and an approximate index is the honest answer.
k
The bias-variance dial. k = 1 has zero training error and enormous variance; large k smooths toward the global majority (3 to 25, cross-validated; odd values for binary classification to avoid ties)
distance metric
The actual model. Euclidean for continuous scaled features, cosine for embeddings and text, Hamming or Gower for mixed and categorical data
weighting
Uniform or inverse-distance. Distance weighting lets a larger k stay locally sensitive, which is usually preferable to a small k
feature scaling
Not optional. An unscaled feature dominates the distance in proportion to its units, so the metric silently encodes your choice of measurement
index structure
Brute force, KD-tree, ball tree, or an approximate index — a latency decision, and above roughly twenty dimensions the only one that matters
No iteration and therefore no convergence, but a real asymptotic guarantee: as n grows with k/n approaching zero, the kNN error rate approaches the Bayes error, and the 1-NN error is bounded by twice it. That guarantee is asymptotic in a way that matters, because the rate depends exponentially on dimension. The characteristic failure mode has a name — the curse of dimensionality — and a precise statement: as dimension grows, the ratio between the nearest and farthest distances in a sample converges to one, so every point becomes equally near and "nearest neighbour" stops carrying information. It does not announce itself; accuracy simply decays toward the base rate while the code keeps working perfectly.
Training: O(1), or O(n log n) to build a tree index. Query, brute force: O(nd) time and O(nd) memory resident for the whole training set. KD-tree query: O(log n) in low dimensions, degrading to O(n) as dimension rises. Approximate indexes trade a recall guarantee for sublinear query time, which is the trade production systems almost always take.
How, where, and why it is used
Time-Series Forecasting
AdaptedAnalog forecasting: embed the series into vectors of the last m observations, find the historical windows most similar to the current one, and forecast by averaging what happened next after each of them. The metric runs over the delay embedding rather than over raw features, so the model asks "when has the recent past looked like this before?" and answers with the empirical continuation.
- Nonlinear and chaotic series where a global functional form fits badly but local structure repeats — hydrology, some load and flow series
- Motif-based short-horizon forecasting on high-frequency data with strong repeating patterns
- A non-parametric sanity baseline against a fitted model: if analog forecasting is competitive, the dynamics are more repetitive than the fitted model assumes
It makes no assumption about the functional form of the dynamics, which is a genuine advantage where the series is nonlinear and the sample is long enough for the past to contain relevant analogues. It fails exactly where those conditions fail: a short history has no analogues, and a trending series has none either, because every recent window is in a region the past never visited. Extrapolation is impossible by construction — the forecast is always a weighted average of things that already happened.
- Delay embedding: stack the last m lags into the query vector, with m chosen from the autocorrelation structure
- Detrend and deseasonalize first, so similarity is judged on shape rather than on level
- Scale the embedding dimensions, since distance across lags is otherwise dominated by the highest-variance ones
- Restrict the candidate pool to matching seasonal phase where the seasonality is strong
Rolling-origin backtesting with an expanding candidate pool, scored with MASE against seasonal-naive. Report how many analogues were actually within a meaningful distance — a forecast built from k distant neighbours is an average of unrelated windows dressed up as a prediction.
- Allowing the neighbour search to see windows from after the forecast origin, which is a leak that looks like skill
- Trend: the current window sits outside the range the past covered, so the nearest analogues are all systematically lower or higher
- Reporting a forecast without reporting neighbour distance, so a no-analogue situation is indistinguishable from a confident one
Anomaly Detection
Primary fitThe distance to the k-th nearest neighbour is itself the anomaly score: a point in a dense region has close neighbours, a point in empty space does not. Local Outlier Factor refines this by comparing a point’s local density to the density of its own neighbours, which is what lets it flag a point that is sparse relative to its cluster even though a global distance threshold would call it normal.
- Fraud and intrusion detection where anomalies are defined by dissimilarity to any known-normal pattern
- Manufacturing and process monitoring with clusters of distinct normal operating regimes
- Data-quality screening, where the anomalies of interest are records unlike any other record
- Novelty detection in embedding space — flagging inputs far from the training distribution of a deployed model
It is the honest default when normal behaviour is multi-modal and you cannot write down its shape. A parametric detector must assume a form for normality; this one only assumes that normal points have neighbours, which is a far weaker commitment. The reasons not to use it are equally concrete: the score costs a full neighbour search per point, and in high dimensions the distances it depends on stop separating, which is the same curse that ends its use as a classifier.
- Scale on a confirmed-clean window and persist those statistics; the metric is the detector
- Reduce dimension first — PCA or a learned embedding — since distance-based scoring degrades well before the dimension count looks alarming
- Choose k above the size of the smallest legitimate cluster, or a small genuine mode gets scored as a set of outliers
Precision@k against confirmed incidents with the threshold taken from score quantiles on clean data, plus PR-AUC rather than ROC-AUC given the class imbalance. Hold out whole anomaly episodes, never individual points from inside one, or neighbouring points from the same episode make each other look normal.
- A cluster of anomalies masks itself: with more than k of them together, they are each other’s neighbours and none scores as isolated
- Contaminated reference data quietly redefines normal to include the thing being hunted
- Query cost grows with the reference set, so the detector gets slower exactly as it accumulates the history that makes it good
Optimization
Not applicableThere is no objective to optimize and no decision to produce under constraints — the method deliberately fits nothing, and its only search problem, finding the k nearest points, is an indexing question rather than a decision one.
Recommendation & Ranking
Primary fitItem-item collaborative filtering is kNN with the interaction matrix as the feature space and cosine similarity as the metric: to recommend for a user, take the items they engaged with, retrieve each item’s nearest neighbours, and rank by aggregated similarity. User-user filtering is the same computation transposed.
- Item-item collaborative filtering, still a production default for "customers who bought this also bought"
- Candidate generation ahead of a learned ranker, where recall matters far more than the ordering
- Cold-start-by-content: nearest neighbours in a content-embedding space for items with no interaction history
Item-item neighbourhoods are stable, cheap to precompute, and give an explanation for free — "because you watched X" is the model’s actual mechanism rather than a post-hoc story. It stops being enough when interactions are sparse and the signal has to be shared across items, which is what matrix factorization and two-tower models exist to do; but as a candidate generator underneath them it remains hard to beat.
- Cosine or adjusted-cosine similarity rather than Euclidean, since interaction vectors differ in magnitude by popularity
- Down-weight popular items, which are otherwise close to everything and dominate every neighbourhood
- Precompute and truncate the item-item neighbour lists offline; only the top few dozen per item are ever read
Recall@k and coverage on a temporal split, never a random one. Measure catalogue coverage explicitly: a neighbourhood model that scores well by recommending the same head items to everyone is a popularity baseline in disguise.
- Popularity bias compounding into a feedback loop — recommended items gain interactions and become closer neighbours to everything
- Cold items have no interaction vector and therefore no neighbours, so the content fallback is not optional
- Neighbour lists computed once and never refreshed, which silently freezes the catalogue as it was months ago
Computer Vision
AdaptedNot on pixels — on embeddings. A pretrained network maps images to vectors in which semantic similarity is approximately distance, and kNN over those vectors gives retrieval, verification, and few-shot classification without training a classifier at all. Face verification in its deployed form is a nearest-neighbour decision on a learned embedding.
- Reverse image search and visual similarity retrieval
- Face and object verification against an enrolled gallery, where classes appear and disappear constantly
- Few-shot classification from a handful of labelled examples per class
- Near-duplicate detection and content moderation against a known-bad reference set
It is the right structure whenever the set of classes is open or changes constantly: adding a new identity means inserting one vector, not retraining a head. Its quality is entirely inherited from the embedding — kNN adds no representational power at all, which is precisely why it works so well on top of a good encoder and so badly on raw pixels, where Euclidean distance measures brightness rather than content.
- L2-normalize embeddings so cosine and Euclidean distance become equivalent and thresholds transfer
- Use the encoder’s own training metric; a model trained with a cosine margin should not be queried with unnormalized Euclidean distance
- Move to an approximate index past a few hundred thousand vectors, where exact search stops meeting latency
Recall@k and mean average precision for retrieval; for verification, the ROC over distance thresholds, with the operating point chosen from the cost of a false accept relative to a false reject.
- A distance threshold tuned on one embedding version silently invalidated when the encoder is retrained
- Gallery drift — enrolled vectors captured under conditions that no longer resemble query conditions
- Treating the top-1 neighbour as a decision without a distance floor, so an out-of-gallery query always matches somebody
What it costs to run
Zero, and that is not a figure of speech — storing the data is the whole of training. Building a tree or approximate index is the only cost, and it is an indexing cost rather than a fitting one.
The expensive half. A brute-force query is O(nd) and grows with every example ever collected, so latency degrades as the dataset improves. Approximate indexes restore sublinear query time at the cost of a recall guarantee, and that trade — not the model — is where the engineering work goes.
There is nothing to retrain; new data is simply inserted. That is a genuine operational advantage, and it hides a genuine operational risk: nobody schedules a review of a model that never retrains, so drift in the metric or in feature scaling goes unnoticed for a long time.
- Track the distance distribution of queries to their neighbours — rising distances mean queries are landing in regions the reference set does not cover, which is the early warning that nothing else gives you
- Monitor reference-set size against p99 latency, since the two move together by construction
- Watch the feature scaling statistics as a first-class artefact; a scaler change silently rewrites the metric and therefore the model
- For approximate indexes, sample exact recall on a schedule — index quality degrades with insertions and nothing surfaces it otherwise
- The training data IS the model and must ship with it, which makes deletion requests, licensing, and access control model problems rather than data problems
- That memorization is literal: neighbours can be inspected and inverted, so a kNN over personal records leaks those records to anyone who can query it enough times
- Feature scaling must be persisted exactly; refitting a scaler changes what "near" means and every prediction with it
- Class imbalance skews the vote — with a 95/5 split, most neighbourhoods are majority-class regardless of the query, so either weight the vote or rebalance the reference set
- A new categorical level has no sensible distance to anything, so the encoding decides the answer; leaving it to a default one-hot puts every unseen level equidistant from everything
Pros and cons, in context
No training, and new data is incorporated by insertion
Decisive where the label set is open or changes constantly — face galleries, product catalogues, evolving fraud patterns — because adding a class costs one write. Worth nothing where the data is static and inference latency is the binding constraint.
Makes no assumption about the shape of the decision boundary
It will fit an arbitrarily complicated boundary given enough data, which is exactly what a linear model cannot do. The cost is that "enough data" grows exponentially with dimension, so the flexibility is real in low dimensions and illusory in high ones.
Predictions come with their own explanation
The neighbours ARE the justification, which is uniquely defensible in review settings — a reviewer can look at the five cases the answer came from. Most interpretability methods approximate the model; here the explanation is the mechanism.
One mechanism serves classification, regression, density estimation, and outlier scoring
The same neighbour search answers all four, which makes it an efficient thing to build infrastructure around. Less compelling if only one of the four is ever needed, where a purpose-built model will beat it.
The curse of dimensionality destroys it, quietly
Above roughly twenty informative dimensions, nearest and farthest distances converge and the neighbourhood stops being informative. Nothing errors — accuracy just decays toward the base rate, which makes this the most under-diagnosed failure in the method.
Inference cost scales with the size of the training set
Backwards from every other model here: it gets slower as it gets better. Past a few hundred thousand references this forces an approximate index, and with it an accuracy budget that has to be monitored.
The whole training set must be stored and shipped
A memory and a governance problem at once. Under a deletion regime the model itself contains the personal data, so "delete the record" means editing the model, which no other entry in this section requires.
Highly sensitive to scaling, irrelevant features, and class imbalance
An irrelevant feature adds noise to every distance and cannot be down-weighted the way a linear model would zero its coefficient. In practice this means feature selection is a prerequisite here, not an improvement.
Intuitive, idiomatic, optimized
The same algorithm in three languages, each progressing through Kent Beck's three stages. Every optimization names the lever it pulls and the trade-off it accepts — cost profiles are illustrative, not measured benchmarks.
The math, transcribed literally. No dependencies.
"""k-NN classification and regression - the prediction rule, transcribed.
There is no fit step: predict() does all the work, every time. The loop below
is the definition of N_k(x) - measure the distance to every stored point, keep
the k smallest, then reduce their labels under the loss.
"""
import math
def euclidean(a, b):
total = 0.0
for j in range(len(a)):
difference = a[j] - b[j]
total += difference * difference
return math.sqrt(total)
def k_nearest(X_train, query, k):
"""Indices of the k closest training points, nearest first."""
distances = []
for i in range(len(X_train)):
distances.append((euclidean(X_train[i], query), i))
# Full sort of all n distances to keep k of them - correct, and wasteful.
distances.sort()
return [index for _, index in distances[:k]]
def predict_label(X_train, y_train, query, k=5):
"""0-1 loss: the minimizer over the neighbourhood is the majority vote."""
votes = {}
for index in k_nearest(X_train, query, k):
label = y_train[index]
votes[label] = votes.get(label, 0) + 1
best_label = None
best_count = -1
for label, count in votes.items():
if count > best_count:
best_label = label
best_count = count
return best_label
def predict_value(X_train, y_train, query, k=5):
"""Squared loss: the minimizer over the neighbourhood is the mean."""
neighbours = k_nearest(X_train, query, k)
total = 0.0
for index in neighbours:
total += y_train[index]
return total / len(neighbours)O(n*d) distances plus an O(n log n) sort per query, in interpreter loops. The sort is the obvious cost and not the largest one.
How this language actually wants it written.
"""k-NN - typed, validated, vectorized distances, partial selection."""
from dataclasses import dataclass
from typing import Literal
import numpy as np
from numpy.typing import NDArray
Vector = NDArray[np.float64]
Matrix = NDArray[np.float64]
Weighting = Literal["uniform", "distance"]
@dataclass(frozen=True)
class NeighbourIndex:
"""The 'fitted' model: the reference set, plus the metric it is read under.
Scaling lives here rather than in the caller's pipeline because it is part
of the metric, and the metric is the model.
"""
points: Matrix
labels: NDArray
centre: Vector
scale: Vector
k: int
weighting: Weighting = "uniform"
def _neighbours(self, queries: Matrix) -> tuple[NDArray, Matrix]:
scaled = (queries - self.centre) / self.scale
# (q, n) matrix of squared distances - the sqrt is monotone, so it is
# never taken for ranking, only where a weight needs the true distance.
squared = ((scaled[:, None, :] - self.points[None, :, :]) ** 2).sum(axis=2)
# argpartition is O(n) selection, not O(n log n) sorting.
selected = np.argpartition(squared, self.k - 1, axis=1)[:, : self.k]
chosen = np.take_along_axis(squared, selected, axis=1)
return selected, chosen
def predict_label(self, queries: Matrix) -> NDArray:
if queries.ndim != 2 or queries.shape[1] != self.points.shape[1]:
raise ValueError(
f"expected (q, {self.points.shape[1]}) queries, got {queries.shape}"
)
selected, squared = self._neighbours(queries)
weights = self._weights(squared)
classes = np.unique(self.labels)
scores = np.empty((queries.shape[0], classes.size), dtype=np.float64)
for position, label in enumerate(classes):
scores[:, position] = (weights * (self.labels[selected] == label)).sum(axis=1)
return classes[scores.argmax(axis=1)]
def predict_value(self, queries: Matrix) -> Vector:
selected, squared = self._neighbours(queries)
weights = self._weights(squared)
return (weights * self.labels[selected]).sum(axis=1) / weights.sum(axis=1)
def _weights(self, squared: Matrix) -> Matrix:
if self.weighting == "uniform":
return np.ones_like(squared)
# Guard the exact-match case explicitly rather than relying on inf.
return 1.0 / (np.sqrt(squared) + 1e-12)
def build(
X: Matrix,
y: NDArray,
k: int = 5,
weighting: Weighting = "uniform",
) -> NeighbourIndex:
"""Store the reference set. Raises ValueError on malformed input."""
if X.ndim != 2:
raise ValueError(f"X must be 2-D, got shape {X.shape}")
if X.shape[0] != y.shape[0]:
raise ValueError(f"X has {X.shape[0]} rows but y has {y.shape[0]}")
if not 1 <= k <= X.shape[0]:
raise ValueError(f"k must lie in [1, {X.shape[0]}], got {k}")
centre = X.mean(axis=0)
scale = X.std(axis=0)
scale[scale == 0.0] = 1.0
return NeighbourIndex(
points=(X - centre) / scale,
labels=y,
centre=centre,
scale=scale,
k=k,
weighting=weighting,
)What changed: Two changes that matter more than the vectorization. The full sort becomes argpartition — selecting the k smallest is an O(n) problem and sorting all n to answer it was the algorithmic mistake in the previous stage. And the square root disappears from ranking entirely, because it is monotone: comparing squared distances gives the same ordering for a fraction of the cost. Structurally, scaling moves inside the index where it belongs (it is part of the metric, not preprocessing), queries are answered in batches rather than one at a time, and malformed input fails at the boundary instead of deep inside a loop.
O(q*n*d) distances plus O(q*n) selection, executed in NumPy rather than the interpreter.
Faster, and the justification for each change.
"""k-NN - squared-norm expansion, one GEMM, chunked over queries."""
import numpy as np
from numpy.typing import NDArray
Vector = NDArray[np.float64]
Matrix = NDArray[np.float64]
class BlockedNeighbours:
"""Exact k-NN, with the distance matrix computed as a matrix product.
||q - x||^2 = ||q||^2 - 2 q . x + ||x||^2
The cross term is the only part that depends on both operands, and it is a
single GEMM. The two squared-norm terms are precomputed once per side, so
the per-query cost collapses onto BLAS instead of onto a broadcast
subtraction that materializes a (q, n, d) intermediate.
"""
def __init__(self, points: Matrix, labels: NDArray, k: int = 5) -> None:
if not 1 <= k <= points.shape[0]:
raise ValueError(f"k must lie in [1, {points.shape[0]}], got {k}")
# One contiguous float64 block: BLAS reads it without an internal copy.
self._points = np.ascontiguousarray(points, dtype=np.float64)
self._labels = labels
self._k = k
self._point_norms = np.einsum("ij,ij->i", self._points, self._points)
def query(self, queries: Matrix, chunk: int = 4096) -> tuple[NDArray, Matrix]:
"""Returns (neighbour indices, squared distances) for every query."""
queries = np.ascontiguousarray(queries, dtype=np.float64)
n_queries = queries.shape[0]
indices = np.empty((n_queries, self._k), dtype=np.int64) # allocated once
distances = np.empty((n_queries, self._k), dtype=np.float64)
# Chunking bounds peak memory at (chunk, n) rather than (q, n), which is
# what makes a million-query batch possible at all.
for start in range(0, n_queries, chunk):
stop = min(start + chunk, n_queries)
block = queries[start:stop]
squared = self._points @ block.T # the GEMM
squared *= -2.0 # in place
squared += self._point_norms[:, None]
# ||q||^2 is constant within a column and does not affect ranking,
# so it is added only because callers read the distances back.
squared += np.einsum("ij,ij->i", block, block)[None, :]
selected = np.argpartition(squared, self._k - 1, axis=0)[: self._k]
chosen = np.take_along_axis(squared, selected, axis=0)
order = np.argsort(chosen, axis=0) # k log k, not n log n
indices[start:stop] = np.take_along_axis(selected, order, axis=0).T
distances[start:stop] = np.take_along_axis(chosen, order, axis=0).T
return indices, distances
def predict_value(self, queries: Matrix) -> Vector:
indices, _ = self.query(queries)
return self._labels[indices].mean(axis=1)What changed: The distance computation stops being a broadcast subtraction and becomes a matrix product, by expanding the square: the only term involving both operands is the inner product, which is a GEMM, and the two norm terms are precomputed once per side. That removes the (q, n, d) intermediate the previous stage materialized, which was the real memory ceiling rather than the arithmetic. Queries are then processed in chunks so peak memory is bounded independently of batch size. The remaining honesty caveat is that this makes exact brute force fast without making it sublinear — past a few hundred thousand references the answer is an approximate index, not a better GEMM.
Vectorize to NumPy/BLAS (@, np.dot, einsum)
The cross term becomes one GEMM, and the squared norms are einsum reductions computed once per side rather than per pair.
Trade-off: The expansion is less numerically stable than a direct subtraction — cancellation between large squared norms can produce small negative distances that must be clipped before any sqrt.
Pre-allocate output arrays and use in-place operations
Index and distance outputs are allocated once for the whole batch, and the distance block is built with in-place multiply and add rather than three temporaries per chunk.
Trade-off: The in-place chain is order-dependent and no longer reads as the formula it implements, so a reordered line changes the result silently instead of failing.
Batch work to amortize interpreter overhead
Chunking over queries turns per-query Python dispatch into one BLAS call per few thousand queries, and bounds peak memory at (chunk, n) instead of (q, n).
Trade-off: Chunk size is a tuning parameter with no good default — too small and dispatch dominates, too large and the distance block no longer fits in cache.
O(q*n*d) in BLAS plus O(q*n) selection, memory bounded by the chunk. Illustrative, not a measured benchmark.