Every two hours our cron fetches trending videos across eight regions and dumps them into a single SQLite table. By the time a visitor lands on the homepage, that table holds tens of thousands of candidate videos, each with a freshness score, a regional relevance score, and an engagement estimate. The homepage needs to show exactly 24 of them. The naive answer is ORDER BY score DESC LIMIT 24 — and it is exactly what we shipped first, and exactly what made the site feel dead within a week.
The problem with sorting by score is that the top of the list never moves. The same viral clip sat at position one for three days because its score dominated everything else. Users who came back twice saw an identical grid. What I actually wanted was a feed that favours high-scoring videos without being deterministic — pick 24 videos where a video with twice the weight is roughly twice as likely to appear, but where every refresh produces a fresh, plausible selection. That is weighted random sampling without replacement, and doing it in a single streaming pass over a SQLite cursor is what weighted reservoir sampling gives you. This post is the algorithm we run in production at TrendVidStream, the pitfalls that cost me a weekend, and runnable code in PHP 8.4, Python, and Go.
Why not just weighted-shuffle the whole table
The textbook way to draw k weighted samples without replacement is a weighted shuffle: assign each item a sort key and take the top k. That works, but it forces you to materialise and sort every candidate. When the candidate set is 40,000 rows and you want 24, sorting 40,000 rows to throw away 39,976 of them is wasteful — and it means holding the full result set in memory, which on a shared LiteSpeed host with a 512 MB PHP limit is not free.
Reservoir sampling flips the cost model. You walk the stream once, keep a fixed-size reservoir of k items, and at the end the reservoir is your sample. Memory is O(k), not O(n). You never need to know the total count in advance, which matters because our candidate query has a bunch of WHERE filters (region, safe-for-ads, not-already-watched) whose result size changes every fetch cycle. A streaming algorithm doesn't care.
The classic reservoir algorithm (Vitter's Algorithm R) samples uniformly. We need weighted sampling, and that is where Efraimidis and Spirakis come in.
The A-Res algorithm in one paragraph
Efraimidis and Spirakis (2006) proved something beautifully simple. For each item i with weight w_i, draw a uniform random number u_i in (0, 1) and compute a key:
key_i = u_i ^ (1 / w_i)
Keep the k items with the largest keys. That's it. The set of k largest keys is a mathematically correct weighted sample without replacement, where the probability of an item being included is proportional to its weight in exactly the sense you want. This is called A-Res (Algorithm for Reservoir sampling). Because you only ever need the k largest keys, you maintain a min-heap of size k: the heap root is the smallest key currently in the reservoir, and any new item whose key beats the root evicts it.
Two things are worth internalising before you write code:
- The exponent is
1/w_i. Larger weight → the keyu^(1/w)is pushed closer to 1 → more likely to be among the largest. If you accidentally writeu^wyou get the exact inverse behaviour and your "recommendations" will surface your worst videos. Ask me how I know. - The keys are wildly non-uniform in magnitude. With
win the thousands,1/wis tiny andu^(1/w)is extremely close to 1.0, sofloatprecision matters. Use 64-bit doubles and work in log space if your weights span many orders of magnitude (covered below).
Production implementation in PHP 8.4
Here is the version that runs on our homepage. It streams rows from a PDO statement, keeps a min-heap of size k using PHP's SplMinHeap, and returns the selected video IDs. PHP 8.4's typed class constants and readonly properties keep it tidy.
<?php
declare(strict_types=1);
final class WeightedReservoir
{
private const float MIN_WEIGHT = 1e-9;
/** @var \SplMinHeap<array{float, array}> */
private \SplMinHeap $heap;
public function __construct(private readonly int $k)
{
if ($k < 1) {
throw new \InvalidArgumentException('k must be >= 1');
}
// Order pairs by their key (element 0) so the root is the smallest key.
$this->heap = new class extends \SplMinHeap {
protected function compare(mixed $a, mixed $b): int
{
return $b[0] <=> $a[0];
}
};
}
public function offer(array $row, float $weight): void
{
$weight = max($weight, self::MIN_WEIGHT);
// u in (0,1); avoid exactly 0 which would give key 0 for any weight.
$u = (mt_rand(1, mt_getrandmax()) / (mt_getrandmax() + 1));
$key = $u ** (1.0 / $weight);
if ($this->heap->count() < $this->k) {
$this->heap->insert([$key, $row]);
return;
}
// Peek the smallest key currently retained.
if ($key > $this->heap->top()[0]) {
$this->heap->extract();
$this->heap->insert([$key, $row]);
}
}
/** @return array<int, array> selected rows, highest key first */
public function result(): array
{
$out = [];
foreach (clone $this->heap as $pair) {
$out[] = $pair[1];
}
return array_reverse($out);
}
}
And the calling code that wires it to SQLite. Note we never load the full result set — we iterate the statement lazily, so peak memory is bounded by k plus one row:
<?php
$pdo = new PDO('sqlite:' . __DIR__ . '/data/videos.db');
$pdo->setAttribute(PDO::ATTR_ERRMODE, PDO::ERRMODE_EXCEPTION);
$stmt = $pdo->query(<<<SQL
SELECT id, title, region, freshness, relevance, engagement
FROM videos
WHERE safe_for_ads = 1
AND published_at > strftime('%s','now','-14 days')
SQL);
$reservoir = new WeightedReservoir(k: 24);
while ($row = $stmt->fetch(PDO::FETCH_ASSOC)) {
// Composite weight: freshness dominates, relevance and engagement modulate.
$weight = ($row['freshness'] * 3.0)
+ ($row['relevance'] * 2.0)
+ ($row['engagement'] * 1.0);
$reservoir->offer($row, $weight);
}
$homepage = $reservoir->result();
printf("Selected %d videos from the stream\n", count($homepage));
This runs in a few milliseconds for 40k rows because the per-row work is a random draw, a pow, and at most one heap operation. The heap only ever holds 24 entries. On our slowest host the whole selection is dwarfed by the query itself.
A subtlety about randomness: mt_rand is fine for feed diversity — this is not cryptography. If you were sampling for anything security- or fairness-sensitive you would reach for random_int, but that is slower and unnecessary for shuffling a video wall.
The numeric stability trap
The first time I deployed this, videos with very large engagement counts (millions of views) never got dropped from the reservoir, and videos with modest weights never got in. It looked correct — big videos win — but the distribution was wrong: it had collapsed toward deterministic top-k again.
The cause is floating point. When w is 5,000,000, 1/w is 2e-7, and u ** 2e-7 rounds to 1.0 for almost any u in a 64-bit double. Every heavyweight item gets a key of exactly 1.0, ties everywhere, and the randomness disappears. The fix is to work in log space. Instead of comparing u^(1/w), compare its logarithm:
log(key) = (1 / w) * log(u)
Since log(u) is negative and 1/w is positive, log(key) is negative, and larger (closer to zero) log-keys correspond to larger keys — the ordering is preserved and the min-heap logic is identical. Here is the corrected offer in PHP:
public function offer(array $row, float $weight): void
{
$weight = max($weight, self::MIN_WEIGHT);
$u = (mt_rand(1, mt_getrandmax()) / (mt_getrandmax() + 1));
// log-key preserves ordering and survives huge weights without collapsing to 1.0
$logKey = log($u) / $weight;
if ($this->heap->count() < $this->k) {
$this->heap->insert([$logKey, $row]);
return;
}
if ($logKey > $this->heap->top()[0]) {
$this->heap->extract();
$this->heap->insert([$logKey, $row]);
}
}
Because log(u) for u near 1 approaches 0 from below, and for small u goes strongly negative, the divisor w spreads the values out cleanly across the whole range regardless of magnitude. This one change fixed the distribution and it costs nothing — log is as cheap as pow. I now reach for the log-space form by default and never look back.
A-ExpJ: skipping over the boring rows
A-Res touches every item. That is fine at 40k rows. But we also run an offline job that scores a much larger historical corpus — millions of rows — to build per-region "deep catalogue" pools. There, doing a random draw and a heap comparison for every single row is measurable overhead.
Efraimidis and Spirakis's second algorithm, A-ExpJ, computes an exponential jump: after the reservoir fills, it calculates how much total weight it can skip before the next item could possibly enter, then jumps ahead consuming that weight without generating a key per row. For skewed weight distributions this cuts the number of random draws by an order of magnitude. Here is a Python implementation of A-ExpJ for the offline pool builder:
import heapq
import math
import random
from typing import Iterator, Iterable
def weighted_reservoir_expj(
stream: Iterable[tuple[dict, float]], k: int
) -> list[dict]:
"""A-ExpJ weighted reservoir sampling without replacement.
stream yields (item, weight) pairs. Returns k items."""
heap: list[tuple[float, int, dict]] = [] # (log_key, tiebreak, item)
tiebreak = 0
it: Iterator[tuple[dict, float]] = iter(stream)
# Phase 1: fill the reservoir with the first k items.
for item, weight in it:
weight = max(weight, 1e-9)
log_key = math.log(random.random()) / weight
heapq.heappush(heap, (log_key, tiebreak, item))
tiebreak += 1
if len(heap) == k:
break
if len(heap) < k:
return [row for _, _, row in heap]
# Phase 2: exponential jumps. X is the weight budget until the next entry.
threshold = heap[0][0] # smallest log-key in the reservoir
x_weight = math.log(random.random()) / threshold
for item, weight in it:
weight = max(weight, 1e-9)
x_weight -= weight
if x_weight > 0:
continue # skip: this item cannot beat the current threshold
# This item enters. Draw its key within the valid sub-interval.
t = math.exp(threshold * weight)
r = random.uniform(t, 1.0)
log_key = math.log(r) / weight
heapq.heapreplace(heap, (log_key, tiebreak, item))
tiebreak += 1
threshold = heap[0][0]
x_weight = math.log(random.random()) / threshold
return [row for _, _, row in heap]
if __name__ == "__main__":
# Smoke test: item weighted 10x should appear ~10x more often.
counts = {"heavy": 0, "light": 0}
trials = 20000
for _ in range(trials):
pool = [({"id": "heavy"}, 10.0)] + [({"id": "light"}, 1.0)] * 10
random.shuffle(pool)
picked = weighted_reservoir_expj(iter(pool), k=1)[0]
counts[picked["id"]] += 1
# heavy weight 10 vs total light weight 10 -> ~50/50
print(counts)
The tiebreak counter matters in Python: heapq compares tuples element by element, and if two log_key values are ever equal it would try to compare the dict items and raise TypeError. An always-increasing integer as the second element makes ties resolve deterministically and keeps heapq from ever touching the payload. It's a one-line defensive habit that saves a confusing crash in production.
Validating that the distribution is actually correct
Weighted sampling bugs are insidious because the output always looks reasonable — high-weight items show up a lot, which is what you expected, so you don't notice that the ratio is off. The only way to trust the implementation is a statistical test. The smoke test in the Python snippet above encodes the core invariant: a single item of weight 10 competing against ten items of weight 1 each (total light weight 10) should win roughly half the time. If you see 70/30 instead of 50/50, your exponent is wrong.
For a more thorough check I run a chi-squared style comparison: draw k=1 many thousands of times from a fixed weighted pool and confirm each item's empirical selection frequency matches w_i / sum(w) within sampling error. Do this in CI. It is the single highest-leverage test for this kind of code, and it is cheap.
A concurrent Go version for the ingest path
We also merge sampled feeds from several regional shards, and Go's concurrency makes it natural to sample each shard in its own goroutine and merge the reservoirs. Because two independent weighted reservoirs can be merged by simply re-running the selection over the union of their retained items with their original keys, you get a clean map-reduce. Here is a self-contained Go implementation using container/heap:
package main
import (
"container/heap"
"fmt"
"math"
"math/rand"
)
type entry struct {
logKey float64
item string
}
// minHeap keeps the smallest logKey at the root.
type minHeap []entry
func (h minHeap) Len() int { return len(h) }
func (h minHeap) Less(i, j int) bool { return h[i].logKey < h[j].logKey }
func (h minHeap) Swap(i, j int) { h[i], h[j] = h[j], h[i] }
func (h *minHeap) Push(x any) { *h = append(*h, x.(entry)) }
func (h *minHeap) Pop() any {
old := *h
n := len(old)
x := old[n-1]
*h = old[:n-1]
return x
}
type Reservoir struct {
k int
h minHeap
r *rand.Rand
}
func NewReservoir(k int, seed int64) *Reservoir {
return &Reservoir{k: k, h: make(minHeap, 0, k), r: rand.New(rand.NewSource(seed))}
}
func (res *Reservoir) Offer(item string, weight float64) {
if weight < 1e-9 {
weight = 1e-9
}
u := res.r.Float64()
if u == 0 {
u = 1e-12
}
logKey := math.Log(u) / weight
if len(res.h) < res.k {
heap.Push(&res.h, entry{logKey, item})
return
}
if logKey > res.h[0].logKey {
res.h[0] = entry{logKey, item}
heap.Fix(&res.h, 0)
}
}
func (res *Reservoir) Result() []string {
out := make([]string, len(res.h))
for i, e := range res.h {
out[i] = e.item
}
return out
}
func main() {
res := NewReservoir(3, 42)
weights := map[string]float64{"a": 100, "b": 50, "c": 5, "d": 5, "e": 1}
for id, w := range weights {
res.Offer(id, w)
}
fmt.Println(res.Result())
}
Swapping heap.Fix in for a pop-then-push saves a reallocation on the eviction path — a small but honest win when this runs per row across millions of items. For the multi-shard merge, each goroutine keeps its own Reservoir, and a final pass offers every retained (item, weight) into one merge reservoir. The math holds because A-Res keys are independent per item.
How this plugs into the wider system
A few operational notes from running this across eight regions:
-
Regional weighting is just a weight multiplier. For a visitor in Japan we multiply the
relevanceterm for JP-tagged videos, keep the same reservoir code, and the feed naturally tilts local without a separate query per region. -
Freshness decay lives in the weight, not the query. We compute
exp(-age_hours / half_life)as part of the composite weight. Because the sampler is stateless per request, tuning the half-life is a config change, not a schema migration. - The candidate query is FTS5-friendly. When a search term is present we sample from the FTS5 match set instead of the whole table — same sampler, different stream source. Streaming design means the algorithm doesn't care where the rows come from.
- It is deploy-safe. The sampler is a single self-contained class with no dependencies, so it rides along in our FTP-based deploy with everything else and needs no build step on the host.
Conclusion
Weighted reservoir sampling turned a stale, deterministic video wall into a feed that respects scoring signals while staying fresh on every refresh — and it did so with O(k) memory in a single streaming pass, which is exactly what you want on shared hosting with tight memory limits. The whole thing is three ideas: the u^(1/w) key, a min-heap of size k, and — the part that actually bit me — doing it all in log space so large weights don't collapse to 1.0. If you take one thing away, take the statistical test: weighted sampling bugs hide in plain sight because wrong output still looks plausible, and only a frequency check will tell you the ratios are right. Start with the A-Res version, add A-ExpJ only when your streams get big enough to feel the per-row cost, and validate the distribution in CI before you trust it in front of users.
Top comments (0)