Byan 13518066

README.md

@@ -31,30 +31,8 @@ Terdapat dua titik capai dalam tugas ini :
 [3] Anda dapat membuat 2 buah algoritma untuk perbanding antar sekuens dan perbandingan antar profil. Coba anda ubah agar anda melakukan keduanya secara langsung, sehingga mengurangi beban kerja anda.<br>
 [4] Anda bisa saja menggunakan Python untuk mengerjakan tugas, namun ingat bahwa **kinerja C++ dan C jauh lebih cepat**. Sebagai pengalaman, asisten menjalankan skoring pensejajaran global 2 DNA dengan panjang ~29000 nukleotida. Algoritma berjalan 100 menit untuk Python, dan algoritma berjalan hanya 3 menit untuk bahasa C++ dengan flag -O3 (optimization) ketika kompilasi (**33 x speedup !**). Sebagai saran (bila anda keukeuh menggunakan Python) , anda bisa menggunakan Python optimizer (misal Numba ataupun Cython) untuk mempercepat eksekusi algoritma anda.
 
-## Pengumpulan
-### Pengerjaan
-Silahkan lakukan *fork* dari *repository* ini.
-
-### Deliverables
-1. File yang berisi hasil pensejajaran sekuens global. Misal , hasil pensejajaran antara sekuens pertama (file1.fasta) dan sekuens kedua (file2.fasta) ditulis dalam folder /result/file1_file2/. Lalu dalam folder file1_file2, tuliksan hasil pensejajaran masing-masing file1.fasta dan file2.fasta sebagai file1.txt dan file2.txt. Tuliskan score dalam sebuah file score.txt. 
-2. Kode sumber. Tuliskan cara kompilasi bila menggunakan *compiled language*, namun lebih baik dilengkapi dengan Makefile.
-3. Ubah Readme ini. Tuliskan pendekatan pensejajaran yang anda lakukan dan cara menjalankan program.
-
-### Teknis Pengumpulan
-&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Kumpulkan dengan membuat *merge request* pada *repository* ini. Batas pengumpulan dan demo adalah 29 Juli 2020.
-
-### Demo
-&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Setelah selesai, jadwalkan demo dengan asisten. Kontak dapat dilihat pada Readme ini. Demo berlangsung 15-30 menit. Demo akan berisi tanya jawab, namun belum tentu akan diisi oleh pengujian, tergantung *runtime* dari algoritma anda.
-
-## Penilaian
-&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Uji cobakan algoritma anda dengan data yang telah disediakan di *repository* ini. Percobaan minimal menghasilkan 3 pensejajaran global untuk kasus 2 sekuens, 2 pensejajaran global kasus 3 sekuens, 1 pensejajaran global untuk kasus 4 sekuens.<br>
-&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Nilai maksimal 1350 untuk *milestone* pertama, dan nilai maksimal 1800 untuk milestone kedua. Nilai maksimal demo adalah 850, sehingga nilai maksimal total adalah 4000. Algoritma anda **wajib optimum** untuk pensejajaran global 2 sekuens. Akan tetapi, algoritma anda tidak harus optimum untuk MSA. Implementasi MSA lebih kepada *proof-of-concept*.
-
-## Kontak
-&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;Silahkan hubungi asisten lewat line @alamhasabiebaru atau lewat email 13517096@std.stei.itb.ac.id dengan subjek  \[SELEKSI IRK - SEQUENCE ALIGNMENT\] . *Note : waktu menjawab bervariasi, namun email biasanya akan dibalas kurang dari sehari. Line mungkin tidak dibalas dalam waktu satu-dua hari. Mohon bersabar :)*. Pertanyaan juga dipersilahkan. Jawaban akan diposting dalam bagian QnA README ini.
-
-## QnA
-null.
+## Run
+python main.py
 
 ## Referensi
 Silahkan gunakan referensi berikut sebagai awal pengerjaan tugas:<br>

main.py

