AlgoScope

Suffix Array

structureadvancedTime O(n log n) build, O(m log n) searchSpace O(n)

Every substring of a text is a prefix of some suffix, so if the suffixes are listed in sorted order, every pattern that occurs sits in one contiguous block of that list and a binary search finds it in O(m log n). The suffix array is just that list, stored as starting positions. Sorting the suffixes directly works but compares long strings; prefix doubling builds the same order by ranking each position on its first character, then its first two, four, eight, each round sorting pairs of ranks from the round before, so no comparison ever looks at more than two numbers.

startsuffix0123455a3ana1anana0banana4na2nana

Does "nan" occur in "banana"? The rows are the suffixes in sorted order, which is the suffix array [5, 3, 1, 0, 4, 2]. Any occurrence of the pattern is the start of some suffix, so the pattern is a prefix of that suffix, and in sorted order every suffix with that prefix sits in one contiguous block. Binary search for the block.

Check your understanding

The player pauses before each decision in this run and asks what happens next. Here are all 3, with their answers.

  1. Row 2 is "anana" and the pattern is "nan". What now?

    • Found: the pattern is a prefix here
    • Pattern is smaller: go left
    • Pattern is larger: go right

    Answer: Pattern is larger: go right. Compare the pattern with the suffix's first 3 characters: equal means found, smaller means look above, larger means look below.

  2. Row 4 is "na" and the pattern is "nan". What now?

    • Found: the pattern is a prefix here
    • Pattern is smaller: go left
    • Pattern is larger: go right

    Answer: Pattern is larger: go right. Compare the pattern with the suffix's first 3 characters: equal means found, smaller means look above, larger means look below.

  3. Row 5 is "nana" and the pattern is "nan". What now?

    • Found: the pattern is a prefix here
    • Pattern is smaller: go left
    • Pattern is larger: go right

    Answer: Found: the pattern is a prefix here. Compare the pattern with the suffix's first 3 characters: equal means found, smaller means look above, larger means look below.

How it runs, step by step

  1. Does "nan" occur in "banana"? The rows are the suffixes in sorted order, which is the suffix array [5, 3, 1, 0, 4, 2]. Any occurrence of the pattern is the start of some suffix, so the pattern is a prefix of that suffix, and in sorted order every suffix with that prefix sits in one contiguous block. Binary search for the block.

    Binary search for nan in the suffix array of banana.

  2. Rows 0 to 5 are still possible; the middle row 2 holds "anana" (start 1). "nan" sorts after it, so the block can only be below: search rows 3 to 5.

    Row 2: anana.

  3. Rows 3 to 5 are still possible; the middle row 4 holds "na" (start 4). "nan" sorts after it, so the block can only be below: search rows 5 to 5.

    Row 4: na.

  4. Rows 5 to 5 are still possible; the middle row 5 holds "nana" (start 2). It begins with "nan": found, at text position 2.

    Row 5: nana.

  5. "nan" occurs at position 2, found in 3 probes. Each probe compares at most 3 characters, so a search costs O(m log n) however long the text: the suffix array is the index that makes a text searchable for any pattern without knowing the patterns in advance.

    Found at 2.

Remember

  • The suffix array lists suffix start positions in sorted order; any pattern is a prefix of some suffix.
  • Prefix doubling: rank by 1 character, then sort pairs (rank at i, rank at i + k) to rank by 2k.
  • Search is a binary search over the rows comparing at most m characters per probe: O(m log n).

Topics covered

Related

Where this is used

Compressionbzip2 and the Burrows-Wheeler transform

bzip2 compresses a block by sorting all of its rotations and keeping only the last column, which is then run-length and Huffman coded. Sorting rotations and sorting suffixes agree once a unique end marker is in play, which is why suffix sorting libraries ship a transform entry point at all: libdivsufsort exposes divbwt next to divsufsort and they are the same sort. bzip2 itself omits the marker and records the row number of the original block instead. The ratio comes entirely from that ordering: rows next to each other share a long prefix, so the characters that precede them arrive in long runs of the same byte.

