Longest common subsequence and edit distance in Python

Why it matters

Comparing two strings is behind a surprising number of everyday tools. A spell checker suggests words that are close to the one you typed. The diff command shows which lines of two files changed. Version control systems find the parts of a document that two people both kept. Bioinformatics compares DNA sequences, which are long strings of four letters, to find regions that match. In each case the program must decide how two strings line up, and the answer comes from a table that compares every prefix of one string with every prefix of the other.

This lesson covers two problems that use that table. The first is the longest common subsequence (LCS), which finds the longest sequence of letters that appears in both strings in the same order. The second is edit distance, which counts the fewest changes that turn one word into the other. Both use a table with one row for each letter of the first string and one column for each letter of the second. Both are standard interview questions, and both teach you how to read a table you have filled in.

You should have finished the earlier lessons in this module. You should know the DP recipe, and you should be comfortable with two-dimensional tables. You should also know how a Python string works. A string holds letters in order, and you read one letter with s[i], where i starts at 0. Every program in this lesson is a plain Python file. Some programs read two lines of input. The lesson tells you exactly what to type.

A subsequence is not a substring

Before the table, fix the vocabulary. A subsequence is a sequence you get by deleting some letters from a string, without changing the order of the rest. You may delete as many letters as you like, including none. The string ace is a subsequence of abcde, because you can delete b and d. The string eca is not a subsequence, because the letters are out of order.

A substring is stricter. It is a block of letters that sits together in the original string, with nothing removed from the middle. The string bcd is a substring of abcde. The string ace is not, because the letters are not next to each other.

This difference decides which algorithm you need. Longest common subsequence allows gaps. Longest common substring does not, and it uses a different recurrence that resets to 0 on a mismatch. Read every problem statement with this difference in mind, because a single word can change the answer.

The longest common subsequence

Take two strings, a = "abcde" and b = "ace". Find the longest string that is a subsequence of both. The letters a, c, e appear in a in that order, and they appear in b in that order too. So ace is a common subsequence of length 3. No common subsequence has length 4, because b has only three letters. The answer is 3.

Think of two people reading two sheets of paper side by side. They walk along each line with a finger, and each time they find a letter that appears on both sheets in the same order, they note it down. They are allowed to skip letters on either sheet. The note they end with is the LCS. The DP table is a systematic way to do that walk for every possible pair of positions at once.

The state

State: dp[i][j] is the length of the LCS of the first i letters of a and the first j letters of b. Both numbers are needed, because the answer depends on how far we have read in each string. Two numbers means a two-dimensional table.

The table has one extra row and one extra column for the empty prefix. dp[0][j] is the LCS of the empty string with the first j letters of b, which is 0. dp[i][0] is the same, with the roles swapped. The empty prefix is a real case, and it is the base of the table.

The recurrence

To fill cell dp[i][j], look at the last letter of each prefix: a[i - 1] and b[j - 1]. Note the - 1, because the prefix of length i ends at index i - 1. Python counts from 0, just as Java does.

• The last letters match. Then that letter can end a common subsequence. Add it to the best answer for the two shorter prefixes: dp[i][j] = dp[i - 1][j - 1] + 1. Taking a matched pair never hurts, so this is the right choice.
• The last letters differ. Then they cannot both end the subsequence. Drop the last letter of one string, or the other, and keep the better answer: dp[i][j] = max(dp[i - 1][j], dp[i][j - 1]).

Order: row by row, from top to bottom, and within each row from left to right. Each cell reads the cell above it, the cell to its left, and the cell diagonally above and to the left. All three are earlier in the order, so they are filled.

A complete trace for abcde and ace

Here is the full table. Rows are the letters of a, with an empty row at the top. Columns are the letters of b, with an empty column at the left. The top row and left column are all 0.

Check three cells by hand. Cell dp[1][1] compares a with a. They match, so the cell is dp[0][0] + 1 = 1. Cell dp[3][2] compares c (the third letter of a) with c (the second letter of b). They match, so it is dp[2][1] + 1 = 1 + 1 = 2. Cell dp[5][3] compares e with e. They match, so it is dp[4][2] + 1 = 2 + 1 = 3, which is the answer.