@@ -0,0 +1,110 @@
+def read_fasta(fp):
+        name, seq = None, []
+        for line in fp:
+            line = line.rstrip()
+            if line.startswith(">"):
+                if name: yield (name, ''.join(seq))
+                name, seq = line, []
+            else:
+                seq.append(line)
+        if name: yield (name, ''.join(seq))
+
+def given_matrices_inserter(filename):
+	try:	
+		data = open(filename, 'r')
+		dimensions = data.readline()
+		try:
+			n = int(dimensions)
+		except ValueError:
+			print("Wrong file")
+			exit()
+		letters = data.readline()
+		letters = letters.replace('\n', '')
+		letters_arr = letters.split(' ')
+		score_matrix = np.zeros((n,n))
+		for i in range(0,n):
+			arr = data.readline().split(" ")
+			for j in range(0,n):
+				score_matrix[i][j] = float(arr[j])
+		return letters_arr, score_matrix
+	except FileNotFoundError:
+		print("There is no matrix file.")
+		exit()		
+
+import numpy as np
+
+def write_to_file(filename, value):
+    with open(filename, "w+") as f:
+        f.write(value)
+
+def match_score(c1, c2, m, mm, letters, score_matrix):
+	mapped_values = {}
+	for i, v in enumerate(letters):
+		mapped_values[v] = i
+	a = mapped_values[c1]
+	b = mapped_values[c2]
+	return score_matrix[a][b]
+
+def global_alignment(first, second, letters, matrix, protein):
+	match, mismatch, gap, ind = 1,-1,-1,-1
+	if ind == 0 and protein == 0:
+		letters, matrix = matrix_chooser(letters, match, mismatch)
+	n = len(first)
+	m = len(second)
+	backtrack = [[(-1,1) for j in range(m+1)] for i in range(n+1)]
+	s = [[0 for j in range(m+1)] for i in range(n+1)]
+
+	for i in range(1, n+1):
+		s[i][0] = s[i-1][0] + gap
+		backtrack[i][0] = (i-1, 0)
+	for j in range(1, m+1):
+		s[0][j] = s[0][j-1] + gap
+		backtrack[0][j] = (0, j-1)
+
+	for i in range(1, n+1):
+		for j in range(1, m+1):
+			s[i][j] = max(s[i-1][j] + gap, s[i][j-1] + gap, s[i-1][j-1] + match_score(first[i-1], second[j-1], match, mismatch, letters, matrix))
+			if s[i][j] == s[i-1][j] + gap:
+				backtrack[i][j] = (i-1, j)
+			elif s[i][j] == s[i][j-1] + gap:
+				backtrack[i][j] = (i, j-1)
+			else:
+				backtrack[i][j] =  (i-1, j-1)
+
+	first_p = ""
+	second_p = ""
+	i = n
+	j = m
+	while (i,j) != (0,0):
+		if backtrack[i][j] == (i-1, j-1):
+			first_p = first[i-1] + first_p
+			second_p = second[j-1] + second_p
+		elif backtrack[i][j] == (i-1, j):
+			first_p = first[i-1] + first_p
+			second_p = '-' + second_p
+		else:
+			first_p = '-' + first_p
+			second_p = second[j-1] + second_p
+		(i, j) = backtrack[i][j]
+
+	return first_p, second_p, s[n][m]
+
+file1 = 'AGP04929.1.fasta'
+file2 = 'AYN64561.1.fasta'
+
+with open(f'./data/protein/{file1}') as fp:
+    A = ''
+    for name, seq in read_fasta(fp):
+        A += seq
+
+with open(f'./data/protein/{file2}') as fp:
+    B = ''
+    for name, seq in read_fasta(fp):
+        B += seq
+
+protein = 1
+letters, matrix = given_matrices_inserter('pam250.txt')
+resA, resB, score = global_alignment(A, B, letters, matrix, protein)
+write_to_file(f'./outputs/{file1}.txt', resA)
+write_to_file(f'./outputs/{file2}.txt', resB)
+print(f'Score: {score}')

