A needle in a haystack

The last chapter left us with a problem: find a word in a text three billion letters long faster than by trying every position. That text is the human genome. A copy of it sits in almost every cell of your body, and its alphabet has four letters: A, C, G and T. The chapter begins in a genetics lab, where substring search is everyday work, and ends with a school committee going through suspiciously similar essays. Both will need the same tools.

The lab. Three billion letters

We’ll start with a genome smaller than ours but larger than φX174’s. Bacteriophage lambda is a workhorse of molecular biology: its DNA is sold in test tubes, and new methods are tried out on it first. Its genome is 48,502 letters long. Labs keep genomes in the FASTA format: the first line starts with a > sign and describes the sequence, and the letters follow, 70 to a line. Files with three genomes from the open NCBI database are in the folder /data/strings-search: phage φX174, phage lambda and the human mitochondrion.

The lab’s first job is to predict an experiment. Bacteria have protein scissors called restriction enzymes: each one recognizes its own short motif in DNA and cuts the molecule there. The enzyme EcoRI, from E. coli, looks for the motif GAATTC and cuts it after the first letter. Know every occurrence of the motif, and you can say in advance which pieces the DNA will fall into. For now we’ll search with the find method from Chapter 7: it returns the position of an occurrence or −1, and its second argument says where to start looking.

Five cuts, six pieces: 21,226, 7,421, 5,804, 5,643, 4,878 and 3,530 letters. Treat lambda DNA with EcoRI and run it through a gel in an electric field, and the pieces separate by length: six bands of these lengths appear on the gel. The mixture is sold ready-made and used as a ruler for measuring other molecules. A dozen lines of Python have predicted an experiment done at the lab bench.

Read the motif GAATTC backward and replace every letter with its partner, A with T and C with G and vice versa. You get GAATTC again. DNA has two strands, and the second reads this way, so the enzyme recognizes the motif on both strands at once. Biologists call such motifs palindromes, though they are palindromes with a twist, not the kind we checked in Chapter 7.

Here find finished instantly. But lambda is more than sixty thousand times shorter than the human genome. To know whether the search will survive three billion letters, we need to know how it works inside. We’ll write it ourselves.

The ruler

The problem is called substring search: given a text of length $n$ and a pattern of length $m$, find every position where the pattern occurs in the text. The first solution suggests itself. Lay the pattern against the start of the text like a ruler and compare letter by letter. On a mismatch, slide the ruler one position along and start over. While we’re at it, we’ll count how many letter comparisons this takes.

On the lambda genome the ruler does well: 1.36 comparisons per letter. The first letter of the pattern almost always fails to match, and the ruler moves straight on. Python would get through three billion letters this way in a few minutes. The second half of the output is another story. In a text of nothing but A, the pattern AAA…AC matches almost to the end at every position, and only its last letter rejects it. Each position costs $m$ comparisons, about $n \cdot m$ in all: a pattern ten times longer takes ten times as long. This is the worst case, and for a genome of three billion letters and a pattern of a thousand it means $3 \cdot 10^{12}$ comparisons, more than a day of work.

A text of nothing but A looks made up, but nature writes texts like it. The ends of every human chromosome carry telomeres: six letters, TTAGGG, repeated hundreds and thousands of times. Repeats, simple and complicated, take up about half of the human genome. For the naive search such a text is a minefield: at every repeat the ruler runs through a long match and then starts over from the next position.

The text is TTAGGG repeated 20,000 times (120,000 letters). The pattern is TTAGGG repeated 200 times with TTAGGC at the end: it almost matches the text but never occurs in it. How many times more comparisons will the naive search make here than on a text of the same length made of random letters?

about twice as many about twenty times as many more than a hundred times as many the same: the pattern won’t be found anyway

About 24 million comparisons against 160 thousand, a hundred and fifty times as many. At every sixth position the ruler goes through more than a thousand matching letters before the final C stops it. We’ll check this below, once we have something to compare it with.

The waste is easiest to see at a position where 999 letters out of a thousand matched. We have read 999 letters of the text and know for certain what they are: the same as the start of the pattern. Then we slide the ruler one position and read them again as if we had never seen them, and the work of all those comparisons is thrown away. To save it, we first have to find something out about the pattern itself.

