Compression in Practice:
Arithmetic Coding, LZ & Beyond

Every compressor is a bet on which of these patterns your data has. When the bet is wrong, the output gets slightly bigger.

import random, zlib

n = 16
strings_of_n_bits = 2 ** n
shorter_strings = sum(2 ** k for k in range(n))   # lengths 0 .. n-1
print(f"{strings_of_n_bits} inputs, only {shorter_strings} shorter outputs")

blob = random.Random(0).randbytes(100_000)
print("random bytes:", len(blob), "->", len(zlib.compress(blob, 9)))
text = b"to be or not to be " * 5000
print("repetitive  :", len(text), "->", len(zlib.compress(text, 9)))

65536 inputs, only 65535 shorter outputs
random bytes: 100000 -> 100041
repetitive  : 95000 -> 271

• Random bytes grow by 41 bytes. That is zlib's header, checksum and "stored block" markers.
• Repetitive text shrinks 350×.
• Good formats cap the loss. DEFLATE falls back to stored blocks, costing about 5 bytes per 64 KiB.

This is why "compress it twice" never helps, and why compressed or encrypted files do not shrink.

The split also lets you swap parts. zstd uses Huffman for literals and tANS for everything else. LZMA uses a range coder with bit-level context models.

Test corpus used on every later slide. It is seeded, so it is the same on each run:

import random
random.seed(42)
WORDS = ("the of and to a in is that for it as was with be by on not he this are "
         "or his from at which but have an they you were her she there been one all "
         "we their has would when if so no will more can who its said about into "
         "data code bits model symbol entropy compress window match stream").split()
weights = [1 / (r + 1) for r in range(len(WORDS))]       # Zipf-like word frequencies
sentences = []
for _ in range(400):
    n = random.randint(5, 14)
    s = " ".join(random.choices(WORDS, weights, k=n))
    sentences.append(s.capitalize() + ".")
TEXT = " ".join(sentences).encode()
print(len(TEXT), "bytes:", TEXT[:60].decode())

15019 bytes: The of of he be we. For the the of it that. To a of there ar

Zipf-weighted words: a few are very common, most are rare. That matches real text, but there are no long repeated phrases.

Empirical conditional entropy H(X_t | previous k bytes):

from collections import Counter, defaultdict
from math import log2

def cond_entropy(data, k):      # bits/byte, given k previous bytes
    ctx = defaultdict(Counter)
    for i in range(k, len(data)):
        ctx[data[i - k:i]][data[i]] += 1
    n = len(data) - k
    return -sum(c * log2(c / sum(cnt.values()))
                for cnt in ctx.values() for c in cnt.values()) / n

for k in range(4):
    print(f"order-{k} model: {cond_entropy(TEXT, k):.3f} bits/byte")

order-0 model: 3.799 bits/byte
order-1 model: 2.219 bits/byte
order-2 model: 1.457 bits/byte
order-3 model: 1.286 bits/byte

Canonical codes are fully set by the code lengths. DEFLATE sends only those lengths. They meet Kraft's inequality with equality: Σ 2^−L = 1. Order-0 entropy here is 3.7995, so Huffman loses 0.03 bits/byte.

import heapq
from collections import Counter

def huffman_lengths(data):
    heap = [(n, i, {s: 0}) for i, (s, n) in enumerate(Counter(data).items())]
    heapq.heapify(heap)
    tie = len(heap)
    while len(heap) > 1:                   # merge the two rarest subtrees
        n1, _, a = heapq.heappop(heap)
        n2, _, b = heapq.heappop(heap)
        merged = {s: L + 1 for s, L in {**a, **b}.items()}
        heapq.heappush(heap, (n1 + n2, tie, merged)); tie += 1
    return heap[0][2]

def canonical(lengths):                    # DEFLATE sends only the lengths
    code, prev, out = 0, 0, {}
    for s, L in sorted(lengths.items(), key=lambda kv: (kv[1], kv[0])):
        code <<= L - prev
        out[s] = format(code, f"0{L}b"); code += 1; prev = L
    return out

