Gaps in one piece
Proteins are compared by aligning them: two amino acid sequences are written one under the other with gaps, dashes, inserted into them. A column of two letters x and y scores S(x, y) from a substitution matrix. One mutation often inserts or deletes several amino acids at once, so a gap is a run of consecutive dashes in one row, and a gap of length k costs O + E·k: expensive to open, cheap to extend. Two runs of dashes back to back in different rows are two different gaps.
Align two proteins end to end, from the first letter to the last, with the largest score: the sum of the scores of two-letter columns minus the gap penalties.
The proteins in the tests are haemoglobins, myoglobins, lysozyme, ribonuclease and proteins of phages λ, HK97 and HK022; the matrix is almost always BLOSUM62.
Input
The first line holds the integers O and E (0 ≤ O ≤ 50, 1 ≤ E ≤ 20). The second holds the 20 amino acid letters separated by spaces. Then 20 lines of the substitution matrix: a letter and 20 integers from −20 to 20, the scores of its pairs with the letters in the order of the second line; the matrix is symmetric. The last two lines hold the proteins: strings of these letters from 1 to 600 long.
Output
The first line holds the largest score. The second and third hold the alignment: the first and second protein, of equal length, with dashes in place of gaps. There are no columns of two dashes. If there are several best alignments, print any.
Examples
Example 1
3 1 A R N D C Q E G H I L K M F P S T W Y V A 4 -1 -2 -2 0 -1 -1 0 -2 -1 -1 -1 -1 -2 -1 1 0 -3 -2 0 R -1 5 0 -2 -3 1 0 -2 0 -3 -2 2 -1 -3 -2 -1 -1 -3 -2 -3 N -2 0 6 1 -3 0 0 0 1 -3 -3 0 -2 -3 -2 1 0 -4 -2 -3 D -2 -2 1 6 -3 0 2 -1 -1 -3 -4 -1 -3 -3 -1 0 -1 -4 -3 -3 C 0 -3 -3 -3 9 -3 -4 -3 -3 -1 -1 -3 -1 -2 -3 -1 -1 -2 -2 -1 Q -1 1 0 0 -3 5 2 -2 0 -3 -2 1 0 -3 -1 0 -1 -2 -1 -2 E -1 0 0 2 -4 2 5 -2 0 -3 -3 1 -2 -3 -1 0 -1 -3 -2 -2 G 0 -2 0 -1 -3 -2 -2 6 -2 -4 -4 -2 -3 -3 -2 0 -2 -2 -3 -3 H -2 0 1 -1 -3 0 0 -2 8 -3 -3 -1 -2 -1 -2 -1 -2 -2 2 -3 I -1 -3 -3 -3 -1 -3 -3 -4 -3 4 2 -3 1 0 -3 -2 -1 -3 -1 3 L -1 -2 -3 -4 -1 -2 -3 -4 -3 2 4 -2 2 0 -3 -2 -1 -2 -1 1 K -1 2 0 -1 -3 1 1 -2 -1 -3 -2 5 -1 -3 -1 0 -1 -3 -2 -2 M -1 -1 -2 -3 -1 0 -2 -3 -2 1 2 -1 5 0 -2 -1 -1 -1 -1 1 F -2 -3 -3 -3 -2 -3 -3 -3 -1 0 0 -3 0 6 -4 -2 -2 1 3 -1 P -1 -2 -2 -1 -3 -1 -1 -2 -2 -3 -3 -1 -2 -4 7 -1 -1 -4 -3 -2 S 1 -1 1 0 -1 0 0 0 -1 -2 -2 0 -1 -2 -1 4 1 -3 -2 -2 T 0 -1 0 -1 -1 -1 -1 -2 -2 -1 -1 -1 -1 -2 -1 1 5 -2 -2 0 W -3 -3 -4 -4 -2 -2 -3 -2 -2 -3 -2 -3 -1 1 -4 -3 -2 11 2 -3 Y -2 -2 -2 -3 -2 -1 -2 -3 2 -1 -1 -2 -1 3 -3 -2 -2 2 7 -1 V 0 -3 -3 -3 -1 -2 -2 -3 -3 3 1 -2 1 -1 -2 -2 0 -3 -1 4 HEAGAWGHEE PAWHEAE
17 HEAGAWGHE-E ---PAW-HEAE
Example 2
10 1 A R N D C Q E G H I L K M F P S T W Y V A 4 -1 -2 -2 0 -1 -1 0 -2 -1 -1 -1 -1 -2 -1 1 0 -3 -2 0 R -1 5 0 -2 -3 1 0 -2 0 -3 -2 2 -1 -3 -2 -1 -1 -3 -2 -3 N -2 0 6 1 -3 0 0 0 1 -3 -3 0 -2 -3 -2 1 0 -4 -2 -3 D -2 -2 1 6 -3 0 2 -1 -1 -3 -4 -1 -3 -3 -1 0 -1 -4 -3 -3 C 0 -3 -3 -3 9 -3 -4 -3 -3 -1 -1 -3 -1 -2 -3 -1 -1 -2 -2 -1 Q -1 1 0 0 -3 5 2 -2 0 -3 -2 1 0 -3 -1 0 -1 -2 -1 -2 E -1 0 0 2 -4 2 5 -2 0 -3 -3 1 -2 -3 -1 0 -1 -3 -2 -2 G 0 -2 0 -1 -3 -2 -2 6 -2 -4 -4 -2 -3 -3 -2 0 -2 -2 -3 -3 H -2 0 1 -1 -3 0 0 -2 8 -3 -3 -1 -2 -1 -2 -1 -2 -2 2 -3 I -1 -3 -3 -3 -1 -3 -3 -4 -3 4 2 -3 1 0 -3 -2 -1 -3 -1 3 L -1 -2 -3 -4 -1 -2 -3 -4 -3 2 4 -2 2 0 -3 -2 -1 -2 -1 1 K -1 2 0 -1 -3 1 1 -2 -1 -3 -2 5 -1 -3 -1 0 -1 -3 -2 -2 M -1 -1 -2 -3 -1 0 -2 -3 -2 1 2 -1 5 0 -2 -1 -1 -1 -1 1 F -2 -3 -3 -3 -2 -3 -3 -3 -1 0 0 -3 0 6 -4 -2 -2 1 3 -1 P -1 -2 -2 -1 -3 -1 -1 -2 -2 -3 -3 -1 -2 -4 7 -1 -1 -4 -3 -2 S 1 -1 1 0 -1 0 0 0 -1 -2 -2 0 -1 -2 -1 4 1 -3 -2 -2 T 0 -1 0 -1 -1 -1 -1 -2 -2 -1 -1 -1 -1 -2 -1 1 5 -2 -2 0 W -3 -3 -4 -4 -2 -2 -3 -2 -2 -3 -2 -3 -1 1 -4 -3 -2 11 2 -3 Y -2 -2 -2 -3 -2 -1 -2 -3 2 -1 -1 -2 -1 3 -3 -2 -2 2 7 -1 V 0 -3 -3 -3 -1 -2 -2 -3 -3 3 1 -2 1 -1 -2 -2 0 -3 -1 4 YFPHFDLSHGSAQVKGH FFESFGDLSTPDAVMGNPKVKAH
23 YFPHF-DLS-----HGSAQVKGH FFESFGDLSTPDAVMGNPKVKAH
Limits
- Time · Python
- 4000 ms
- Time · Go
- 1000 ms
- Memory
- 256 MB
- Judging
- 16 tests
Walkthrough
Open it step by step: each step gives away a little more. You can jump straight to the solution, but a problem solved with one hint stays with you longer.
1Hint 1
If every dash cost the same, one table would be enough, as for the edit distance. What do you need to know about an alignment of the first letters to tell what the next dash will cost?
2Hint 2
The next dash costs E if the last column is a dash in the same row, and O + E otherwise. So it matters which column an alignment ends with. Keep a table for each kind of last column: a pair of letters, a letter of the first protein over a dash, a dash over a letter of the second.
3Hint 3
The edges: aligning the first i letters with an empty string is one gap of length i, scoring −O − E·i, and it ends with a dash in the second row. Everything else on the edges is impossible: minus infinity, but one from which subtraction does not overflow.
4The idea
Let M(i, j), X(i, j), Y(i, j) be the best scores of alignments of the first i letters of protein a with the first j letters of protein b that end in a pair of letters, in the letter ai over a dash, and in a dash over bj. Then
M(i, j) = S(ai, bj) + max(M, X, Y at (i − 1, j − 1)),
X(i, j) = max(M(i − 1, j) − O − E, X(i − 1, j) − E, Y(i − 1, j) − O − E),
Y(i, j) = max(M(i, j − 1) − O − E, Y(i, j − 1) − E, X(i, j − 1) − O − E).
Why this is right. Remove the last column from the best alignment ending in a letter over a dash. The rest is an alignment of the first i − 1 letters with the first j, and the last column cost E if the rest also ends in a letter over a dash, and O + E otherwise. The cost of the last column depends on the rest only through the kind of its last column, so the rest must be the best among alignments of its kind, otherwise it could be replaced by a better one of the same kind and the whole alignment improved. Trying the three kinds of rest gives the formula. The argument for the other two tables is the same. The answer is the largest of the three values in cell (n, m).
5Walkthrough
Fill the three tables row by row, then walk back from cell (n, m), remembering which table you are in. In M write out a pair of letters and find which of the three tables in the diagonal cell gave the maximum; in X write out a letter over a dash and check whether it is a continuation of the gap or its opening, and so on. Reverse what was written out.
The first example shows it at once. If each dash costs 4, it pays to scatter them: --P-AW-HEAE. With an opening of 3 and extension of 1, three dashes in a row cost 3 + 3 = 6, while the same three dashes in two pieces cost 5 + 4 = 9, and the best alignment ---PAW-HEAE scores 17. In the second example there is a corner of the haemoglobin α and β chains: five dashes in a row stand opposite the β chain’s helix D, which the α chain lacks.
5+Solution
6Complexity
Time and memory are O(n·m): three tables of 361,201 cells in the worst test. If only the score is needed, two rows of each table are enough. The alignment too can be obtained in linear memory, by the method of Myers and Miller: compute the tables from both ends to the middle row, find where the best path crosses it, and solve the two halves separately.
7Common mistakes
- Opening for every dash. One table with a penalty of O + E per letter is a linear penalty, and a long gap costs as much as as many short ones. In the first example this gives an alignment scoring 15, in the second 14 instead of 23.
- Back to back means one piece. If you move from X to Y without an opening penalty, a dash over a letter followed right away by a letter over a dash merge into one gap, while by the statement they are two.
- Edges without opening. A gap at the start of an alignment is a gap too: −O − E·i, not −E·i.
- A traceback without state. If in each cell you take whichever table has the larger value, you can jump onto a path that came into the cell differently and write out an alignment with a different score. You must remember which table you are in.
- Overflow. A “minus infinity” such as the smallest 32-bit number becomes a huge positive one after subtracting a penalty. Leave a margin: −260 in 64 bits.