What the pattern knows about itself

Take the pattern abracadabra and imagine that the ruler, laid at position $i$, matched the text on its first ten letters, abracadabr, and failed on the eleventh. Which shifts of the ruler make sense? We know those ten letters of the text without rereading anything: they are abracadabr. A shift by one puts the pattern’s first a against b, which is useless. The only useful shift is one after which the start of the pattern coincides with the end of the piece already read. A shift by seven will do: then the abr at the start of the pattern lands on the abr at the end of what we’ve read, and the comparison carries on from the fourth letter of the pattern. The text needs no rereading, since we already know those three letters match.

Everything, then, hinges on a question with no text in it at all: what is the longest piece that the pattern both begins and ends with? Such a piece is called a border of the string: a beginning of the string that is also its end, but not the whole string. The longest border of abracadabr is abr, that of the whole word abracadabra is abra, and abc has no borders. The lengths of the longest borders of all the beginnings of the pattern make up a table, the prefix function $\pi$: $\pi[i]$ is the length of the longest border of the beginning p[:i + 1].

Computing it by brute force, taking each beginning and trying every border length, is slow: on the order of $m^3$ steps. But the table can be built out of itself. Suppose we know $k = \pi[i - 1]$: the beginning of length $k$ coincides with the end of p[:i]. If the next letter p[i] equals p[k], the border grows by one letter, and $\pi[i] = k + 1$. If not, we need a shorter border that will grow. The next longest border of a string is the longest border of its border, that is, $\pi[k - 1]$, and that is already in the table. We try it, then the border of the border’s border, until one grows or $k$ drops to zero.

Under each letter is the length of the longest border of the word up to and including that letter. At the end of abracadabra stands 4, the border abra. In the telomere pattern a second TTAGGG begins at the seventh position, and the border grows 1, 2, 3, 4, 5 until the final C cuts it back to zero. Press “Steps” to watch k fall back in the while loop.

Never a step back

Now the search. We walk through the text from left to right and remember one number, $k$: how many of the pattern’s first letters match the end of what we have read. We read the next letter. If it equals pattern[k], $k$ grows. If not, the pattern shifts according to the border table: $k$ becomes $\pi[k - 1]$, and the same letter of the text is compared again, now with a different letter of the pattern. When $k$ reaches $m$, the pattern is found. The pointer into the text never moves back: each letter of the text is read once, and all the backing up happens inside the pattern, by the table. This is the Knuth–Morris–Pratt algorithm, KMP for short.

On the random text both searches make 1.2–1.3 comparisons per letter. On the telomere the naive search makes twenty-four million, and KMP a hundred and forty thousand, even fewer than on the random text. That settles the bet above: a hundred and fifty times. What remains is to prove that KMP is this fast on any text.

The last line of the output is the built-in in. On the same telomere it takes a fraction of a millisecond: it is written in C, and since version 3.10 it uses yet another linear algorithm for long patterns, the Two-Way algorithm published by Maxime Crochemore and Dominique Perrin in 1991. Before 3.10, in and find had a quadratic worst case, like our ruler.

In everyday work, write in and find. Knowing how linear search works pays off where the built-in methods can’t help: in a stream that arrives one letter at a time and doesn’t fit in memory, in a search for many patterns at once, in a list of numbers or events instead of a string.

The number $k$ in KMP is a state: “how many letters of the pattern have matched so far.” A new letter moves the state to another one, and the algorithm remembers nothing else. That is how a finite automaton works; automata get a whole chapter later on. The automaton is built from a single pattern and serves any text. Below you can build one for a word of your own.

The fallbacks can be worked out in advance

If the alphabet is small, as DNA’s is, you can record in advance, for every state $k$ and every letter, which state the automaton will land in once all the fallbacks are done. Then each letter of the text costs one table lookup and no while loop at all. The table takes $(m + 1) \cdot |\Sigma|$ cells, where $|\Sigma|$ is the size of the alphabet: nothing for DNA, a great deal for Unicode. This is how programs like grep search for regular expressions: the expression is turned into an automaton, and the text is run through it letter by letter.

Fingerprints: Rabin and Karp