L = huffman_lengths(TEXT)
codes = canonical(L)
for s in b"e t.":
    print(repr(chr(s)), codes[s])
bits = sum(L[b] for b in TEXT)
print(f"Huffman: {bits / len(TEXT):.4f} bits/byte")

'e' 1001
' ' 00
't' 011
'.' 10110
Huffman: 3.8289 bits/byte

Six letters in a 100-letter file (the classic CLRS example). Rule: take the two smallest counts, join them, repeat.

Now walk down from the root. Left is 0, right is 1. A letter's code is its path.

Fractions of a bit per symbol are fine. They add up across the whole message.

The decoder gets the number 0.625 and the same model: A = [0, .5), B = [.5, .8), C = [.8, 1).

Each step: rescale the number into the current interval, t = (x − low) / (high − low). See which slot t falls in. Output that symbol and shrink.

Steps 1–3 give back BAC, the message from the last slide.

Overhead is 2 bits per message, not per symbol. Huffman can lose up to 1 bit per symbol.

This is the Witten–Neal–Cleary coder (CACM 1987). Rissanen and Pasco found the idea in 1976. IBM's patents slowed its use for about 20 years.

FULL, HALF, QTR = (1 << 32) - 1, 1 << 31, 1 << 30

class Encoder:
    def __init__(self):
        self.low, self.high, self.pending, self.bits = 0, FULL, 0, []

    def emit(self, b):
        self.bits.append(b)
        self.bits.extend([1 - b] * self.pending)   # flush straddle bits
        self.pending = 0

    def encode(self, lo, hi, total):               # symbol owns [lo, hi) of total
        rng = self.high - self.low + 1
        self.high = self.low + rng * hi // total - 1
        self.low = self.low + rng * lo // total
        while True:
            if self.high < HALF:
                self.emit(0)
            elif self.low >= HALF:
                self.emit(1); self.low -= HALF; self.high -= HALF
            elif self.low >= QTR and self.high < HALF + QTR:
                self.pending += 1; self.low -= QTR; self.high -= QTR
            else:
                break
            self.low, self.high = 2 * self.low, 2 * self.high + 1

    def finish(self):
        self.pending += 1
        self.emit(0 if self.low < QTR else 1)
        return self.bits

About 30 lines. The same design, with 64-bit words and byte output, sits inside LZMA, CABAC (H.264/HEVC) and many PAQ-family compressors.

class Decoder:
    def __init__(self, bits):
        self.bits, self.pos = bits, 0
        self.low, self.high, self.value = 0, FULL, 0
        for _ in range(32):
            self.value = 2 * self.value + self.next_bit()

    def next_bit(self):
        b = self.bits[self.pos] if self.pos < len(self.bits) else 0
        self.pos += 1
        return b

    def target(self, total):                       # where does value fall?
        rng = self.high - self.low + 1
        return ((self.value - self.low + 1) * total - 1) // rng

    def consume(self, lo, hi, total):              # mirror the encoder exactly
        rng = self.high - self.low + 1
        self.high = self.low + rng * hi // total - 1
        self.low = self.low + rng * lo // total
        while True:
            if self.high < HALF:
                pass
            elif self.low >= HALF:
                self.low -= HALF; self.high -= HALF; self.value -= HALF
            elif self.low >= QTR and self.high < HALF + QTR:
                self.low -= QTR; self.high -= QTR; self.value -= QTR
            else:
                break
            self.low, self.high = 2 * self.low, 2 * self.high + 1
            self.value = 2 * self.value + self.next_bit()

Reading past the end returns 0 bits. The encoder's two final bits make any padding decode the same.

from collections import Counter
from math import log2

counts = Counter(TEXT)                        # static order-0 model
syms = sorted(counts)
cum, c = {}, 0
for s in syms:
    cum[s] = (c, c + counts[s]); c += counts[s]
total = c

