Back to String Processing Dsa

FM-Index Backward Search

September 17, 2026 Wasil Zafar 7 min read

LF-mapping, rank queries, and compressed pattern intervals

Table of Contents

  1. BWT Columns and LF-Mapping
  2. Backward Search
  3. Rank, Select, and Wavelet Structures
  4. Building a Practical FM-Index

BWT Columns and LF-Mapping

Append a unique sentinel $ that is lexicographically smaller than every text symbol. Sorting the suffixes of banana$ gives suffix-array order $[6,5,3,1,0,4,2]$. The first column $F$ contains the first symbols of those sorted suffixes, while the Burrows-Wheeler column $L$ contains the symbol immediately preceding each suffix in the original cyclic text; here, $L=\texttt{annb\$aa}$.

Stable sorting preserves the relative order of equal symbols, so the $k$-th occurrence of a symbol in $L$ corresponds to the $k$-th occurrence of that symbol in $F$. Let $C[c]$ be the number of text symbols strictly smaller than $c$, and let $\operatorname{rank}_c(i)$ count occurrences of $c$ in the half-open prefix $L[0\dots i)$. Then

$$ \operatorname{LF}(i)=C[L[i]]+\operatorname{rank}_{L[i]}(i). $$

Moving from row $i$ to $\operatorname{LF}(i)$ prepends $L[i]$ to that row's suffix. Repeating LF steps therefore walks backward through the text and is also the basis of BWT inversion.

Row $i$$\text{SA}[i]$Sorted suffix$F$ occurrence$L$ occurrence$\operatorname{LF}(i)$
06$$ #1a #11
15a$a #1n #15
23ana$a #2n #26
31anana$a #3b #14
40banana$b #1$ #10
54na$n #1a #22
62nana$n #2a #33

All suffix-array rows initially match the empty pattern, so backward search starts with the half-open interval $[l,r)=[0,n)$. To prepend character $c$ to the portion of the pattern already processed, restrict the interval to rows whose preceding BWT character is $c$:

$$ l'=C[c]+\operatorname{rank}_c(l),\qquad r'=C[c]+\operatorname{rank}_c(r). $$

The new interval is empty when $l'\ge r'$. Otherwise, after processing the pattern from right to left, every suffix-array row in $[l,r)$ begins with the complete pattern. Its occurrence count is therefore $r-l$; locating occurrences additionally requires recovering the suffix-array values for those rows.

Worked Example

Backward Search

Search for ana in banana$, whose BWT is annb$aa. The cumulative counts are $C[\$]=0$, $C[a]=1$, $C[b]=4$, and $C[n]=5$.

  1. Begin with every row: $[0,7)$.
  2. Read the pattern's final a: $[0,7)\rightarrow[1,4)$.
  3. Prepend n: $[1,4)\rightarrow[5,7)$.
  4. Prepend the first a: $[5,7)\rightarrow[2,4)$.

Rows 2 and 3 are the suffixes ana$ and anana$. Their suffix-array values are 3 and 1, so the reported text positions are $[1,3]$ after sorting. Notice that backward search performs one pair of rank queries per pattern symbol and never compares the pattern against individual suffixes.

Rank, Select, and Wavelet Structures

A dense occurrence table stores $\operatorname{rank}_c(i)$ for every character and position, giving constant-time rank but using $O(n|\Sigma|)$ space. Checkpointed counts reduce space and scan a short block per query. For small alphabets, one bitvector per symbol plus word-level popcount is often both compact and fast.

A wavelet tree recursively partitions the alphabet and stores one routing bitvector at each node. A rank query follows the queried symbol's root-to-leaf path, translating the prefix length with rank at each level, for $O(\log |\Sigma|)$ time and approximately $n\log |\Sigma|$ bits before succinct-bitvector overhead. Wavelet matrices improve layout; Huffman-shaped trees favor frequent symbols; and run-length-compressed BWT structures exploit highly repetitive texts.

Wavelet-tree rank query for n in the BWT annb$aa The first six BWT symbols contain three characters routed right at the root. In the right child, two of those route to n, so rank n of six equals two. rankₙ(6) in L = annb$aa Root: {$, a} → 0; {b, n} → 1 0 1 1 1 0 0 0 Left subsequence: a $ a a not followed for symbol n Right subsequence: n n b 1 1 0 root rank₁(6) = 3; child rank₁(3) = 2; therefore rankₙ(6) = 2

Building a Practical FM-Index