KMP compares letters, but we can compare something bigger. In Chapter 16 we turned a string into a number with a polynomial hash: multiply what you have by the base, add the code of the letter, take the remainder. We compute such a hash for the pattern, and then for every window of length $m$ in the text. Where the hashes differ, the strings certainly differ. Where they agree, we have most likely found the pattern, and a comparison will confirm it. The hash works here as a fingerprint: short, and yet it almost always identifies its owner.

The catch is that there are $n - m + 1$ windows, and if each hash is computed from scratch, we are back to the ruler’s $n \cdot m$ steps. But neighboring windows are almost the same: when the window slides by one letter, one letter leaves on the left and one arrives on the right. By Horner’s method the hash of the window $c_i c_{i+1} \ldots c_{i+m-1}$ is $c_i B^{m-1} + c_{i+1} B^{m-2} + \ldots + c_{i+m-1}$ modulo $p$, where $B$ is the base. To get the hash of the next window, subtract the share of the departing letter, multiply by $B$ and add the arriving one:

$$h_{i+1} = \bigl( (h_i - c_i \cdot B^{m-1}) \cdot B + c_{i+m} \bigr) \bmod p.$$

Three arithmetic operations per slide, however long the pattern. A hash updated this way is called rolling. The power $B^{m-1}$ is computed once, in advance: the built-in pow(base, m - 1, mod) can raise to a power modulo a number directly. We’ll pick a random base; the theorem about the secret base in Chapter 16 explained why.

The five EcoRI sites turn up with every modulus: equal strings have equal hashes, so the hash never misses a true occurrence. The difference lies in the false alarms, when the hashes agree and the strings don’t. With modulus 7 the hash takes only seven values, and thousands of windows raise an alarm. With 101 it’s hundreds; with 10,007, anywhere from none to a few dozen; with $2^{61} - 1$, none. Run the cell a few times: the base is random, and the numbers will change. There are about $n/p$ false alarms, because each window hits the pattern’s hash with a probability of about $1/p$.

The sliding formula has a subtraction in it, and h - ord(...) * top can be negative. In Python the remainder of division by a positive number is never negative, -5 % 7 is 2, and the formula works as written. In C and Java the remainder of a negative number is negative, -5 % 7 is −5 there, and you have to add the modulus to it. Anyone who ports hashing code from Java to Python or back trips over this sooner or later.

The algorithm was published by Richard Karp and Michael Rabin in 1987. Both won the Turing Award, and Rabin is also the second name in the Miller–Rabin primality test from the last chapter. Their search is randomized too, and it shows both kinds of randomized algorithm. Our version checks every match of the hashes, so the answer is always right and only the time is a matter of chance: with an unlucky base there are more false alarms and the search is slower. That makes it a Las Vegas algorithm. Drop the check and the search saves time on every alarm but is now and then wrong: a Monte Carlo algorithm. By the secret-base theorem, the probability that one given window raises a false alarm is at most $(m - 1)/p$. For lambda, a six-letter pattern and $p = 2^{61} - 1$, an error anywhere in the search will happen on average no more than once in nine trillion runs.

KMP guarantees linear time with no probability involved, yet fingerprints have a job of their own: they can do what an automaton for one pattern can’t. Given a thousand patterns, all of the same length, put their hashes into a set, and each window of the text will cost one check, h in hashes, however many patterns there are. And compare the fingerprints of all the pieces of one text with those of all the pieces of another, and you have a plagiarism checker. We’ll come to that in the second half of the chapter.

A catalog of scissors: the trie

Many restriction enzymes are known, and their motifs differ in length: AluI and HaeIII have four letters each, EcoRI six, NotI eight. Rabin–Karp fingerprints compare only windows of one length, so the catalog needs a different idea. We put all the motifs into one tree. Arrows leave the root, each labeled with a letter that some motif begins with; from every node below come arrows labeled with the letters that can follow. Motifs with a common beginning share a path: GATC and GATATC part ways only at the fourth letter. The node where a motif ends is marked with the enzyme’s name.

Such a tree is called a trie. René de la Briandais described it in 1959, and in 1960 Edward Fredkin named it trie, from the middle of the word retrieval. Fredkin pronounced it “tree”; others say “try,” to tell the two apart. In Python a trie is easy to build out of dictionaries: a node is a dictionary “letter → next node,” and the end of a motif is marked with a special key.