The code for LCS

Reading the input needs one small point. Python's input() returns the next line as a string, without the newline at the end. If a line is blank, input() returns an empty string, which is what we want when the problem allows empty words. You do not need to import anything for this.

The program below reads two lines. Type these two lines when it asks for input:

abcde
ace

def lcs(a, b):
    m, n = len(a), len(b)
    dp = [[0] * (n + 1) for _ in range(m + 1)]
    for i in range(1, m + 1):
        for j in range(1, n + 1):
            if a[i - 1] == b[j - 1]:
                dp[i][j] = dp[i - 1][j - 1] + 1
            else:
                dp[i][j] = max(dp[i - 1][j], dp[i][j - 1])
    return dp[m][n]
 
 
a = input()
b = input()
print("LCS length:", lcs(a, b))

The output is:

LCS length: 3

Look at the loops. The table has m + 1 rows and n + 1 columns, so the indexes i and j run from 1 to m and 1 to n, and the row 0 and column 0 keep their zeros. Inside the loop, the letters are read with a[i - 1]. Reading a[i] is the most common bug in this kind of code, and the common mistakes section shows its error.

Recovering the common letters

Knowing the length is often not enough. Interviewers usually ask for the actual common subsequence, which means the letters themselves. To get them, walk back from the bottom-right cell, the same way you walked back through the knapsack table in the last lesson. At each cell, look at the letters.

• If the letters match, that letter is in the LCS. Record it, then move diagonally up and to the left, because both prefixes shrink by one.
• If the letters do not match, move to whichever neighbour, up or left, has the larger value. If they are equal, either one works. Move there, and do not record a letter.

The recorded letters come out in reverse order, because the walk starts at the end. Reverse them at the end.

Here is the walk for AGGTAB and GXTXAYB, which is the classic example. The length is 4. Start at the bottom-right corner with prefixes of length 6 and 7. The last letters are B and B, which match, so record B and move diagonally. Now compare A with Y. They differ, and the cell to the left has the larger value, so move there. Here A and A match, so record A. Continue this process until you reach the top-left. Recorded in order: B, A, T, G. Reversed, the common subsequence is GTAB.

The program below does this walk. It also prints the length from the table, so you can check that the string it returns is as long as the number says. It does not read input.

The output is:

length: 4
letters: GTAB

The walk checks the match first. If the letters match, it must take the diagonal, because taking a match is always part of some longest answer. Only when the letters differ does the code choose a direction. The >= in the comparison means that when the two neighbours tie, the walk goes up. Either choice gives a correct common subsequence, so the tie rule does not change correctness, but it does change which of the valid answers you get.

The list found collects the letters as the walk adds them. Each append adds one letter to the end. The call reversed(found) walks the list backward, and "".join(...) glues the letters into one string. Joining once at the end is faster than adding strings with + inside a loop, because each + creates a new string.

Edit distance: counting the changes

The second problem measures how different two words are. The edit distance between two words is the fewest single-letter changes that turn one word into the other. The allowed changes are three.

• Insert a letter.
• Delete a letter.
• Replace one letter with another.

For example, horse becomes ros in three changes. Replace h with r, giving rorse. Delete the second r, giving rose. Delete the final e, giving ros. No sequence of two changes does it, so the edit distance is 3. Spell checkers use this measure to rank suggestions: a word one change away is a better suggestion than a word three changes away.

The state

State: dp[i][j] is the edit distance between the first i letters of a and the first j letters of b. Again the table is two-dimensional, and again it has an extra row and column for the empty prefix.

Base cases: the edit distance from a prefix to the empty string is the number of letters in the prefix, because you must delete each one. So dp[i][0] = i. The edit distance from the empty string to a prefix is the number of letters in the prefix, because you must insert each one. So dp[0][j] = j. These base cases are not 0, which is the most common surprise for beginners.

The recurrence

Look at the last letters again, a[i - 1] and b[j - 1].

