Suffix trie
#
Find similar titles
- (rev. 3)
- JehongLee
Structured data
- Category
- Algorithm
접미사 트라이(Suffix trie)의 활용 #
NGS sequencing의 결과 적게는 수 천만 read 에서 수 억개의 read 가 생산되게 된다.
이 수많은 read를 genome에 align하기 위해 툴들이 가장 많이 사용하는 방법이 genome을 suffix-trie 자료구조로 변환하는 것이다. 변환된 genome의 자료구조를 바탕으로 read에 matching 시켜 위치를 찾는다.
아래 코드는 "" AATCGGGTTCAATCGGGGT ATCG GGGT "" 를 input으로 받을 경우 첫번째 string을 genome으로 인식하여 suffix-trie 자료구조를 만들고 두번째 행 이후로의 string을 pattern 혹은 read라고 생각하여 genome 상의 위치를 출력하는 코드이다.
class Node:
def __init__(self, name):
self.name = name
self.posterior_node = {}
self.pos = None
def set_posterior_node(self, node, edge):
self.posterior_node[edge] = node
def get_posterior_node(self, edge):
try:
return self.posterior_node[edge]
except KeyError:
return -1
def __str__(self):
return str(self.name)
def set_pos(self, pos):
self.pos = pos
def get_pos(self,pos):
return self.pos
def depth_first(node):
for edge in node.posterior_node:
posterior_node = node.posterior_node[edge]
print('{}->{}:{}'.format(node.name, posterior_node.name, edge))
depth_first(posterior_node)
def check(node, text):
next_node = node.get_posterior_node(text[0])
if next_node == -1:
return False
elif next_node.get_posterior_node('$') != -1:
return True
else:
return check(next_node, text[1:])
def main(file_):
global number
number = 0
with open(file_) as f:
lines = f.readlines()
text = lines[0]
patterns = [line.strip() for line in lines[1:]]
first_node = Node(number)
number += 1
## make trie
suffix_pos = 0
for pattern in patterns:
for i in range(len(pattern)):
edge = pattern[i]
if i == 0:
prior_node = first_node
else:
prior_node = next_node
next_node = prior_node.get_posterior_node(edge)
if next_node == -1:
next_node = Node(number)
number += 1
prior_node.set_posterior_node(next_node, edge)
last_node = Node(number)
next_node.set_posterior_node(last_node, '$')
number += 1
#depth_first(first_node)
## find pattern match position
result_match_pos = []
for i in range(len(text)):
if check(first_node, text[i:]):
result_match_pos.append(i)
print ' '.join([str(pos) for pos in result_match_pos])
if __name__ == '__main__':
main(sys.argv[1])