The setdefault(ch, {}) method returns the value for a key, and if the key is missing, it first puts an empty dictionary there: that way the path through the trie gets laid without extra checks. Fourteen enzymes are looked up in one pass over the genome, at fewer than two trie steps per letter. The four-letter motifs occur in lambda more than a hundred times each, where a random text would give about $48\,502 / 4^4 \approx 190$, and the eight-letter NotI doesn’t occur at all. That is why geneticists love NotI: it cuts a big genome into a few long pieces.

A trie plus a prefix function: Aho–Corasick

Our walk through the trie starts over from every position, like the ruler, and in the worst case costs $n \cdot L$ steps, where $L$ is the length of the longest motif. The cure is the one KMP uses: work out in advance, for every node of the trie, where to fall back on a mismatch, to the node of the longest border, only now the border is sought among all the motifs at once. The result is an automaton that passes through the text once without going back and finds every occurrence of every motif in $O(n + \text{total length of the motifs} + \text{number of matches})$. Alfred Aho and Margaret Corasick invented it at Bell Labs in 1975, and it powered the fgrep utility, which searches files for many words at once.

The same trie lives in every phone. You type “nat,” and the keyboard offers “Natasha,” so somewhere there is a dictionary that quickly returns every word with a given beginning. In a sorted list or a hash table you would have to go through the words; in a trie you walk down the three letters of the prefix, and everything you need lies in one subtree. We’ll build a trie from all 562 thousand words of War and Peace and count in each node how many times its word occurred, since it makes sense to order the suggestions by frequency. A node is now more convenient as a class, as in Chapter 17.

The 20,478 different words hold 163,930 letters between them, yet the trie has only 63,617 nodes, two and a half times fewer: common beginnings are stored once. A query costs as many steps as the prefix has letters, plus the walk through the subtree, and the size of the dictionary doesn’t lengthen the descent. The suggestions are right but naive. To “bor” the novel answers with Boris, Borodino and the borzois of the wolf hunt; to “war,” with “warm.” And to “kutu,” after Kutuzov and his possessive, it offers “kutuzov—having” and “kutuzov—the”: this edition sets its dashes without spaces, and our word cutter takes such a pair for a single word. Phone keyboards also take the previous word into account, like the word chains in Chapter 8, and your own habits as well. Try it yourself.

In Chapter 48 the trie comes back as part of a search engine for our textbooks. For now the lab closes, and we move to the staff room.

The committee. Three essays

The literature teacher has brought in three essays on the same topic, “The oak in War and Peace.” They were written by Amy, Ben and Clara. Amy’s and Ben’s essays strike her as suspiciously alike, but Ben says it’s a coincidence: one topic, one novel. The committee needs a measure of similarity that can be computed rather than felt.

The first thing that comes to mind is diff from Chapter 22: the longest common subsequence of words. It copes with insertions and replacements but handles rearrangement badly. Ben moved the sentence “The oak in the novel is…” from the beginning to the middle, and a common subsequence can take it in only one of the two places. Besides, the table costs $n \cdot m$ for each pair, and a school with a thousand essays has half a million pairs. We need a measure that ignores the order of the pieces and is quick to compute.

Cut each text into shingles: runs of $k$ consecutive words, overlapping like shingles on a roof, which is where the name comes from. A text of $N$ words gives $N - k + 1$ shingles. The similarity of two texts is the share of common shingles among all of them: $J(A, B) = \frac{|A \cap B|}{|A \cup B|}$. The measure was introduced in 1901 by the Swiss botanist Paul Jaccard, who compared the flora of different parts of the Alps and the Jura this way: how many plant species two areas have in common. It is called the Jaccard index. Andrei Broder proposed shingles and this measure for documents in 1997: the AltaVista search engine used them to find nearly identical pages on the web.

