Существует ли функция, которая может вычислить оценку для выровненных последовательностей с учетом параметров выравнивания? - PullRequest
6 голосов
/ 16 апреля 2011

Я пытаюсь оценить уже выровненные последовательности. Пусть скажут

seq1 = 'PAVKDLGAEG-ASDKGT--SHVVY----------TI-QLASTFE'
seq2 = 'PAVEDLGATG-ANDKGT--LYNIYARNTEGHPRSTV-QLGSTFE'

с заданными параметрами

substitution matrix : blosum62
gap open penalty : -5
gap extension penalty : -1

Я просмотрел кулинарную книгу по биопифону, но все, что я могу получить, это подстановочная матрица blogsum62, но я чувствую, что кто-то уже реализовал такую ​​библиотеку.

Так кто-нибудь может предложить какие-либо библиотеки или кратчайший код, который может решить мою проблему?

Спасибо заранее

Ответы [ 2 ]

9 голосов
/ 18 апреля 2011

Jessada,

Матрица Blosum62 (обратите внимание на орфографию;) находится в Bio.SubsMat.MatrixInfo и представляет собой словарь с кортежами, преобразующимися в оценки (поэтому ('A', 'A') стоит 4 балла).Это не имеет пробелов, и это только один треугольник матрицы (так что это может быть ('T', 'A'), но не ('A', 'T'). В Biopython есть некоторые вспомогательные функции,включая некоторые из них в Bio.Pairwise, но вот что я придумал в качестве ответа:

from Bio.SubsMat import MatrixInfo

def score_match(pair, matrix):
    if pair not in matrix:
        return matrix[(tuple(reversed(pair)))]
    else:
        return matrix[pair]

def score_pairwise(seq1, seq2, matrix, gap_s, gap_e):
    score = 0
    gap = False
    for i in range(len(seq1)):
        pair = (seq1[i], seq2[i])
        if not gap:
            if '-' in pair:
                gap = True
                score += gap_s
            else:
                score += score_match(pair, matrix)
        else:
            if '-' not in pair:
                gap = False
                score += score_match(pair, matrix)
            else:
                score += gap_e
    return score

seq1 = 'PAVKDLGAEG-ASDKGT--SHVVY----------TI-QLASTFE'
seq2 = 'PAVEDLGATG-ANDKGT--LYNIYARNTEGHPRSTV-QLGSTFE'

blosum = MatrixInfo.blosum62

score_pairwise(seq1, seq2, blosum, -5, -1)

, который возвращает 82 для вашего выравнивания. Есть почти наверняка более красивые способы сделать все это, но это должно бытьхорошее начало.

6 голосов
/ 19 апреля 2011

blosum62 - это диктонар из 276 наименований.

Я предпочел дополнить недостающие элементы, потому что он представляет собой итерацию всего 276 ходов, в то время как анализируемые последовательности, вероятно, будут иметь более 276 элементов. Следовательно, если вы найдете оценку каждой пары с помощью функции Score_match (), эта функция должна будет выполнить тест if pair not in matrix для каждого из элементов последовательностей, то есть, безусловно, гораздо больше, чем 276 раз .

Еще одна вещь, которая занимает много времени: каждый score += something создает новое целое число и связывает имя оценка с этим новым объектом. Каждое связывание занимает количество времени, которое не существует с потоком целых чисел от генератора, которые немедленно добавляются к текущему количеству.

from Bio.SubsMat.MatrixInfo import blosum62 as blosum
from itertools import izip

blosum.update(((b,a),val) for (a,b),val in blosum.items())

def score_pairwise(seq1, seq2, matrix, gap_s, gap_e, gap = True):
    for A,B in izip(seq1, seq2):
        diag = ('-'==A) or ('-'==B)
        yield (gap_e if gap else gap_s) if diag else matrix[(A,B)]
        gap = diag

seq1 = 'PAVKDLGAEG-ASDKGT--SHVVY----------TI-QLASTFE'
seq2 = 'PAVEDLGATG-ANDKGT--LYNIYARNTEGHPRSTV-QLGSTFE'

print sum(score_pairwise(seq1, seq2, blosum, -5, -1))

Этот Score_pairwise () является функцией генератора, потому что есть yield вместо return .

...