7

すでにアラインメントされたシーケンスをスコアリングしようとします。言いましょう

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

与えられたパラメータで

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

私はbiopythonクックブックを調べましたが、取得できるのは置換マトリックスblogsum62だけですが、誰かがこの種のライブラリをすでに実装している必要があると感じています。

それで、誰かが私の問題を解決できるライブラリまたは最短のコードを提案できますか?

事前にThx

4

2 に答える 2

11

ジェサダ、

Blosum62マトリックス(スペルに注意してください;)はBio.SubsMat.MatrixInfoにあり、タプルがスコアに解決される辞書です(つまり、('A', 'A')4ポイントの価値があります)。ギャップはなく、マトリックスの三角形は1つだけです(したがって、('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を返します。これらすべてを行うには、ほぼ確実にきれいな方法がありますが、それは良いスタートになるはずです。

于 2011-04-18T05:38:56.207 に答える
9

blosum62は276アイテムの辞書です。

分析対象のシーケンスには276を超える要素が含まれている可能性が高いのに対し、これは276ターンの反復のみを表すため、不足している項目で完了することを好みました。したがって、関数score_match()を使用して各ペアのスコアを見つけると、この関数は if pair not in matrixシーケンスの各要素に対してテストを実行する必要があります。つまり、確かに276回をはるかに超えます。

時間がかかるもう1つのことは、それぞれが新しい整数を作成し、名前スコアをこの新しいオブジェクトに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()は、returnの代わりにyieldがあるため、ジェネレーター関数です。

編集:Python 3のコードを更新しました:

from Bio.SubsMat.MatrixInfo import blosum62 as blosum

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

def score_pairwise(seq1, seq2, matrix, gap_s, gap_e, gap = True):
    for A,B in zip(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)))
于 2011-04-18T20:49:38.640 に答える