Our normalize is the words function from Chapter 7 with one addition: the curly apostrophe ’ becomes a straight one. Otherwise “doesn’t” in Amy’s essay, typed in a word processor that curls apostrophes, and “doesn't” in Ben’s, typed in a plain editor, would count as different words, and the edition of the novel we’ll check against later uses curly ones throughout. With $k = 1$ the shingles are single words, and even the independent essays of Amy and Clara share a fifth of theirs: “oak,” “prince,” “andrew,” “tolstoy,” and of course “the” and “and.” By $k = 3$ the matches between independent pieces of work are gone, while Amy and Ben hold on: at $k = 4$ they share 71 shingles, almost half of all they have. One matching word is chance; four words in a row matching seventy-one times is not. The moved sentence produced shared shingles too, since the measure ignores the order of the pieces.

Winnowing

Three essays can be compared pair by pair. But the committee wants more: to check each new essay against all the school’s past essays, and against the novel itself, which can be copied from too. War and Peace has more than half a million shingles, and ten years of the archive add millions more. Storing them all is expensive. We would like to keep a small sample, the document’s fingerprints, chosen so that no noticeable borrowing gets lost.

The obvious ways fail. Taking every fourth shingle won’t do: let Ben add one word at the beginning, and all the numbers shift, so the copied piece yields different shingles from Amy’s. We could choose shingles by the hash itself, say only those whose hash is divisible by four. That choice survives shifts but gives no guarantee: a long shared piece may contain no hash divisible by four at all. The solution was found in 2003 by Saul Schleimer, Daniel Wilkerson and Alex Aiken. Compute the hashes of all the shingles in order and pass a window of width $w$ over them. In each window, record the smallest hash as a fingerprint. Neighboring windows often pick the same minimum, so there are few fingerprints. The method is called winnowing, the old word for separating grain from chaff.

A borrowing shorter than this threshold may slip through, but that does no harm: matches of three or four words happen by chance. The share of fingerprints can be estimated in advance. If the hashes behave like random numbers, a window of width $w$ puts on average about $\frac{2}{w + 1}$ of all shingles into the fingerprints: 40 percent at $w = 4$, four percent at $w = 50$. We’ll test it all at once on a harder case. A fourth student, Dan, has handed in a description of the oak. Where did it come from?

More than half a million words of the novel turned into about 223 thousand fingerprints, forty percent, as promised. The index takes under a second to build, and a lookup in it is a few dictionary queries: once again the hash table does all the work. Dan added a word, reworded one sentence, threw out another and cut half of a third, but more than a dozen fingerprints matched and pointed straight at the source: Book Six, Chapter I, the spring of 1809, Prince Andrew on his way to his Ryazan estates.

The hashes here come from the built-in hash, and as we saw in Chapter 16, it has a salt. Within one run this makes no difference: the index and the query are computed by the same function. But such an index can’t be stored on disk, because after a restart the salt will be different. A plagiarism checker that builds up an archive over the years uses a hash function without salt, such as a polynomial hash with a fixed base.

The warning applies to our committee too. If Clara had put the description of the oak in quotation marks and credited Tolstoy, the fingerprints would have found it all the same: a quotation is a shared piece as well. Only a person who reads both texts can tell it from copying.

Back to the lab: all the suffixes at once

We return to the geneticists with one last question. So far the pattern was known in advance and the text was read once. In the lab it’s the other way around: there is one genome, and it doesn’t change, while patterns arrive by the hundreds every day. So it pays to work hard on the text once and build an index that finds any pattern quickly, like the index at the back of a book.

Every occurrence of a pattern is the beginning of some suffix of the text, its tail from some position to the end. Write out all $n$ suffixes and sort them alphabetically, keeping only the positions where they begin. This is the suffix array. All the suffixes that begin with the pattern stand next to one another in it, and binary search from Chapter 20 finds that stretch in $O(m \log n)$, however many patterns come in. Repeats show up in it too: if a piece of text occurs twice, two suffixes begin with it, and sorting puts them side by side. The longest repeat is the longest common beginning of two neighbors in the array.

Sorting the suffixes head-on, sorted(range(n), key=lambda i: s[i:]), works, but every key is a copy of a tail: for the lambda genome that is more than a billion letters in memory. The doubling technique serves better. First order the suffixes by their first letter and give each one a rank, the number of its group. Then by the first two letters: that is a pair of ranks, “the first letter, the letter one position on.” Then by four: a pair of ranks, “the first two letters, the next two.” Each time the length doubles, and after $\log n$ such sorts the order is complete.

