# BioInfa 9.1

Mar 1st, 2021 (edited)
446
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. import numpy
2.
3.
4. def match_mismatch(first_dna, second_dna, i, j):
5.     first_dna = " " + first_dna
6.     second_dna = " " + second_dna
7.     return 1 if first_dna[j] == second_dna[i] else 0
8.
9. def print_norm_first(first):
10.     for i in range(1, len(first)):
11.         if first[i] == "-":
12.             first[i] = first[i-1]
13.             first[i-1] = "-"
14.     print(first)
15.
16. def find_alignment(p, first_dna, second_dna, n, m, new_first, new_second, d, match, mismatch):
17.     i = n
18.     j = m
19.     while i > 0 or j > 0:
20.         score = p[i, j]
21.         score_up = p[i - 1, j]
22.         score_left = p[i, j - 1]
23.         score_diag = p[i - 1, j - 1]
24.         if score == score_up + d:
25.             new_first.append("-")
26.             new_second.append(second_dna[i - 1])
27.             i -= 1
28.         elif score == score_left + d:
29.             new_first.append(first_dna[j - 1])
30.             new_second.append("-")
31.             j -= 1
32.
33.         elif score == score_diag + match or score == score_diag + mismatch:
34.             new_first.append(first_dna[j - 1])
35.             new_second.append(second_dna[i - 1])
36.             i -= 1
37.             j -= 1
38.
39.     if i > 0:
40.         new_second.append("_")
41.         i -= 1
42.     if j > 0:
43.         new_first.append("_")
44.         j -= 1
45.     return new_first, new_second
46.
47.
48. def find_sum(match, mismatch, gap, first, second, local):
49.     n = len(first)
50.     l = 0  # count of matches
51.     f = 0  # count of mismatches
52.     u = 0  # count of indels
53.     for k in range(0, n - 1):
54.         if first[k] == second[k]:
55.             l += 1
56.         elif first[k] == "-" or second[k] == "-":
57.             u += 1
58.         else:
59.             f += 1
60.     if local:
61.         print("Final score=", match * (l + 1) + mismatch * f + gap * u, " match=", l + 1, "mismatch=", f, "gap=", u)
62.     else:
63.         print("Final score=", match * l + mismatch * f + gap * u, " match=", l, "mismatch=", f, "gap=", u)
64.
65.
66. def global_alignment(n, m, match, mismatch, gap, first_dna, second_dna):
67.     local = 0
68.     p = numpy.zeros((n + 1, m + 1), dtype=numpy.int64)
69.     new_first = []
70.     new_second = []
71.     for j in range(1, m + 1):
72.         p[0, j] = p[0, j - 1] + gap
73.     for i in range(1, n + 1):
74.         p[i, 0] = p[i - 1, 0] + gap
75.     for i in range(1, n + 1):
76.         for j in range(1, m + 1):
77.             g = match_mismatch(first_dna, second_dna, i, j)
78.             if g:
79.                 first_func = p[i - 1, j - 1] + match
80.             else:
81.                 first_func = p[i - 1, j - 1] + mismatch
82.             second_func = p[i - 1, j] + gap
83.             third_func = p[i, j - 1] + gap
84.             p[i, j] = max(first_func, second_func, third_func)
85.     find_alignment(p, first_dna, second_dna, n, m, new_first, new_second, gap_fine, match, mismatch)
86.     find_sum(match, mismatch, gap, new_first, new_second, local)
87.     print_norm_first(new_first[::-1])
88.     print(new_second[::-1])
89.
90.
91. def max_in_matrix(matrix, n, m):
92.     for i in range(1, n + 1):
93.         for j in range(1, m + 1):
94.             if max(map(max, matrix)) == matrix[i, j]:
95.                 return i, j
96.
97.
98. def find_alignment_local(p, first_dna, second_dna, n, m, new_first, new_second, d, match, mismatch):
99.     i, j = max_in_matrix(p, n, m)
100.     score = p[i, j]
101.     while score != 0:
102.         score_up = p[i - 1, j]
103.         score_left = p[i, j - 1]
104.         score_diag = p[i - 1, j - 1]
105.         if score == score_up + d and i > 0:
106.             new_first.append("-")
107.             new_second.append(second_dna[i - 1])
108.             i = i - 1
109.             score = p[i, j]
110.         elif score == score_left + d:
111.             new_first.append(first_dna[j - 1])
112.             new_second.append("-")
113.             j = j - 1
114.             score = p[i, j]
115.         elif (score == score_diag + match or score == score_diag + mismatch) and i > 0 and j > 0:
116.             new_first.append(first_dna[j - 1])
117.             new_second.append(second_dna[i - 1])
118.             i = i - 1
119.             j = j - 1
120.             score = p[i, j]
121.
122.     return new_first, new_second
123.
124.
125. def local_alignment(n, m, match, mismatch, gap, first_dna, second_dna):
126.     p = numpy.zeros((n + 1, m + 1), dtype=numpy.int64)
127.     local = 1
128.     new_first = []
129.     new_second = []
130.     for i in range(1, n + 1):
131.         for j in range(1, m + 1):
132.             if first_dna[j - 1] == second_dna[i - 1]:
133.                 first_func = p[i - 1, j - 1] + match
134.             else:
135.                 first_func = p[i - 1, j - 1] + mismatch
136.             second_func = p[i - 1, j] + gap
137.             third_func = p[i, j - 1] + gap
138.             p[i, j] = max(0, first_func, second_func, third_func)
139.     find_alignment_local(p, first_dna, second_dna, n, m, new_first, new_second, gap, match, mismatch)
140.     find_sum(match, mismatch, gap, new_first, new_second, local)
141.     print(new_first[::-1])
142.     print(new_second[::-1])
143.
144.
145. if __name__ == '__main__':
146.     first_dna = input("Enter first DNA\n")
147.     second_dna = input("Enter second DNA\n")
148.     match_bonus = int(input("Enter match bonus\n"))
149.     mismatch_fine = int(input("Enter mismatch fine\n"))
150.     gap_fine = int(input("Enter gap fine\n"))
151.     alignmentType = input("Enter type of alignment(G or L)\n")
152.     m = len(first_dna)
153.     n = len(second_dna)
154.     if alignmentType == "G":
155.         global_alignment(n, m, match_bonus, mismatch_fine, gap_fine, first_dna, second_dna)
156.     else:
157.         local_alignment(n, m, match_bonus, mismatch_fine, gap_fine, first_dna, second_dna)
158.
RAW Paste Data