Реализация алгоритма Смита-Уотермана для локального выравнивания в python

Я создал инструмент выравнивания последовательностей для сравнения двух цепочек ДНК (X и Y), чтобы найти наилучшее выравнивание подстрок из X и Y. Алгоритм кратко изложен здесь (https://en.wikipedia.org)./wiki/Smith%E2%80%93Waterman_algorithm). Я смог создать списки списков, заполнив их нулями, чтобы представить мою матрицу. Я создал алгоритм подсчета очков, чтобы возвращать числовой балл для каждого типа выравнивания между базами (например, плюс 4 для совпадения). Затем я создал алгоритм выравнивания, который должен ставить оценку каждой координате моей матрицы. Однако, когда я иду печатать матрицу, она возвращает только оригинал со всеми нулями (а не фактические баллы).

Я знаю, что есть другие методы реализации этого метода (например, с помощью numpy), поэтому не могли бы вы сказать мне, почему этот конкретный код (ниже) не работает? Есть ли способ изменить его, чтобы он работал?

код:

def zeros(X: int, Y: int):
    lenX = len(X) + 1
    lenY = len(Y) + 1
    matrix = []
    for i in range(lenX):
        matrix.append([0] * lenY)
    
    def score(X, Y):
        if X[n] == Y[m]: return 4
        if X[n] == '-' or Y[m] == '-': return -4
        else: return -2

    def SmithWaterman(X, Y, score):
        for n in range(1, len(X) + 1):
            for m in range(1, len(Y) + 1):
                align = matrix[n-1, m-1] + (score(X[n-1], Y[m-1]))
                indelX = matrix[n-1, m] + (score(X[n-1], Y[m]))
                indelY = matrix[n, m-1] + (score(X[n], Y[m-1]))
        matrix[n, m] = max(align, indelX, indelY, 0)
    print(matrix)

zeros("ACGT", "ACGT")

выход:

[[0, 0, 0, 0, 0], [0, 0, 0, 0, 0], [0, 0, 0, 0, 0], [0, 0, 0, 0, 0], [0, 0, 0, 0, 0]]

python algorithm sequence-alignment dna-sequence
person LampoonBonanza    schedule 25.02.2021    source источник

Ответы (1)


arrow_upward
0
arrow_downward

Причина, по которой он просто выводит обнуленную матрицу, заключается в том, что функция SmithWaterman никогда не вызывается, поэтому матрица никогда не обновляется.

Вам нужно будет сделать что-то вроде

# ...
SmithWaterman(X, Y, score)
print(matrix)
# ...

Однако, если вы сделаете это, вы обнаружите, что этот код на самом деле не работает во многих других отношениях. Я просмотрел и аннотировал некоторые синтаксические ошибки и другие проблемы с кодом:

def zeros(X: int, Y: int):
    #        ^       ^  incorrect type annotations. should be str
    lenX = len(X) + 1
    lenY = len(Y) + 1
    matrix = []
    for i in range(lenX):
        matrix.append([0] * lenY)
    # A more "pythonic" way of expressing the above would be:
    # matrix = [[0] * len(Y) + 1 for _ in range(len(x) + 1)]
    
    def score(X, Y):
        #     ^  ^ shadowing variables from outer scope. this is not a bug per se but it's considered bad practice
        if X[n] == Y[m]: return 4
        #    ^       ^  variables not defined in scope
        if X[n] == '-' or Y[m] == '-': return -4
        #    ^              ^  variables not defined in scope
        else: return -2

    def SmithWaterman(X, Y, score): # this function is never called
        #                   ^ unnecessary function passed as parameter. function is defined in scope
        for n in range(1, len(X) + 1):
            for m in range(1, len(Y) + 1):
                align = matrix[n-1, m-1] + (score(X[n-1], Y[m-1]))
                #                 ^ invalid list lookup. should be: matrix[n-1][m-1]
                indelX = matrix[n-1, m] + (score(X[n-1], Y[m]))
                #                                          ^ out of bounds error when m == len(Y)
                indelY = matrix[n, m-1] + (score(X[n], Y[m-1]))
                #                                  ^ out of bounds error when n == len(X)
        matrix[n, m] = max(align, indelX, indelY, 0)
        # this should be nested in the inner for-loop. m, n, indelX, and indelY are not defined in scope here
    print(matrix)

zeros("ACGT", "ACGT")
person George D.    schedule 26.02.2021
comment
Я исправил это с помощью ваших комментариев и вместо этого добавил print(*matrix, sep = \n), чтобы он читался в более знакомом формате матрицы (каждый список в новой строке). Конечно, это менее идеальный способ сделать это, поэтому я также пошел дальше и создал его с помощью NumPy. - person LampoonBonanza; 02.03.2021