Since Python 3.10 the bisect functions accept a key: here it cuts the first six letters off each suffix, and the search compares only those with the pattern. The five EcoRI sites are found again, this time by two binary searches over a ready-made array.

The repeats, though, are short: 12 letters in φX174, 15 each in the mitochondrion and in lambda. Whether that is a lot or a little, the birthday paradox from Chapter 16 will tell us. A random text of $n$ letters has about $n^2/2$ pairs of windows of length $L$, and each pair matches with probability $4^{-L}$. The expected number of matches drops to about one when $4^L \approx n^2/2$, that is, at $L \approx 2\log_4 n$. For φX174 that gives 12.4, for the mitochondrion 14, for lambda 15.6. In these genomes the repeats are no longer than chance would make them: a phage’s genome is thrifty and carries nothing spare.

Our own genome is another matter. About a tenth of it is more than a million copies of a single piece some three hundred letters long. The piece is called Alu, because it was first cut out by the enzyme AluI from our catalog. And near the centromeres, the waists of the chromosomes, repeats stretch for millions of letters. That is why the draft of 2001 could be assembled while the last few percent of the genome had to wait until 2022: a fragment a few hundred letters long from such a place fits a thousand places at once. New instruments that read pieces tens of thousands of letters long helped finish the job.

Tasks

Four tasks: two linear searches, keyboard suggestions and repeats. In each one the starter code already gives the right answers but runs slowly, and the tests with big inputs won’t let it through.

Write prefix_function(p), the list $\pi$ for a sequence p, and find_all(text, pattern), the list of all positions where pattern occurs in text, in increasing order, overlapping occurrences included. The pattern is never empty. The functions must work on lists as well as strings: the tests look, for example, for the tune of “Frère Jacques” in a list of notes, the numbers of a synthesizer’s keys. So the string methods find and in are no help here. The tests include sequences of two to three hundred thousand items, with two seconds for everything.

The starter is right but slow: prefix_function tries every border length at every position, and find_all is the ruler, which on a list of zeros with a one in the middle compares twenty thousand items at each position. Bring over prefix_function from the section on borders: there is nothing string-specific in it, since indices and != work the same for strings and lists.

In find_all, keep one number, k: how many items of the pattern match the end of what has been read. For each item of the text: while k > 0 and the item isn’t equal to pattern[k], fall back with k = pi[k - 1]; if it is equal, k += 1.

When you find an occurrence (k == len(pattern)), don’t reset k to zero: the next occurrence may overlap this one, like aa in aaaa. Fall back along the border: k = pi[k - 1].

Both functions are linear by the coin-counting theorem: k grows by at most one per step, and every fallback lowers it. KMP knows nothing about letters; all it needs is a test for equality. That is why it can find a tune in a list of notes, a sequence in an event log, or a pattern in a stream that arrives one item at a time and doesn’t fit in memory.

Write window_hashes(s, k, base, mod), the list of the polynomial hashes of all windows of length k in the string s, in order. A window’s hash is computed by Horner’s method, like poly_hash in the chapter: start from zero, and for each letter multiply by base, add ord(letter) and take the remainder modulo mod. If the window is longer than the string, the answer is an empty list. Then write rk_find(text, pattern, base, mod), the list of all positions where the pattern occurs. The tests give a string of a million letters with a window of ten thousand and three seconds, and they check rk_find with modulus 7, which raises false alarms at every turn.

The starter has two problems. window_hashes computes the hash of every window from scratch, $n \cdot k$ steps, ten billion on the big test. And rk_find trusts every match of the hashes.

Compute the hash of the first window by Horner’s method, then slide: h = ((h - ord(s[i - 1]) * top) * base + ord(s[i + k - 1])) % mod, where top = pow(base, k - 1, mod) is the weight of the window’s first letter, computed once.

In rk_find a match of the hashes is only a reason to check: compare the slice text[i:i + m] with the pattern.

