Analysis on the basic algorithm and Dynamic Programming Algorithm for String Matching
The new phase I of the algorithm class is about DNA matching. DNA can be considered as a string, so it can be converted into a string matching problem. String Matching is not easy, because it is a computing workload and a huge project. For example, the best match between "abc" and "ac" should be "abc" and "a-c", but to calculate such a short string, you must consider many situations.
Therefore, a dynamic planning algorithm is introduced. Dynamic Programming is an algorithm optimized for repeated computing. For example, a problem is split into several small problems and then divided into several small problems. If the calculation of a small problem is repeated, it can be directly substituted. For example, when calculating the matching matrix:
def compute_alignment_matrix(seq_x, seq_y, scoring_matrix, global_flag): length_x = len(seq_x) length_y = len(seq_y) alignment_matrix = [[0 for dummy_col in range(length_y + 1)] for dummy_row in range(length_x + 1)] for idx in range(1, length_x+1): alignment_matrix[idx][0] = alignment_matrix[idx-1][0] + scoring_matrix[seq_x[idx-1]]['-'] for idx in range(1, length_y+1): alignment_matrix[0][idx] = alignment_matrix[0][idx-1] + scoring_matrix['-'][seq_y[idx-1]] for idx_x in range(1, length_x+1): for idx_y in range(1, length_y+1): alignment_matrix[idx_x][idx_y] = max( alignment_matrix[idx_x-1][idx_y-1] + scoring_matrix[seq_x[idx_x-1]][seq_y[idx_y-1]], alignment_matrix[idx_x-1][idx_y] + scoring_matrix[seq_x[idx_x-1]]['-'], alignment_matrix[idx_x][idx_y-1] + scoring_matrix['-'][seq_y[idx_y-1]], ) return alignment_matrix
Scoring_matrix is A score matrix. For example, A and C match A score of-5, A and A match A score of 10, and so on.
This Code uses enumeration to calculate all possible combinations and place the obtained scores on the corresponding positions of the matrix. The previously calculated value can be used for subsequent calculations.
Note:
for idx in range(1, length_x+1): alignment_matrix[idx][0] = alignment_matrix[idx-1][0] + scoring_matrix[seq_x[idx-1]]['-']
For example, the value at the position (1, 0) is actually (0, 0) plus the score at the current position. Then (2, 0) is (1, 0) plus the score of the current position. This greatly reduces the computing workload. Because all sub-problems are solved. The calculation in the y direction is the same as that in the x direction. The calculation in the middle area is beyond the calculation in all directions, and then the maximum value is obtained.
This process is a bit like sliding from the upper left corner to the lower right corner when selecting the image area.
What we get is a score matrix. We select the maximum value from it, which is actually the best match. We can push back the above process to get two matching strings respectively.
def compute_global_alignment(seq_x, seq_y, scoring_matrix, alignment_matrix): score_optimal = max([max(item) for item in alignment_matrix]) for item in alignment_matrix: if score_optimal in item: idx_y = item.index(score_optimal) idx_x = alignment_matrix.index(item) #initial the align_x and align_y, reserve it align_x = seq_x[idx_x:] align_y = seq_y[idx_y:] align_x = align_x[::-1] align_y = align_y[::-1]
while idx_x and idx_y: if alignment_matrix[idx_x][idx_y] == alignment_matrix[idx_x-1][idx_y-1] + scoring_matrix[seq_x[idx_x-1]][seq_y[idx_y-1]]: align_x += seq_x[idx_x-1] align_y += seq_y[idx_y-1] idx_x -= 1 idx_y -= 1 else: if alignment_matrix[idx_x][idx_y] == alignment_matrix[idx_x-1][idx_y] + scoring_matrix[seq_x[idx_x-1]]['-']: align_x += seq_x[idx_x-1] align_y += '-' idx_x -= 1 #import pdb; pdb.set_trace() else: align_x += '-' align_y += seq_y[idx_y-1] idx_y -= 1 while idx_x: align_x += seq_x[idx_x-1] align_y += '-' idx_x -= 1 while idx_y: align_x += '-' align_y += seq_y[idx_y-1] idx_y -= 1 #reserve the string align_x = align_x[::-1] align_y = align_y[::-1] # align_x and align_y have same length length_x = len(align_x) length_y = len(align_y) if length_x < length_y: align_x += '-' * (length_y - length_x) else: align_y += '-' * (length_x - length_y) #The score need to recalculate, cause there may be '-' more score_updated = sum([scoring_matrix[align_x[idx]][align_y[idx]] for idx in range(len(align_x))]) return score_updated, align_x, align_y
It is worth noting that the matched strings are in reverse order, so I reversed the strings twice in this module. This module looks very long. The core of this module is the while loop, which returns the entire process to obtain the matching string.
Before while, first obtain the Index of the maximum and maximum values of the matrix, and then take out the characters that do not enter the while LOOP first.
The characters that do not enter the while loop are the final parts of the match and all remain the original characters.
After a while loop, add the header with deficiencies (the end of the generated header) "-".
Finally, because "-" is added, the score must be recalculated. This completes the global match.
The global match corresponds to a local match, which indicates that the maximum value is formed in the local match. For example, the best match between "tttabcttt" and "ac" is
tttabcttt a-c
The principle of local matching is basically the same as that of global matching, except that when calculating the score matrix, the minimum value is 0 and the first '-' is removed '-'. For more information, see the source code.
The application scope of Dynamic Programming is not as common as iteration, but the efficiency improvement is obvious. String Matching, similarity calculation, and spelling check are all usable. In the last part of this blog series, I will create an algorithm for calculating the similarity between people.
For source code, refer to: Github