This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Need help with my Needleman-Wunsch Python program (traceback) (error)

Good Afternoon Biostar Community,

I have a problem with my python program. I tried to implement the Needleman-Wunsch algorithm but I get the wrong alignment due to a mistake in the traceback section I assume. I got the correct matrix but not the correct sequence alignment. With this NW("SCHULFACH","SCHLUESSEL","2","-1","-1") input values I should get this solution "SCH_ULFACH || SCHLUESSEL" but I get this "----SCHULFACH and SCHL---UESSEL". It would help me a lot if you have a look at my program and try to figure out what I've done wrong.

Thanks in advance!

def NW(Seq1,Seq2,match,mismatch,gap):

#Initialization of the matrix
matrix=[[x for x in range(len(Seq1)+1)] for i in range(len(Seq2)+1)]

#set matrix to zero
for x in range(len(Seq1)+1):
    for i in range(len(Seq2)+1):
        matrix[i][x]=0



#Calculation of the matrix
traceback=[[x for x in range(len(Seq1)+1)] for i in range(len(Seq2)+1)]
a=0
b=0
if len(Seq1)==len(Seq2):
    a=len(Seq1)+1
    b=len(Seq2)+1
else:
    if len(Seq1)>len(Seq2):
        a=len(Seq2)+1
        b=len(Seq1)+1

    else:
        b=len(Seq1)+1
        a=len(Seq2)+1





for i in range(a):
    for x in range(b):

     if x==0 and i==0:
        matrix[0][0]=0

     else:

        if x==0:
            f1=float("-inf")
            f2=matrix[i-1][x]
            f3=float("-inf")

        elif i==0:
            f1=float("-inf")
            f2=float("-inf")
            f3=matrix[i][x-1]

        else:
            f1=matrix[i-1][x-1]
            f2=matrix[i][x-1]
            f3=matrix[i-1][x]



        if(x>0 and i>0):
            if(Seq2[i-1]==Seq1[x-1]):
                m=int(match)
            else:
                m=int(mismatch)
        else:
            m=int(mismatch)



        w1=f1+m
        w2=f2+int(gap)
        w3=f3+int(gap)
        erg=(w1,w2,w3)
        w4=max(erg)


        if w4==w1:
            traceback[i][x]=0

        elif w4==w2:
            traceback[i][x]=1 


        else:#w4==w3
            traceback[i][x]=-1


        try:
            matrix[i][x]=w4

        except IndexError:
            print("ERROR")

#traceback

alnseqA=""#horizontal
alnseqB=""#vertical
j=len(Seq1)#horizontal
i=len(Seq2)#vertical
print(j,i)

while(i>0 or j>0):
    if traceback[i][j]==0:
        alnseqA += Seq1[j-1]
        alnseqB += Seq2[i-1]
        j=j-1
        i=i-1
    elif traceback[i][j]==1:
        alnseqA += "-"
        alnseqB += Seq2[i-1]
        i=i-1
    else:
        alnseqB += "-"
        alnseqA += Seq1[j-1]
        j=j-1
print(matrix)
print(alnseqA[::-1])
print(alnseqB[::-1])


return (alnseqA,alnseqB)

NW("GCTCAA","CTAAG","2","-1","-2")

NW("SCHULFACH","SCHLUESSEL","2","-1","-1")

alignment rna-seq software error

It seems you're using the wrong scores:

if x==0:
            f1=float("-inf")
            **f2=matrix[i-1][x]** 
            f3=float("-inf")

        elif i==0:
            f1=float("-inf")
            f2=float("-inf")
            **f3=matrix[i][x-1]**

        else:
            f1=matrix[i-1][x-1]
            **f2=matrix[i][x-1]** # should be f2=matrix[i-1][x]
            **f3=matrix[i-1][x]** # should be f3=matrix[i][x-1]

edit: fixed mis-annotation

I would like to know, why your code return an error,

return (alnseqA,alnseqB) ^ SyntaxError: 'return' outside function

Please, can you explain it,

There is some function define somewhere?

NW("GCTCAA","CTAAG","2","-1","-2") NW("SCHULFACH","SCHLUESSEL","2","-1","-1")

Thanks for your help

0 answers

No answers yet.

Log in to answer this question.