Traceback in der dynamischen Programmierung die Umsetzung des Needleman-Wunsch-Algorithmus
Ich habe fast mein needleman-wunsch-Implementierung arbeiten, aber ich bin verwirrt über wie man mit den traceback auf einen bestimmten Fall.
Die Idee ist, dass, um zu re konstruieren die Folge (längster Pfad), die wir neu zu berechnen, um zu bestimmen, die matrix der Partitur kam. Die edge-Fall, ich habe ein problem mit ist, wenn der rechten unteren Kerbe ist nicht in der match-matrix, sondern ist in die Spalte einfügen-matrix (was bedeutet, dass die resultierende zurückgeführt Folge haben sollte, einfügen.
Diese Sequenzen sind aufgezeichnet in der a2m-format, wobei Einsätze in der Sequenz aufgezeichnet werden als Kleinbuchstaben. Also in der letzten Ausgabe der Ausrichtung der ZZ zu AAAC sollte AAac. Wenn ich mit trace back von hand habe ich am Ende mit AAAc weil ich nur besuchen Sie die Ic-matrix einmal. Hier ist ein Bild von meinem whiteboard. Wie Sie sehen können, habe ich drei schwarze Pfeile und einen grünen Pfeil, der ist, warum meine traceback gibt mir AAAc. Sollte ich das zählen der ersten Zelle, dann halt an position 1,1? Ich bin nicht sicher, wie ich möchte ändern die Weise, die ich umgesetzt habe, dies zu tun.
Beachten Sie, dass die substitution matrix verwendet wird, hier ist BLOSUM62. Die Wiederholung Beziehungen sind
M(i,j) = max(Ic(i-1,j-1)+subst, M(i-1,j-1)+subst, Ir(i-1,j-1)+subst)
Ic(i,j) = max(Ic(i,j-1)-extend, M(i,j-1)-open, Ir(i,j-1)-double)
Ir(i,j) = max(Ic(i-1,j)-double, M(i-1,j)-open, Ir(i-1,j)-extend)EDIT: Hier ist der traceback_col_seq Funktion neu geschrieben werden sauberer. Beachten Sie, dass score_cell jetzt gibt thisM,thisC,thisR statt der max dieser. Diese version erzielt die Ausrichtung als AaAc, immer noch das gleiche problem, und nun mit einem anderen problem, warum geht es dann in die Ic wieder bei 1,2. Diese version ist viel besser lesbar ist jedoch.
def traceback_col_seq(self):
i, j = self.maxI-1, self.maxJ-1
self.traceback = list()
matrixDict = {0:'M',1:'Ic',2:'Ir',3:'M',4:'Ic',5:'Ir',6:'M',7:'Ic',8:'Ir'}
while i > 0 or j > 0:
chars = self.col_seq[j-1] + self.row_seq[i-1]
thisM, thisC, thisR = self.score_cell(i, j, chars)
cell = thisM + thisC + thisR
prevMatrix = matrixDict[cell.index(max(cell))]
print(cell, prevMatrix,i,j)
if prevMatrix == 'M':
i -= 1; j -= 1
self.traceback.append(self.col_seq[j])
elif prevMatrix == 'Ic':
j -= 1
self.traceback.append(self.col_seq[j].lower())
elif prevMatrix == 'Ir':
i -= 1
self.traceback.append('-')
return ''.join(self.traceback[::-1])Hier ist die python-Klasse, die generiert die dynamische Programmierung matrix und Spuren zurück, die Ausrichtung. Es gibt auch einen score-Funktion verwendet, um zu überprüfen, wenn die Ausrichtung hergestellt ist richtig.
class global_aligner():
def __init__(self, subst, open=12, extend=1, double=3):
self.extend, self.open, self.double, self.subst = extend, open, double, subst
def __call__(self, row_seq, col_seq):
#add alphabet error checking?
score_align(row_seq, col_seq)
return traceback_col_seq()
def init_array(self):
"""initialize three numpy arrays, set values in 1st column and row"""
self.M = zeros((self.maxI, self.maxJ), float)
self.Ic = zeros((self.maxI, self.maxJ), float)
self.Ir = zeros((self.maxI, self.maxJ), float)
for i in xrange(1,self.maxI):
self.M[i][0], self.Ic[i][0], self.Ir[i][0] = \
-float('inf'), -float('inf'), -(self.open+self.extend*(i-1))
for j in xrange(1,self.maxJ):
self.M[0][j], self.Ir[0][j], self.Ic[0][j] = \
-float('inf'), -float('inf'), -(self.open+self.extend*(j-1))
self.Ic[0][0] = self.Ir[0][0] = -float('inf')
def score_cell(self, i, j, chars):
"""score a matrix cell based on the 9 total neighbors (3 each direction)"""
thisM = [self.M[i-1][j-1]+self.subst[chars], self.Ic[i-1][j-1]+ \
self.subst[chars], self.Ir[i-1][j-1]+self.subst[chars]]
thisC = [self.M[i][j-1]-self.open, self.Ic[i][j-1]-self.extend, \
self.Ir[i][j-1]-self.double]
thisR = [self.M[i-1][j]-self.open, self.Ic[i-1][j]-self.double, \
self.Ir[i-1][j]-self.extend]
return max(thisM), max(thisC), max(thisR)
def score_align(self, row_seq, col_seq):
"""build dynamic programming matrices to align two sequences"""
self.row_seq, self.col_seq = list(row_seq), list(col_seq)
self.maxI, self.maxJ = len(self.row_seq)+1, len(self.col_seq)+1
self.init_array() #initialize arrays
for i in xrange(1, self.maxI):
row_char = self.row_seq[i-1]
for j in xrange(1, self.maxJ):
chars = row_char+self.col_seq[j-1]
self.M[i][j], self.Ic[i][j], self.Ir[i][j] = self.score_cell(i, j, chars)
def traceback_col_seq(self):
"""trace back column sequence in matrices in a2m format"""
i, j = self.maxI-1, self.maxJ-1
self.traceback = list()
#find which matrix to start in
#THIS IS WHERE THE PROBLEM LIES I THINK
cell = (self.M[i][j], self.Ic[i][j], self.Ir[i][j])
prevMatrix = cell.index(max(cell))
while i > 1 and j > 1:
if prevMatrix == 0: #M matrix
i -= 1; j -= 1 #step up diagonally
prevChars = self.row_seq[i-1]+self.col_seq[j-1]
diag = self.score_cell(i, j, prevChars) #re-score diagonal cell
prevMatrix = diag.index(max(diag)) #determine which matrix that was
self.traceback.append(self.col_seq[j])
elif prevMatrix == 1: #Ic matrix
j -= 1
prevChars = self.row_seq[i-1]+self.col_seq[j-1]
left = self.score_cell(i, j, prevChars)
prevMatrix = left.index(max(left))
self.traceback.append(self.col_seq[j].lower())
elif prevMatrix == 2: #Ir matrix
i -= 1
prevChars = self.row_seq[i-1]+self.col_seq[j-1]
up = self.score_cell(i, j, prevChars)
prevMatrix = up.index(max(up))
self.traceback.append('-')
for j in xrange(j,0,-1): #hit top of matrix before ending, add chars
self.traceback.append(self.col_seq[j-1])
for i in xrange(i,0,-1): #hit left of matrix before ending, add gaps
self.traceback.append('-')
print(''.join(self.row[::-1]))
return ''.join(self.traceback[::-1])
def score_a2m(self, s1, s2):
"""scores an a2m alignment of two sequences. I believe this function correctly
scores alignments and is used to test my alignments. The value produced by this
function should be the same as the largest value in the bottom right of the three
matrices"""
s1, s2 = list(s1.strip('.')), list(s2.strip('.'))
s1_pos, s2_pos = len(s1)-1, len(s2)-1
score, gap = 0, False
while s1_pos >= 0 and s2_pos >= 0:
if s1[s1_pos].islower() and gap is False:
score -= self.open; s1_pos -= 1; gap = True
elif s1[s1_pos].islower() and gap is True:
score -= self.extend; s1_pos -= 1
elif s2[s2_pos].islower() and gap is False:
score -= self.open; s2_pos -= 1; gap = True
elif s2[s2_pos].islower() and gap is True:
score -= self.extend; s2_pos -= 1
elif s1[s1_pos] == '-' and gap is False:
score -= self.open; s1_pos -= 1; s2_pos -= 1; gap = True
elif s1[s1_pos] == '-' and gap is True:
score -= self.extend; s1_pos -= 1; s2_pos -= 1
elif s2[s2_pos] == '-' and gap is False:
score -= self.open; s1_pos -= 1; s2_pos -= 1; gap = True
elif s2[s2_pos] == '-' and gap is True:
score -= self.extend; s1_pos -= 1; s2_pos -= 1
elif gap is True:
score += self.subst[s1[s1_pos].upper() + s2[s2_pos].upper()]
s1_pos -= 1; s2_pos -= 1; gap = False
else:
score += self.subst[s1[s1_pos].upper() + s2[s2_pos].upper()]
s1_pos -= 1; s2_pos -= 1
if s1_pos >= 0 and gap is True:
score -= self.extend*s1_pos
elif s1_pos >= 0 and gap is False:
score -= self.open+s1_pos*self.extend
if s2_pos >= 0 and gap is True:
score -= self.extend*s2_pos
elif s2_pos >= 0 and gap is False:
score -= self.open+s2_pos*self.extend
return score
test = global_aligner(blosumMatrix)
s1,s2 = 'ZZ','AAAC'
test.score_align(s1, s2)
align = test.traceback_col_seq()
print('This score: ', test.score_a2m(s1,align)
print('Correct score: ', test.score_a2m(s1,'AAac'))Blosum-parser
def parse_blosum(blosumFile):
blosumMatrix, commentFlag = dict(), False
for line in blosumFile:
if not line.startswith('#') and not commentFlag:
alphabet = line.rstrip().split()
commentFlag = True
elif commentFlag:
line = line.rstrip().split()
thisChar, line = line[0], line[1:]
for x in xrange(len(line)):
alphaChar, thisValue = alphabet[x], line[x]
blosumMatrix[thisChar+alphaChar] = int(thisValue)
return blosumMatrixDu musst angemeldet sein, um einen Kommentar abzugeben.
Als Sie spielte auf Ihre Bemerkung über den Pfeil Farbe, das grundlegende problem ist, dass die zweite horizontale Pfeil ist immer gekennzeichnet als eine Bewegung, die in M (schwarzer Pfeil), wenn es wirklich eine Bewegung im Ic (grüner Pfeil). Dies scheint passiert zu sein, weil die
prevMatrixvariable gibt die beste matrix (i-1, j-1), (i-1, j) oder (i, j-1). Das ist die matrix, in die Vorherige Zelle von den vorausschauenden Perspektive. Da du tust trace*zurück* an diesem Punkt die Zelle, die Sie "aus" in der Spur ist wirklich (i, j) - die, die Sie sind momentan in. Bestimmt die Richtung, in die Sie sich bewegen und so, ob Sie sind, idem Sie eine Lücke oder nicht. Ihre beste option ist wahrscheinlich eine variable, verfolgt die aktuelle matrix iteration von Schleife zu Schleife, iteration, und verwenden Sie diese, um zu bestimmen, welche Zeichen Anhängen, um Ihre Ausrichtung zu string. Bei jeder iteration können Sie aktualisieren Sie es nach der Festlegung, welche Zelle du gehst zum nächsten.Ich denke, dass diese änderungen klären das problem, dass Sie mit Ihrem neu geschrieben traceback_col_sequence Funktion, aber vorläufig sieht es aus wie Sie nicht, unter Berücksichtigung der Informationen, die Sie gewonnen haben durch die Verfolgung zurück. Du bist Wiederverwendung
score_cell(), die sieht aus wie es entworfen wurde, zu berechnen, erzielt der nach vorne geht. Wenn Sie die Ablaufverfolgung zurück, Sie wissen bereits, das Ende, so dass Sie nicht wollen, um die maximale Punktzahl für die Vorherige Zelle. Sie wollen wählen Sie die maximale Punktzahl, die führen konnte, um die Lage in der matrix, die Sie derzeit bei. Zum Beispiel, wenn Sie wissen, dass die matrix, die Sie mit Ihren traceback ist M, dann wissen Sie, Sie nicht bekommen, da durch die Erweiterung eine Lücke, so sollte man nicht überlegen, diese Optionen.