Demo entry 6789519

python3

Submitted by anonymous on Apr 30, 2019 at 04:17
Language: Python 3. Code size: 1.9 kB.

# -*- coding: utf-8 -*-
# @Author: orres
# @Email:  jasonjin22@gmail.com
# @Date:   2019-04-29 11:08:22
import copy
import math

# used to calcualte the score
def score(a, b):
if (a == b):
return (1)
else:
return (-1)

P = 'GGT'
T = 'ACCGGTGCTGTAA'
penalty = -1

global_max = 0
best_k = 1
best_seq = ''

for k in range(1, math.floor(len(T) / len(P) + 1)):
# iterate k from 1 to maximum
string = P
for i in range(k-1):
string += P
P = string
# init the dp matrix
dp = []
trace = []
for i in range(k*len(P) + 1):
line = []
for j in range(len(T) + 1):
line.append(0)
dp.append(line)
trace.append(copy.deepcopy(line))
for i in range(k*len(P) + 1):
dp[i] = -i

# calculate the remaining elements in matrix
for i in range(1, len(P) + 1):
for j in range(1, len(T) + 1):
dp[i][j] = max(dp[i-1][j-1] + score(P[i-1], T[j-1]), dp[i-1][j] + penalty, dp[i][j-1]+penalty)
# the trace matrix used to record when the element comes from
if (dp[i][j] == dp[i-1][j] + penalty):
trace[i][j] = 2
elif (dp[i][j] == dp[i-1][j-1] + score(P[i-1], T[j-1])):
trace[i][j] = 1
else:
trace[i][j] = 3
# find the highest score in the matrix
big = 0
pos = (0, 0)
for i in range(k*len(P) + 1):
for j in range(len(T) + 1):
if dp[i][j] >= big:
big = dp[i][j]
pos = (i, j)

# get the aligned sequence
x = pos
y = pos
seq = ""
for i in range(len(P)):
if (x < 0 or y < 0):
break
if (trace[x][y] == 1):
seq = T[y-1] + seq
x -= 1
y -= 1
elif (trace[x][y] == 2):
x -= 1
seq = '-' + seq
else:
y -= 1

# choose the best k
if (big > global_max):
global_max = big
best_k = k
best_seq = seq

# print the best k and the best aligned sequence
print(best_k)
print(best_seq)

This snippet took 0.01 seconds to highlight.

Back to the Entry List or Home.