• If they match, the last letters need no change. The problem for the shorter prefixes is already solved, so dp[i][j] = dp[i - 1][j - 1].
• If they differ, one change is needed, and there are three ways to make it. Each one leads to a smaller problem. We take the cheapest of the three and add 1 for the change.
  • Replace the last letter of a with the last letter of b. The prefixes shrink by one each. This uses dp[i - 1][j - 1].
  • Delete the last letter of a. Then a shrinks and b does not. This uses dp[i - 1][j].
  • Insert the last letter of b into a. Then b shrinks and a does not. This uses dp[i][j - 1].

So for a mismatch: dp[i][j] = 1 + min(dp[i - 1][j - 1], dp[i - 1][j], dp[i][j - 1]).

Order: row by row, left to right, the same as LCS.

A complete trace for horse and ros

Here is the full table for horse (rows) and ros (columns). The first row and column are the base cases, counting up from 0.

Check the cell dp[1][1]. The letters are h and r, which differ. The three options are: replace (dp[0][0] = 0), delete (dp[0][1] = 1), and insert (dp[1][0] = 1). The smallest is 0, so the cell is 1 + 0 = 1. Check dp[2][2]. The letters are o and o, which match, so the cell is dp[1][1] = 1. Check dp[5][3]. The letters are e and s, which differ. The options are replace (dp[4][2] = 3), delete (dp[4][3] = 2), and insert (dp[5][2] = 4). The smallest is 2, so the cell is 1 + 2 = 3, the answer.

The table also shows where the three edits come from. Walk back from the bottom-right corner. Each step tells you which operation was used: a diagonal step is a replace or a match, a step up is a delete, and a step left is an insert. The walk from dp[5][3] goes up to dp[4][3] (delete), then diagonal to dp[3][2] (a letter pair that matches or is replaced), and so on.

The code for edit distance

The next program reads two words, one per line, and prints the edit distance. Empty lines are valid input here, so the program handles them correctly. Type these two lines when it asks for input:

horse
ros

def edit_distance(a, b):
    m, n = len(a), len(b)
    dp = [[0] * (n + 1) for _ in range(m + 1)]
    for i in range(m + 1):
        dp[i][0] = i
    for j in range(n + 1):
        dp[0][j] = j
    for i in range(1, m + 1):
        for j in range(1, n + 1):
            if a[i - 1] == b[j - 1]:
                dp[i][j] = dp[i - 1][j - 1]
            else:
                dp[i][j] = 1 + min(dp[i - 1][j - 1], dp[i - 1][j], dp[i][j - 1])
    return dp[m][n]
 
 
a = input()
b = input()
print("edit distance:", edit_distance(a, b))

The output is:

edit distance: 3

The function min in Python takes any number of arguments, so min(x, y, z) picks the smallest of three. In Java the same call needs two nested calls, which is a common source of mistakes there. Forgetting one of the three options gives a wrong answer that is sometimes right. The common mistakes section shows how to spot it.

Try an empty first line. Replace the input with a blank line followed by abc. The output is edit distance: 3, because turning the empty string into abc takes three inserts. The base case dp[0][3] = 3 handles it.

A batch of edit distances

The next program does not read input. It runs the function on a set of word pairs, including the empty string and some well-known examples. It also prints the whole table for horse and ros, so you can compare it with the trace above. Study the results, because a good test set is the fastest way to find the edge cases before an interviewer does.

The output is:

[0, 1, 2, 3]
[1, 1, 2, 3]
[2, 2, 1, 2]
[3, 2, 2, 2]
[4, 3, 3, 2]
[5, 4, 4, 3]
"" -> "abc": 3
"abc" -> "abc": 0
"cat" -> "cut": 1
"kitten" -> "sitting": 3
"sunday" -> "saturday": 3
"intention" -> "execution": 5

Each row teaches something. The empty string needs one insert per letter. Identical words need no changes. cat and cut differ by one replacement. kitten and sitting need two replacements and one insert. sunday and saturday need three changes. The pair intention and execution is a classic from textbooks, and its distance of 5 is worth checking by hand once.