enc = Encoder()
for b in TEXT:
    enc.encode(*cum[b], total)
bits = enc.finish()

dec = Decoder(bits)
out = bytearray()
for _ in range(len(TEXT)):
    t = dec.target(total)
    s = next(s for s in syms if cum[s][0] <= t < cum[s][1])
    dec.consume(*cum[s], total)
    out.append(s)

H = -sum(n / len(TEXT) * log2(n / len(TEXT)) for n in counts.values())
print("roundtrip ok:", bytes(out) == TEXT)
print(f"entropy bound : {H * len(TEXT):.0f} bits")
print(f"arith. coded  : {len(bits)} bits ({len(bits) / len(TEXT):.4f} bits/byte)")

roundtrip ok: True
entropy bound : 57063 bits
arith. coded  : 57064 bits (3.7995 bits/byte)

class AdaptiveModel:
    """Order-k context model: counts start at 1 and grow as we code."""
    def __init__(self, k):
        self.k, self.tables = k, {}

    def freqs(self, ctx):
        return self.tables.setdefault(ctx[-self.k:] if self.k else b"", [1] * 256)

    def update(self, ctx, sym):
        self.freqs(ctx)[sym] += 32              # big step = learn fast

def interval(f, s):
    lo = sum(f[:s])
    return lo, lo + f[s], sum(f)

def compress(data, k):
    enc, m = Encoder(), AdaptiveModel(k)
    for i, b in enumerate(data):
        enc.encode(*interval(m.freqs(data[:i]), b))
        m.update(data[:i], b)
    return enc.finish()

Uses the Encoder/Decoder from the two slides before. The O(256) scan per symbol is for clarity. Real coders use Fenwick trees or binary decomposition.

def decompress(bits, n, k):
    dec, m, out = Decoder(bits), AdaptiveModel(k), bytearray()
    for _ in range(n):
        f = m.freqs(bytes(out))
        t, acc, s = dec.target(sum(f)), 0, 0
        while acc + f[s] <= t:
            acc += f[s]; s += 1
        dec.consume(*interval(f, s))
        m.update(bytes(out), s); out.append(s)
    return bytes(out)

for k in range(3):
    bits = compress(TEXT, k)
    assert decompress(bits, len(TEXT), k) == TEXT
    print(f"order-{k} adaptive: {len(bits) / len(TEXT):.3f} bits/byte")

order-0 adaptive: 3.822 bits/byte
order-1 adaptive: 2.389 bits/byte
order-2 adaptive: 1.970 bits/byte

• LIFO: encode in reverse, decode forward.
• Real coders keep x in a machine word and renormalize by streaming out bytes.
• tANS precomputes a table (FSE in zstd, LZFSE). Hundreds of MB/s, no multiplies.

from collections import Counter