Sliding the window costs three operations, and a million windows take a fraction of a second, however long the window is. The check by slicing makes the search a Las Vegas algorithm: the answer is always right, and the time depends on luck with the hash. With modulus 7 there are many checks; with $2^{61} - 1$, almost none beyond the true finds. The check if k > len(s) is needed: without it s[:k] would quietly return the whole string, and a short string would get the hash of a window that doesn’t exist.

Write a class Trie with three methods. insert(word) adds one use of a word. count(word) says how many times the word was inserted (zero if never; the beginning of a word doesn’t count as a word). complete(prefix, k) returns at most k words that begin with prefix, the most frequent first and, among equally frequent ones, in alphabetical order; the word prefix itself, if it was inserted, qualifies too. The tests load all 562 thousand words of War and Peace and make twenty thousand queries, with three seconds for the queries.

On every query the starter goes through all twenty thousand different words. Twenty thousand queries make four hundred million startswith checks. A trie walks down the letters of the prefix and looks only into the subtree it needs.

A node is an object with the fields children (a dictionary “letter → node”) and count. The walk down a prefix is needed in both count and complete: move it into a separate method that returns the node or None.

In complete, walk the subtree with a stack, as in the chapter, collecting pairs (-count, word). Sorted in increasing order, such pairs line up by themselves in the order you need: higher frequencies first, and alphabetical among equals.

A query costs the length of the prefix plus the size of the subtree. For prefixes of four to six letters the subtree has a few dozen nodes, and twenty thousand queries fit in a fraction of a second. For a one-letter prefix the subtree is huge, so keyboards keep a ready list of the best continuations in each node and update it on every insertion. You can also do without a trie, with a sorted list of words: all the words with a given prefix stand together in it, and the binary search of Chapter 20 finds that stretch.

Write longest_repeat(s), the longest piece of the string that occurs in it at least twice; the occurrences may overlap, like ana in banana. If there are several such pieces, return any of them; if there are no repeats, return the empty string. The tests take the φX174 genome, a random DNA sequence of a hundred thousand letters with a repeat planted in it, and strings made of nothing but repeats: a telomere and a run of As. The genome gets four seconds, the other big strings five each.

If a piece of length $L$ repeats, so does its beginning of length $L - 1$. So the answer to “is there a repeat of length $L$?” is first “yes,” then “no,” and the boundary can be found by binary search on $L$, in $\log n$ checks instead of $n$.

One check, “is there a repeat of length $L$?”, is a rolling hash over all the windows of length $L$ plus a dictionary “hash → where it was seen.” When two hashes match, compare the pieces themselves, so as not to fall for a false alarm.

The suffix array from the chapter will also do, but carefully: the common beginnings of neighbors, computed by comparing letters, take quadratic time on a string of nothing but As. If you take this road, common beginnings can be computed in linear time with Kasai’s algorithm; look it up yourself.

The binary search on length makes about seventeen checks for a hundred thousand letters, each one a single pass of the rolling hash, $O(n \log n)$ in all. The starter begins with the longest pieces and spends on the order of $n^3$ steps on them for a random string. The test with nothing but As has a different catch: almost everything there is a repeat, and solutions that compare pieces letter by letter, without hashes, turn quadratic. A hash compares a piece of any length in one step, and the slice for the check is needed only once per check.

What next

KMP, fingerprints, tries and suffix arrays all quietly rely on one thing: if letters look the same, the characters of the string are the same. We can test that on the same novel.

In the edition of the novel in the sandbox, Hélène’s name occurs 165 times, but a search for the same name in its other form finds nothing, though both look identical on screen. The difference is that in one string each accented letter is a single character, and in the other it is two: a plain e and a separate accent mark that attaches itself to the letter before it. The name has six characters in one form and eight in the other. The comparison says False, and neither KMP nor a trie nor in will find one Hélène in the other.

The last line of the output shows what the strings are made of. The encode method turns a string into bytes, the stuff that lies in memory and in files. The same six-letter name takes eight bytes in one form and ten in the other: a plain letter takes one byte, and an accented one takes two, or three when the accent is a separate mark. For a whole chapter we have been searching for letters without once asking what a letter is inside the machine. Why does “é” take two numbers and “e” one? What happens if the numbers are read differently from the way they were written? These are questions for the next chapter, the first in the part about how the machine itself works.