The f-string f'"{a}" -> "{b}": {edit_distance(a, b)}' uses single quotes outside, so the double quotes inside are plain characters. In a normal string you would write \" to put a double quote inside. Either way the printed text has the quotes.

What happens inside: memory for a two-dimensional table

A two-dimensional list in Python is a list of lists. [[0] * (n + 1) for _ in range(m + 1)] creates one outer list with m + 1 slots, and each slot holds a reference to an inner list with n + 1 integers. Reading dp[i][j] first finds the inner list at position i, then reads position j inside it.

Walking row by row touches each inner list in order, which is kind to the processor's cache, the fast memory that holds recently used data. Walking column by column jumps between lists and is slower on large inputs. Our loops go row by row, so they are the fast order.

The memory for the full table is about (m + 1) * (n + 1) integers. Each Python integer takes more space than a Java int, because a Python int is an object that can grow. For two words of 1000 letters each, the table has about a million cells, which is a few tens of megabytes in Python. Fine. For two strings of a million letters each, the table would have a trillion cells, which is impossible. For such inputs, use a different algorithm, or accept an approximate answer. Interviewers may ask about this, and you should be able to say it in one sentence.

Shrinking to two rows

Each row of the table reads only the row above it and the cell to its left. So you can keep just two rows, the current one and the previous one. The memory drops from m * n to 2 * n, which is O(n) in the length of the second word. The table in the code is easier to read, so the lessons keep it. Mention the two-row version in an interview as a memory improvement. For edit distance, the diagonal cell must be saved in a variable before the current row is overwritten, which is a small but common mistake.

Reading a finished table

A finished table holds more than the answer. Read it the way you would read a map, looking for patterns that tell you whether the numbers make sense.

The answer always sits in the bottom-right corner. Start there when you check a result, and walk back through the table to see how the answer was built. The walk shows you which letters were matched (for LCS) or which operations were used (for edit distance). This is the same walk the recovery programs use, so reading a table by hand and writing the recovery code are one skill.

Two sanity checks catch most wrong tables. First, in the LCS table, no value ever goes down as you move right or down the table. A longer prefix cannot share a shorter common subsequence than a shorter prefix, so each row and each column only grows or stays the same. In the abcde and ace table, each row starts with 0 and then only goes up or stays level, like 0, 1, 2, 2. If you see a value drop, the table is wrong.

Second, in both tables, two neighbouring cells differ by at most 1. Adding one letter to a string can raise the LCS by at most 1, and it can change the edit distance by at most 1, because one extra change is enough to fix one extra letter. If two neighbours differ by 2 or more, a cell was computed from the wrong inputs. Find it by checking the row above and the cell to the left.

Edit distance does not have the first check. Its values can go down as you move through the table, as the horse and ros table shows in the column for o, where the value drops from 2 to 1. Use only the neighbour check there.

Common mistakes

Mistake 1: using a[i] instead of a[i - 1]. The table has an extra empty row and column, so the prefix of length i ends at index i - 1. If the code reads a[i], it compares the wrong letters, and at the last row it reads past the end of the string. Here is what happens when the code uses a[i] with the string abc. The loop reaches i = 3, and a[3] asks for a letter that does not exist:

a = "abc"
for i in range(1, len(a) + 1):
    print(a[i])

The program prints b and c, and then stops with this error:

Traceback (most recent call last):
  File "main.py", line 3, in <module>
    print(a[i])
          ~^^^
IndexError: string index out of range

Read the last line: index 3 is not allowed for a string of length 3, because the valid indexes are 0 to 2. The line above it shows the expression that failed. The fix is a[i - 1].

Mistake 2: base cases of 0 for edit distance. Writing dp[i][0] = 0 is wrong. Turning abc into the empty string takes three deletes, not zero. A table that starts with zeros will report a distance of 0 for anything against the empty string, which is clearly false.