BioinformaticsAligning sequencing reads to a genome

Bowtie and BWA index a three billion character reference genome once, then place hundreds of millions of short reads into it. A suffix tree over that genome would need tens of bytes per base and would not fit in memory, while the suffix array is four bytes per base and the FM-index derived from it, which keeps the Burrows-Wheeler transform plus a sampled suffix array, is smaller still. Each read is a pattern whose occurrences form one contiguous block of the suffix order, so the aligner narrows a range of rows instead of scanning the genome.

Operationsbsdiff and Chrome's Courgette

Shipping an update as a patch needs, for every offset in the new binary, the longest run of bytes that already appears somewhere in the old one. bsdiff builds a suffix array of the old file using Larsson and Sadakane's qsufsort and answers each of those queries with a binary search, which is how it produces patches many times smaller than compressing the new file alone. Courgette disassembles both Chrome binaries first so that shifted jump targets stop looking like edits, then hands the result to the same suffix-array diff.

Developer toolsTraining a zstd dictionary

Compressing thousands of small records separately throws away every pattern the records share, which is what zstd's dictionary trainer exists to recover. Its original algorithm, still reachable as zstd --train-legacy, concatenates the samples, calls divsufsort on the whole buffer, and walks the resulting suffix array for segments that recur across samples: in sorted suffix order a repeat of any length is just a run of adjacent rows sharing that prefix. The current default, --train-fastcover, gave that up for speed and counts fixed-length dmers instead, so it only sees repeats at the length it samples.

Why it works this way

How to get the m out of O(m log n) search

The m is work repeated: two probes that land on suffixes sharing a long prefix each re-compare that prefix from its first character. Manber and Myers showed you can do better by storing, beside the array, the longest common prefix each suffix shares with the one before it. Carrying the count of characters already known to match into the next probe means no character is ever compared twice, which brings search down to O(m + log n). That LCP array is also what turns a suffix array into a substitute for a suffix tree rather than just a search index, and Kasai's algorithm builds it in O(n) from the suffix array and its inverse.

Why a sentinel character that sorts below everything else

Without one, a suffix can be a prefix of another suffix - "a" and "ana" in banana - so their order is settled by a rule about running out of characters rather than by a comparison. Appending a character smaller than every real one makes every suffix end differently, so no suffix is a prefix of another and every comparison is decided by an actual character. The doubling code here does the same job with its -1 for positions past the end of the string, since a rank of -1 sorts below every real rank. It matters downstream too: inverting the Burrows-Wheeler transform requires knowing which sorted row was the original string, and a sentinel gives that away for free, since exactly one row ends with it. Implementations that leave the sentinel out, bzip2 among them, have to store that row number beside the transformed block instead.

The rank array has to be rebuilt, not updated in place

Round k derives each position's new rank from the pair (rank[i], rank[i + k]), and both halves of that pair are ranks from the previous round. Writing new ranks straight back into rank while the pass is still reading it means some positions are compared at 2k characters and some at k, and the order that comes out is quietly wrong rather than crashing. The code here fills a fresh next array and only swaps it in once the round is finished, which is the entire reason the swap exists.

Why the stated build cost is O(n log n) when this code is not

The loop runs log n rounds and the version here calls a comparison sort inside each one, which is O(n log^2 n). Both halves of the key are ranks smaller than n, so a two-pass radix sort - by the second element, then a stable pass by the first - sorts a round in O(n) and brings the total to O(n log n). Genuinely linear constructions exist as well: DC3 from Karkkainen and Sanders, and SA-IS from Nong, Zhang and Chan, whose induced sorting is what Chromium ships in Zucchini. Linear does not mean fastest, though: libdivsufsort, the sorter most production code reaches for, is O(n log n) in the worst case and still beats the linear methods on real inputs. Doubling is still worth knowing because the rank array it produces on the way is the inverse suffix array, which Kasai's LCP algorithm and many string problems need anyway.

Read more

Next up