1. import numpy as np
2.  
3. #seqProtein1
4. sekwencja1 = []
5. plik1 = open('seqProtein1.fasta', "r").readlines()
6. for linia in plik1[1:]:
7.     sekwencja1 += linia.rstrip()
8.     print(sekwencja1)
9.  
10. print()
11.  
12. sekw1 = sekwencja1
13.  
14. #seqProtein2
15. sekwencja2 = []
16. plik2 = open('seqProtein2.fasta', "r").readlines()
17. for linia in plik2[1:]:
18.     sekwencja2 += linia.rstrip()
19.     print(sekwencja2)
20.  
21. sekw2 = sekwencja2
22.  
23. #Macierz BLOSUM
24. S = []
25. plik3 = open('blosum62.txt', "r")
26. S = eval(plik3.read())
27.  
28. print(S)
29.  
30. # Sekwencje dlugosci
31. sekw1_len = len(sekw1)
32. sekw2_len = len(sekw2)
33.  
34. # kara
35. kara = 4
36.  
37. # Macierz scoringowa
38. F = np.zeros((sekw2_len + 1, sekw1_len + 1), dtype=int)
39.  
40. #wypelnianie macierzy
41. for i in range(1, F.shape[0]):
42.     for j in range(1, F.shape[1]):
43.  
44.         literka_gura = sekw1[j - 1]
45.         literka_bok = sekw2[i - 1]
46.  
47.         pszekontna = F[i - 1][j - 1] + S[literka_gura][literka_bok]
48.         gura = F[i - 1][j] - kara
49.         dul = F[i][j - 1] - kara
50.  
51.         F[i][j] = max(pszekontna, gura, dul, 0)
52.  
53. print(F)
54.  
55. maksymalny = F[0][0]
56. indeks_i = list()
57. indeks_j = list()
58.  
59.  
60. #indeksy maksymalnego elementu
61. for i in range(1, F.shape[0]):
62.     for j in range(1, F.shape[1]):
63.         if F[i][j] >= maksymalny:
64.             maksymalny = F[i][j]
65.  
66.  
67. #znalezienie indeksu maksymalnych
68. for i in range(1, F.shape[0]):
69.     for j in range(1, F.shape[1]):
70.         if F[i][j] == maksymalny:
71.             indeks_i.append(i)
72.             indeks_j.append(j)
73.  
74.  
75. #sekwencje
76. sekw1_wynikowa = list()
77. sekw2_wynikowa = list()
78.  
79. #traceback
80. i = indeks_i[0]
81. j = indeks_j[0]
82.  
83. while i != 0 or j != 0:
84.         literka_gura = sekw1[j - 1]
85.         literka_bok = sekw2[i - 1]
86.         literka = literka_gura + literka_bok
87.  
88.         pszekontna = F[i - 1][j - 1]
89.         gura = F[i - 1][j] - kara
90.         dul = F[i][j - 1] - kara
91.  
92.         F[i][j] = max(pszekontna, gura, dul)
93.  
94.         if F[i][j] == pszekontna:
95.             i = i - 1
96.             j = j - 1
97.             sekw1_wynikowa.append(literka_gura)
98.             sekw2_wynikowa.append(literka_bok)
99.         elif F[i][j] == gura:
100.             i = i - 1
101.             sekw1_wynikowa.append("-")
102.             sekw2_wynikowa.append(literka_bok)
103.         else:
104.             j = j - 1
105.             sekw1_wynikowa.append(literka_gura)
106.             sekw2_wynikowa.append("-")
107.  
108. print()
109.  
110. # print(sekw1_2_wynikowa)
111. tt = 0
112. sekw1_2_wynikowa = list()
113. for i in reversed(sekw1_wynikowa):
114.     sekw1_2_wynikowa.append(i)
115.     tt = tt + 1
116.  
117. print("".join(sekw1_2_wynikowa))
118.  
119. # print(sekw2_2_wynikowa)
120. tt = 0
121. sekw2_2_wynikowa = list()
122. for i in reversed(sekw2_wynikowa):
123.     sekw2_2_wynikowa.append(i)
124.     tt = tt + 1
125.  
126. print("".join(sekw2_2_wynikowa))