Mistake 3: taking the minimum of only two options. Edit distance has three options: replace, delete and insert. Leave one out and the program still prints a number, but the number is sometimes too big. For kitten and sitting the correct answer is 3, and it uses two replacements and one insert. Without the replace option the same table gives 5, and without the insert option it gives 6. A test with the expected answer is the only way to catch this, which is why the batch program above includes that pair.

Mistake 4: confusing LCS with the longest common substring. The substring version resets to 0 on a mismatch. The subsequence version keeps the best answer from the neighbours. If you write dp[i][j] = 0 on a mismatch for LCS, you compute the wrong quantity. The next program does this, and it prints 1 for abcde and ace, where the correct LCS length is 3. The program prints the largest value in the table, which is the longest common substring of the two strings.

The output is 1. The longest block of letters that sits together in both strings is a single letter.

Mistake 5: an empty line that is not read as empty. Python's input() returns the empty string for a blank line, so the program in this lesson handles empty words without any extra step. Some other ways of reading input skip blank lines and read the next word instead. If your answer on a blank input seems to come from the wrong word, check how the program reads its lines. In a Python file, input() is the reliable choice for this lesson.

Mistake 6: a missing line in the input. If the program asks for a second line and the input has none, Python stops with an error. The program below reads two lines, but the input has only one.

a = input()
b = input()
print(a, b)

With the input abc and then nothing, the last line of the error is:

EOFError: EOF when reading a line

The message means the program asked for another line, and the input had none left. Python does not check this in advance, so you find out only when the program runs. Give the input the number of lines the program expects.

Debugging a string table

When a string answer is wrong, print the table, exactly as you did for knapsack. For abcde and ace, the printed table must match the one in this lesson. Compare it row by row. The first differing cell tells you which pair of letters was handled wrongly.

For edit distance, a useful debugging trick is to check the base row and column first. They must be 0, 1, 2, 3 across the top and 0, 1, 2, 3, 4, 5 down the side for horse and ros. If they are wrong, nothing after them can be right. Print only the first row and first column. That single check catches most base-case bugs.

Another trick is to test a pair that you know by heart. cat and cut differ by one letter, so the answer must be 1. If the answer is 2, look at the replace option first, because the single change is a replacement.

Edge cases you should test

• Both strings empty: LCS is 0, edit distance is 0.
• One string empty: LCS is 0, edit distance is the length of the other string.
• Identical strings: LCS is the full length, edit distance is 0.
• Completely different letters, such as abc and xyz: LCS is 0, edit distance is 3.
• Repeated letters, such as aaaa and aa: LCS is 2, edit distance is 2.

Write these tests before you code, and check each by hand against the table. A test set that covers the empty string and a pair with no common letters catches almost every base-case bug.

Interview questions you should expect

"What is the time and space complexity of LCS?" Time is O(m * n), because the table has that many cells and each takes constant work. Space is O(m * n) for the full table, or O(n) with the two-row version.

"How would you return the actual LCS string?" Keep the full table and walk back from the bottom-right corner, as the recovery program does. The two-row version cannot do this, because the earlier rows are gone. That trade-off is a strong point to mention.

"Why are the base cases for edit distance i and j, and not 0?" Because turning a prefix into the empty string takes one delete per letter, and turning the empty string into a prefix takes one insert per letter. The base case is the cost of the obvious route.

"What is the difference between LCS and edit distance?" LCS counts matches and ignores differences. Edit distance counts the changes needed and ignores matches. They can be related: for the insert and delete operations alone, the number of edits equals the total length minus twice the LCS. Give an example if the interviewer asks for proof, such as abc and abd, where the LCS is ab. With only inserts and deletes, the count is 3 + 3 - 2 * 2 = 2, which is one delete and one insert. With replace allowed, the edit distance is 1, because one replacement turns abc into abd.

"How is this used in practice?" Spell checkers rank suggestions by edit distance. The diff tool uses an LCS-style method to show which lines changed. Version control and DNA alignment use the same table with different scores for each operation.

Practice

Try these before reading the solutions. Work each one by hand first, then check it with the programs.

Exercise 1 (easy). Find the LCS length of abc and abc. Then find the LCS length of abc and def.

