Python Arena
Step 5 of 6

Longest common subsequence and edit distance

Cross out some characters of a string, and what's left is a subsequence of it: the characters that stay keep their order, with gaps where others were. "ace" is a subsequence of "abcde". Take two strings, and the longest string you can get from both by crossing out characters is their longest common subsequence (LCS). That's how diff and Git compare two versions of a file, line by line: the lines in the LCS are the ones nobody touched. Trying every subsequence is hopeless, since a string of n characters has 2ⁿ of them. DP works on prefixes instead, the first few characters of each string.

The table

length[i][j] holds the length of the LCS of a[:i] and b[:j], the first i characters of a and the first j of b. Each entry depends on the last characters of those prefixes, a[i - 1] and b[j - 1]:

  • When they match, the LCS of the two prefixes may as well finish with that character, after the LCS of the shorter prefixes: length[i][j] = length[i - 1][j - 1] + 1.
  • When they don't, a common subsequence can't finish with both, so one of the two can go. Try leaving out each and keep the longer result: length[i][j] = max(length[i - 1][j], length[i][j - 1]).

An empty prefix has nothing in common with anything, so row 0 and column 0 are all 0. That's why the table has len(a) + 1 rows and len(b) + 1 columns, and why the characters are a[i - 1] and b[j - 1]: entry i is about the first i characters, and the last of them is at index i - 1.

a, b = "PYTHON", "TYPHOON"
length = [[0] * (len(b) + 1) for _ in range(len(a) + 1)]
for i in range(1, len(a) + 1):
    for j in range(1, len(b) + 1):
        if a[i - 1] == b[j - 1]:
            length[i][j] = length[i - 1][j - 1] + 1
        else:
            length[i][j] = max(length[i - 1][j], length[i][j - 1])

print("   " + "  ".join(b))
for i in range(1, len(a) + 1):
    print(a[i - 1], length[i][1:])
Output
   T  Y  P  H  O  O  N
P [0, 0, 1, 1, 1, 1, 1]
Y [0, 1, 1, 1, 1, 1, 1]
T [1, 1, 1, 1, 1, 1, 1]
H [1, 1, 1, 2, 2, 2, 2]
O [1, 1, 1, 2, 3, 3, 3]
N [1, 1, 1, 2, 3, 3, 4]

The bottom-right entry is the answer: the LCS is 4 characters long. The table has (n + 1) × (m + 1) entries for strings of n and m characters, each O(1): O(n × m) time and memory. For two strings of 1,000 characters, that's a million entries instead of 2¹⁰⁰⁰ subsequences.

The classic slip is to compare a[i] and b[j], out of habit:

a, b = "cat", "act"
length = [[0] * (len(b) + 1) for _ in range(len(a) + 1)]
for i in range(1, len(a) + 1):
    for j in range(1, len(b) + 1):
        if a[i] == b[j]:              # should be a[i - 1] and b[j - 1]
            length[i][j] = length[i - 1][j - 1] + 1
        else:
            length[i][j] = max(length[i - 1][j], length[i][j - 1])
Python stops with
Traceback (most recent call last):
  File "main.py", line 5, in <module>
    if a[i] == b[j]:              # should be a[i - 1] and b[j - 1]
               ~^^^
IndexError: string index out of range

Row 1 is about "c", but it compares a[1], the "a", and at j = 3 it asks for b[3], one past the end of "act".

Getting the string back

The table holds lengths, but it also remembers which choice made each entry, so you can walk back from the bottom-right corner, i = len(a) and j = len(b), until i or j is 0:

  • When a[i - 1] and b[j - 1] are the same character, it's part of the answer. Save it in a list, then step back in both strings: i - 1 and j - 1.
  • Otherwise go to the neighbor the entry came from: up (i - 1) when length[i - 1][j] >= length[i][j - 1], else left (j - 1).

The walk finds the characters from the end backward, so found[::-1] turns the list around. In the PYTHON table it starts at N, where the letters match, and finds "PHON". Two strings often have several longest common subsequences ("THON" and "YHON" are just as long), and a different rule for ties finds a different one.

Edit distance

The edit distance (or Levenshtein distance) between two words is the fewest one-character changes that turn one into the other, where a change inserts, deletes or replaces a character. "brain" becomes "bread" in 3, for example by inserting e after the r, deleting the i and replacing the n with a d. A spell checker offers the words at the smallest distance from your typo.

dist[i][j] holds the edit distance from a[:i] to b[:j]:

  • The edges come first. Emptying i characters takes i deletes, so dist[i][0] = i, and building j characters from nothing takes j inserts, so dist[0][j] = j. Leaving these edges at 0 is the common mistake: the distance from "abc" to "" then comes out 0.
  • When a[i - 1] and b[j - 1] are the same, the ends need no change, and what's left costs dist[i - 1][j - 1].
  • When they're different, one change fixes the end, on top of what came before: a replace of a[i - 1] by b[j - 1] (dist[i - 1][j - 1] + 1), a delete of a[i - 1] (dist[i - 1][j] + 1) or an insert of b[j - 1] (dist[i][j - 1] + 1). dist[i][j] is the smallest of the three.

Here is the table from "cat" (down) to "act" (across):

"" a c t
"" 0 1 2 3
c 1 1 1 2
a 2 1 2 2
t 3 2 2 2

Swapping two neighbors costs 2 edits. Like the LCS, the table costs O(n × m) time and memory, while the plain recursion tries up to three edits at every step: exponential.

Your task

Challenge 1 of 3

This program fills the edit distance table from "an" (down) to "can" (across), and prints it row by row. What does it print? Type each line it prints.

Challenge 2

Write edit_distance(a, b): the smallest number of one-character inserts, deletes and replaces that turn the string a into the string b. Use the lesson's table dist: set its edges, then work out each other entry from the three next to it. The tests have strings of 600 and 450 characters, where trying every edit never ends.

For example, print(edit_distance("brain", "bread"), edit_distance("cat", "act")) prints:

3 2

Challenge 3

Write lcs(a, b), which returns one longest common subsequence of the strings a and b, built with the lesson's table and its walk back, going up on a tie. The tests have strings of 600 and 450 characters, so fill the table without recursion.

For example, lcs("DYNAMIC", "PROGRAMMING") returns:

'AMI'