pam250.txt

@@ -0,0 +1,26 @@
+23
+G A V L I P S T D E N Q K R H F Y W M C B Z X
+5 1 -1 -4 -3 0 1 0 1 0 0 -1 -2 -3 -2 -5 -5 -7 -3 -3 0 0 -1 -8
+1 2 0 -2 -1 1 1 1 0 0 0 0 -1 -2 -1 -3 -3 -6 -1 -2 0 0 0 -8
+-1 0 4 2 4 -1 -1 0 -2 -2 -2 -2 -2 -2 -2 -1 -2 -6 2 -2 -2 -2 -1 -8 
+-4 -2 2 6 2 -3 -3 -2 -4 -3 -3 -2 -3 -3 -2 2 -1 -2 4 -6 -3 -3 -1 -8
+-3 -1 4 2 5 -2 -1 0 -2 -2 -2 -2 -2 -2 -2 1 -1 -5 2 -2 -2 -2 -1 -8 
+0 1 -1 -3 -2 6 1 0 -1 -1 0 0 -1 0 0 -5 -5 -6 -2 -3 -1 0 -1 -8
+1 1 -1 -3 -1 1 2 1 0 0 1 -1 0 0 -1 -3 -3 -2 -2 0 0 0 0 -8
+0 1 0 -2 0 0 1 3 0 0 0 -1 0 -1 -1 -3 -3 -5 -1 -2 0 -1 0 -8
+1 0 -2 -4 -2 -1 0 0 4 3 2 2 0 -1 1 -6 -4 -7 -3 -5 3 3 -1 -8
+0 0 -2 -3 -2 -1 0 0 3 4 1 2 0 -1 1 -5 -4 -7 -2 -5 3 3 -1 -8
+0 0 -2 -3 -2 0 1 0 2 1 2 1 1 0 2 -3 -2 -4 -2 -4 2 1 0 -8
+-1 0 -2 -2 -2 0 -1 -1 2 2 1 4 1 1 3 -5 -4 -5 -1 -5 1 3 -1 -8
+-2 -1 -2 -3 -2 -1 0 0 0 0 1 1 5 3 0 -5 -4 -3 0 -5 1 0 -1 -8
+-3 -2 -2 -3 -2 0 0 -1 -1 -1 0 1 3 6 2 -4 -4 -2 0 -4 -1 0 -1 -8
+-2 -1 -2 -2 -2 0 -1 -1 1 1 2 3 0 2 6 -2 0 -3 -2 -3 1 2 -1 -8
+-5 -3 -1 2 1 -5 -3 -3 -6 -5 -3 -5 -5 -4 -2 0 7 0 0 -4 -4 -5 -2 -8
+-5 -3 -2 -1 -1 -5 -3 -3 -4 -4 -2 -4 -4 -4 0 7 10 0 -2 0 -3 -4 -2 -8
+-7 -6 -6 -2 -5 -6 -2 -5 -7 -7 -4 -5 -3 -2 -3 0 0 17 -4 8 -5 -6 -4 -8
+-3 -1 2 4 2 -2 -2 -1 -3 -2 -2 -1 0 0 -2 0 -2 -4 6 -5 -2 -2 -1 -8
+-3 -2 -2 -6 -2 -3 0 -2 -5 -5 -4 -5 -5 -4 -3 -4 0 8 -5 12 -4 -5 -3 -8
+0 0 -2 -3 -2 -1 0 0 3 3 2 1 1 -1 1 -4 -3 -5 -2 -4 3 2 -1 -8
+0 0 -2 -3 -2 0 0 -1 3 3 1 3 0 0 2 -5 -4 -6 -2 -5 2 3 -1 -8
+-1 0 -1 -1 -1 -1 0 0 -1 -1 0 -1 -1 -1 -1 -2 -2 -4 -1 -3 -1 -1 -1 -8
+-8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 -8 1