Show the answer

abc and abc share all three letters in order, so the LCS length is 3. abc and def share no letter, so the LCS length is 0.

Exercise 2 (easy). Find the edit distance from cat to cut, and from cat to at.

Show the answer

cat to cut replaces a with u, which is one change. cat to at deletes the c, which is one change. Both answers are 1.

Exercise 3 (easy). Write the base row and base column of the edit distance table for abc (rows) and xy (columns).

Show the answer

The top row is 0, 1, 2 for the empty prefix and the prefixes of xy. The left column is 0, 1, 2, 3 for the empty prefix and the prefixes of abc. So dp[0][0] = 0, dp[0][1] = 1, dp[0][2] = 2, dp[1][0] = 1, dp[2][0] = 2, and dp[3][0] = 3.

Exercise 4 (medium). Fill the full LCS table for ab (rows) and ba (columns). What is the LCS length?

Show the answer

The table is [[0, 0, 0], [0, 0, 1], [0, 1, 1]]. The length is 1. Either a or b is a common subsequence, but no string of two letters is common, because the letters are in opposite order in the two words.

Exercise 5 (medium). Using the walk-back rule, recover a common subsequence of abcd and acbd.

Show the answer

The LCS length is 3, and one answer is abd. Walk back from dp[4][4]. The last letters are d and d, which match, so record d. Then compare c with b. They differ, and the cell above is dp[2][3], which is 2, and the cell to the left is dp[3][2], which is 2. Take the up move on a tie, to dp[2][3], where the letters are b and b. Record b. Then a and c differ. The cell above holds 0 and the cell to the left holds 1, so move left, where a and a match, so record a. The recorded letters in order are a, b, d, so the common subsequence is abd. A quick run of the walk-back function on these two strings prints (3, 'abd').

Exercise 6 (medium). Explain why dp[i][j] = dp[i - 1][j - 1] is correct when the last letters of the prefixes match in edit distance.

Show the answer

When the last letters match, the best way to fix the two prefixes is to fix the shorter prefixes and leave the matched letter alone. Adding a change for a matched letter would only add cost, and no other option is cheaper. So the cost equals the cost for the two prefixes without their last letters.

Exercise 7 (medium). Count the edit operations for changing sunday into saturday by listing them in order.

Show the answer

One route with three changes: insert a after the first s to get saunday, insert t after that a to get satunday, then replace the n with r to get saturday. The answer is 3. Each change is one of the three allowed operations. Two changes are not enough: the words differ in length by 2, so both changes would have to be inserts, and then sunday would have to appear inside saturday in order, which it does not, because saturday has no n.

Exercise 8 (hard). Explain why the LCS table fills correctly row by row, and what would go wrong if you filled it column by column but kept the same formula.

Show the answer

Filling column by column still works, because the formula reads only the cell above, the cell to the left and the diagonal cell, and all three are filled before the current cell in either order. The order only matters if the code reads a cell that is not yet filled. The formula never reads a cell further along, so either row-by-row or column-by-column order is correct. What would fail is a formula that reads the cell below or to the right.

Exercise 9 (hard). Show that the edit distance can never be larger than max(m, n), where m and n are the lengths of the two words.

Show the answer

Use replace operations for the first min(m, n) letters, then insert or delete the rest. That uses at most min(m, n) replacements and |m - n| more operations. The total is min(m, n) + |m - n| = max(m, n). So the true distance is at most max(m, n), and the table, which finds the smallest cost, can never report more.

What you now know

• A subsequence keeps the order of letters and may skip some. A substring must be continuous.
• The LCS table compares the last letters of two prefixes. A match adds 1 to the diagonal cell, and a mismatch takes the larger neighbour.
• To recover the LCS letters, walk back from the bottom-right corner and record the matches.
• Edit distance has three operations: insert, delete and replace. Its base cases are the lengths of the prefixes, not zeros.
• Use a[i - 1] for the prefix of length i, and check that the input has as many lines as the program reads.
• Print the table and check the base row and column first when an answer is wrong.