def build(data, M_BITS=12):               # quantize counts to sum M = 4096
    M, cnt = 1 << M_BITS, Counter(data)
    freq = {s: max(1, c * M // len(data)) for s, c in cnt.items()}
    top = max(freq, key=freq.get)
    freq[top] += M - sum(freq.values())    # fix rounding so the sum is exactly M
    cum, c = {}, 0
    for s in sorted(freq):
        cum[s] = c; c += freq[s]
    return M, freq, cum

def rans_encode(data, M, freq, cum):
    x = 1
    for s in reversed(data):              # ANS is last-in, first-out
        f = freq[s]
        x = (x // f) * M + cum[s] + x % f
    return x

def rans_decode(x, n, M, freq, cum):
    slot_to_sym = [s for s in sorted(freq) for _ in range(freq[s])]
    out = bytearray()
    for _ in range(n):
        s = slot_to_sym[x % M]
        x = freq[s] * (x // M) + x % M - cum[s]
        out.append(s)
    return bytes(out)

M, freq, cum = build(TEXT)
x = rans_encode(TEXT, M, freq, cum)
assert rans_decode(x, len(TEXT), M, freq, cum) == TEXT
print(f"rANS state: {x.bit_length()} bits "
      f"({x.bit_length() / len(TEXT):.4f} bits/byte)")

rANS state: 57076 bits (3.8003 bits/byte)

Python big integers let us skip renormalization. 3.8003 vs 3.7995 bits/byte: the tiny gap comes from rounding the frequencies to 12 bits. JPEG XL, zstd and Apple's LZFSE all use ANS.

Range coding (Martin, 1979) is arithmetic coding that outputs whole bytes. It is faster, and its patent status was clearer. LZMA/xz uses it.

def rle_encode(data):
    out, i = [], 0
    while i < len(data):
        j = i
        while j < len(data) and data[j] == data[i] and j - i < 255:
            j += 1
        out.append((j - i, data[i]))
        i = j
    return out

def rle_decode(pairs):
    return b"".join(bytes([v]) * n for n, v in pairs)

row = b"\x00" * 40 + b"\xff" * 12 + b"\x00" * 76     # one scanline of a fax
pairs = rle_encode(row)
print(pairs)
assert rle_decode(pairs) == row
print(len(row), "bytes ->", 2 * len(pairs), "bytes")

[(40, 0), (12, 255), (76, 0)]
128 bytes -> 6 bytes

The < 255 cap keeps each count in one byte. Longer runs just split.

RLE is LZ77 limited to distance 1. Generalize the distance and you get LZ77.

def lz77_encode(data, window=4096, min_len=3):
    i, out = 0, []
    while i < len(data):
        best_len = best_dist = 0
        for j in range(max(0, i - window), i):     # naive O(n * window) search
            L = 0
            while i + L < len(data) and data[j + L] == data[i + L]:
                L += 1                             # may run past i: overlap is OK
            if L > best_len:
                best_len, best_dist = L, i - j
        if best_len >= min_len:
            out.append((best_dist, best_len)); i += best_len
        else:
            out.append(data[i]); i += 1
    return out

['a', 'b', 'r', 'a', 'c', 'a', 'd', (7, 4), ' ', (12, 13), (1, 9)]

Decoding is just copying, so it runs at GB/s. The encoder does all the searching.

def lz77_decode(tokens):
    out = bytearray()
    for t in tokens:
        if isinstance(t, tuple):
            dist, length = t
            for _ in range(length):     # byte by byte: overlap OK
                out.append(out[-dist])
        else:
            out.append(t)
    return bytes(out)

s = b"abracadabra abracadabra aaaaaaaaaa"
toks = lz77_encode(s)
print([chr(t) if isinstance(t, int) else t for t in toks])
assert lz77_decode(toks) == s

Encode abababcaba (min match 3)

Output: ['a', 'b', (2, 4), 'c', (7, 3)]. Ten bytes became five tokens.

More traces from the deck's encoder

def lzw_encode(data):
    table = {bytes([i]): i for i in range(256)}
    w, out = b"", []
    for b in data:
        wc = w + bytes([b])
        if wc in table:
            w = wc
        else:
            out.append(table[w])
            table[wc] = len(table)              # learn a new phrase
            w = bytes([b])
    if w:
        out.append(table[w])
    return out

def lzw_decode(codes):
    table = {i: bytes([i]) for i in range(256)}
    w = table[codes[0]]
    out = [w]
    for c in codes[1:]:
        entry = table[c] if c in table else w + w[:1]   # the "cScSc" case
        out.append(entry)
        table[len(table)] = w + entry[:1]
        w = entry
    return b"".join(out)

s = b"TOBEORNOTTOBEORTOBEORNOT"
codes = lzw_encode(s)
print(codes)
assert lzw_decode(codes) == s
big = lzw_encode(TEXT)
print(f"TEXT: {len(TEXT)} bytes -> {len(big)} codes "
      f"(~{len(big) * 12 / 8:.0f} bytes at 12 bits/code)")

[84, 79, 66, 69, 79, 82, 78, 79, 84, 256, 258, 260, 265, 259, 261, 263]
TEXT: 15019 bytes -> 3520 codes (~5280 bytes at 12 bits/code)

def bwt(s):
    s += b"\x00"                                   # unique end marker
    rots = sorted(range(len(s)), key=lambda i: s[i:] + s[:i])
    return bytes(s[i - 1] for i in rots)

def ibwt(last):
    first = sorted(range(len(last)), key=lambda i: (last[i], i))  # stable LF map
    i, out = first[0], bytearray()
    for _ in range(len(last) - 1):
        i = first[i]
        out.append(last[i])
    return bytes(out)

def mtf(data):
    table, out = list(range(256)), []
    for b in data:
        k = table.index(b)
        out.append(k)
        table.insert(0, table.pop(k))              # move to front
    return out

s = b"banana_bandana_banana"
L = bwt(s)
print(L)
print(mtf(L))
assert ibwt(L) == s

b'aaannnndnbbb_\x00_naaaaaa'
[97, 0, 0, 110, 0, 0, 0, 101, 1, 100, 0, 0, 99, 5, 1, 3, 5, 0, 0, 0, 0, 0]

Real BWT uses suffix arrays (SA-IS, O(n)). The same trick powers the FM-index for DNA alignment (BWA, Bowtie).

Add an end marker $ that sorts first. Write all 7 rotations of banana$. Sort them.

Output: L = annb$aa. First column F = $aaabnn is just L sorted, so it is free.

See the grouping? Rows 1–3 all start with "a", and two of them end in "n". "n" comes before "a" in banana. On real text these runs get long.

The price: CM runs at about 1–100 KB/s, and decoding is as slow as encoding. Use it for archives, not for the web.

import bz2, lzma, zlib, random

samples = {
    "text (TEXT)": TEXT,
    "random": random.Random(1).randbytes(15000),
    "zeros": bytes(15000),
}
codecs = {
    "zlib -9": lambda d: zlib.compress(d, 9),
    "bz2 -9": lambda d: bz2.compress(d, 9),
    "lzma": lambda d: lzma.compress(d, preset=9),
}
print(f"{'':12}" + "".join(f"{c:>9}" for c in codecs))
for name, d in samples.items():
    print(f"{name:12}" + "".join(f"{len(f(d)):>9}" for f in codecs.values()))

              zlib -9   bz2 -9     lzma
text (TEXT)      4616     3450     4364
random          15011    15481    15060
zeros              37       46      116

Always benchmark on your data. Rankings flip between text, logs, binaries and images.

import zlib

msgs = [f'{{"user": {i}, "event": "click", "page": "/home", "ok": true}}'.encode()
        for i in range(100)]
dictionary = b'{"user": , "event": "click", "view", "page": "/home", "ok": true}'

plain = sum(len(zlib.compress(m, 9)) for m in msgs)

def with_dict(m):
    c = zlib.compressobj(9, zdict=dictionary)
    return c.compress(m) + c.flush()

primed = sum(len(with_dict(m)) for m in msgs)
d = zlib.decompressobj(zdict=dictionary)
assert d.decompress(with_dict(msgs[7])) == msgs[7]
print(f"raw {sum(map(len, msgs))}  zlib {plain}  zlib+dict {primed}")

raw 5890  zlib 5890  zlib+dict 2090

Each JSON event is about 59 bytes. Alone, zlib cannot shrink it at all. With a 65-byte shared dictionary, the total drops by 65%.

import zlib
from math import sin

signal = bytes(int(128 + 100 * sin(i / 40)) for i in range(20000))  # smooth wave

def delta(d):                                   # PNG "Sub" filter idea
    return bytes([d[0]] + [(d[i] - d[i - 1]) % 256 for i in range(1, len(d))])

def undelta(d):
    out = bytearray([d[0]])
    for b in d[1:]:
        out.append((out[-1] + b) % 256)
    return bytes(out)

assert undelta(delta(signal)) == signal
print("zlib raw  :", len(zlib.compress(signal, 9)))
print("zlib delta:", len(zlib.compress(delta(signal), 9)))

zlib raw  : 1962
zlib delta: 1337

Neighboring samples of a smooth signal are close. Their differences are small numbers, near 0. Those have a skewed distribution that the entropy coder loves. Same data, 32% smaller.

Pattern: predict, then code the residual. A residual close to 0 has low entropy. This is the same "model + coder" split again.

Lossy codecs push this further. Their transform (DCT, wavelet) packs energy into a few coefficients, so the rest can be thrown away.

import random
from math import log2

random.seed(3)
xs = [random.gauss(0, 1) for _ in range(100_000)]      # sigma = 1

def quantize(x, step):
    return step * round(x / step)                     # uniform mid-tread quantizer

print(" step   MSE     entropy(bits)  Shannon D(R)")
for step in (2.0, 1.0, 0.5, 0.25):
    q = [quantize(x, step) for x in xs]
    mse = sum((a - b) ** 2 for a, b in zip(xs, q)) / len(xs)
    counts = {}
    for v in q:
        counts[v] = counts.get(v, 0) + 1
    R = -sum(c / len(q) * log2(c / len(q)) for c in counts.values())
    print(f"{step:5}  {mse:.4f}   {R:.3f}          {2 ** (-2 * R):.4f}")

 step   MSE     entropy(bits)  Shannon D(R)
  2.0  0.3314   1.244          0.1784
  1.0  0.0832   2.104          0.0541
  0.5  0.0209   3.061          0.0144
 0.25  0.0052   4.050          0.0036

A simple uniform quantizer plus ideal entropy coding lands at a fixed gap from the bound. At high rate, MSE / D(R) → πe/6 ≈ 1.42. That is 1.53 dB, or 0.25 bits per sample. At step 0.25 we measure 0.0052 / 0.0036 = 1.44. Vector quantization and trellis coding close some of that gap.

from math import cos, pi, sqrt

def dct(x):                                    # orthonormal DCT-II
    N = len(x)
    return [sqrt((1 if k == 0 else 2) / N) *
            sum(x[n] * cos(pi * (n + 0.5) * k / N) for n in range(N))
            for k in range(N)]

def idct(X):                                   # its inverse (DCT-III)
    N = len(X)
    return [sum(sqrt((1 if k == 0 else 2) / N) * X[k] * cos(pi * (n + 0.5) * k / N)
                for k in range(N)) for n in range(N)]

row = [52, 55, 61, 66, 70, 61, 64, 73]         # 8 pixels from a JPEG test block
X = dct(row)
print("DCT :", [round(v, 1) for v in X])
Q = [round(v / 10) for v in X]                 # quantize with step 10
print("Q   :", Q)
back = idct([q * 10 for q in Q])
print("back:", [round(v) for v in back])
print("max error:", round(max(abs(a - b) for a, b in zip(row, back)), 2))

DCT : [177.5, -14.4, -5.7, -6.7, 7.1, -3.1, -0.7, 2.4]
Q   : [18, -1, -1, -1, 1, 0, 0, 0]
back: [53, 55, 64, 74, 70, 60, 61, 72]
max error: 7.6

MP3/AAC use the MDCT, a lapped version that avoids block edges. JPEG 2000 uses wavelets. H.264 and later use small integer DCTs, so encoder and decoder match bit for bit.

1. No free lunch: most inputs cannot shrink. Compression is a bet on structure.
2. Model + coder. Arithmetic coding and ANS spend within a few bits of Σ −log₂ p. So compression = prediction.
3. Huffman rounds to whole bits. Arithmetic coding does not.
4. Adaptive models learn as they go. Learning costs real bits.
5. LZ points back at repeats. It is universal but converges slowly, so it is paired with entropy coding: DEFLATE, zstd, xz.
6. BWT sorts contexts together. CM mixes many models. Both beat LZ on text.
7. Lossy: transform → quantize → entropy code, bounded by R(D).