Construction needs the suffix array, BWT, cumulative counts, and a rank structure. Counting stops after backward search and can remain highly compressed. Locating must recover suffix-array positions: sample selected SA values, repeatedly apply LF from an unsampled row until reaching a sample, and return $(\text{sample}+\text{steps})\bmod n$. Sampling every $d$ text positions uses roughly $n/d$ samples but can require up to $d-1$ LF steps per located occurrence.

The compact implementation below uses integer encoding with sentinel 0, a dense prefix-count table for clarity, and samples suffix positions divisible by the chosen rate. Production indexes replace the naive suffix sort and dense rank table with scalable suffix-array construction and succinct rank/select structures, while preserving the same search equations.

from dataclasses import dataclass
from typing import Dict, List


@dataclass(frozen=True)
class SearchInterval:
    left: int
    right: int

    @property
    def count(self) -> int:
        return self.right - self.left


class FMIndex:
    """Small educational FM-index with count and sampled locate support."""

    def __init__(self, text: str, sample_rate: int = 4) -> None:
        if sample_rate <= 0:
            raise ValueError("sample_rate must be positive")

        alphabet = sorted(set(text))
        self.to_int = {character: index + 1 for index, character in enumerate(alphabet)}
        self.to_char = {value: key for key, value in self.to_int.items()}
        encoded = [self.to_int[character] for character in text] + [0]
        self.n = len(encoded)

        # Clear but intentionally naive construction; replace for large inputs.
        self.sa = sorted(range(self.n), key=lambda start: encoded[start:])
        self.bwt = [encoded[(start - 1) % self.n] for start in self.sa]

        counts: Dict[int, int] = {}
        for symbol in self.bwt:
            counts[symbol] = counts.get(symbol, 0) + 1

        self.cumulative: Dict[int, int] = {}
        total = 0
        for symbol in sorted(counts):
            self.cumulative[symbol] = total
            total += counts[symbol]

        self.occ: Dict[int, List[int]] = {
            symbol: [0] * (self.n + 1) for symbol in counts
        }
        for index, symbol in enumerate(self.bwt):
            for current, prefix in self.occ.items():
                prefix[index + 1] = prefix[index] + (current == symbol)

        self.sa_samples = {
            row: suffix_start
            for row, suffix_start in enumerate(self.sa)
            if suffix_start % sample_rate == 0
        }

    def rank(self, symbol: int, end: int) -> int:
        """Occurrences of symbol in BWT[0:end]."""
        prefix = self.occ.get(symbol)
        return 0 if prefix is None else prefix[end]

    def lf(self, row: int) -> int:
        symbol = self.bwt[row]
        return self.cumulative[symbol] + self.rank(symbol, row)

    def interval(self, pattern: str) -> SearchInterval:
        left, right = 0, self.n
        for character in reversed(pattern):
            symbol = self.to_int.get(character)
            if symbol is None:
                return SearchInterval(0, 0)
            left = self.cumulative[symbol] + self.rank(symbol, left)
            right = self.cumulative[symbol] + self.rank(symbol, right)
            if left >= right:
                return SearchInterval(0, 0)
        return SearchInterval(left, right)

    def count(self, pattern: str) -> int:
        return self.interval(pattern).count

    def _locate_row(self, row: int) -> int:
        steps = 0
        while row not in self.sa_samples:
            row = self.lf(row)
            steps += 1
        return (self.sa_samples[row] + steps) % self.n

    def locate(self, pattern: str) -> List[int]:
        interval = self.interval(pattern)
        return sorted(self._locate_row(row) for row in range(interval.left, interval.right))

    def decoded_bwt(self) -> str:
        return "".join("$" if symbol == 0 else self.to_char[symbol] for symbol in self.bwt)


def naive_occurrences(text: str, pattern: str) -> List[int]:
    if pattern == "":
        return list(range(len(text) + 1))
    return [
        start
        for start in range(len(text) - len(pattern) + 1)
        if text.startswith(pattern, start)
    ]


if __name__ == "__main__":
    text = "banana"
    index = FMIndex(text, sample_rate=2)

    assert index.decoded_bwt() == "annb$aa"
    assert index.interval("ana") == SearchInterval(2, 4)
    assert index.count("ana") == 2
    assert index.locate("ana") == [1, 3]
    assert index.locate("na") == [2, 4]
    assert index.count("xyz") == 0

    for pattern in ["", "a", "ana", "banana", "nana", "x"]:
        assert index.locate(pattern) == naive_occurrences(text, pattern)
        print(f"{pattern!r}: count={index.count(pattern)}, positions={index.locate(pattern)}")