문제해결(PS)/ROSALIND

Character-Based Phylogeny

곰탱이장 2024. 12. 17. 10:00

https://rosalind.info/problems/chbp/

 

ROSALIND | Character-Based Phylogeny

It appears that your browser has JavaScript disabled. Rosalind requires your browser to be JavaScript enabled. Character-Based Phylogeny solved by 229 2012년 7월 2일 12:00:00 오전 by Rosalind Team Topics: Phylogeny Introduction to Character-Based Phyl

rosalind.info

 

 이 문제는 서로 모순되지 않는 복수의 charactor table이 주어진 후, newick 형식으로 tree를 제출하는 문제이다.

필자는 본래 이 문제를 말 그대로 나눈다에 집착하여 tree의 가장 말단에 있는 taxa들을 추려 이 중 지금의 charator table의 taxa가 말단의 taxa의 subset인 것을 찾아 그대로 두개로 나누어 자식 clade로 이어주는 그런 형식을 생각했다. 

from Bio import Phylo
from Bio.Phylo.Newick import Clade,Tree

def making_split(table,species):
    new_split_0 = []
    new_split_1 = []
    for idx,i in enumerate(table):
        if i=='0':
            new_split_0.append(species[idx])
        elif i=='1':
            new_split_1.append(species[idx])
    if len(new_split_0) < len(new_split_1):
        return new_split_0
    else:
        return new_split_1

def find_terminal(new_split,tree):
    for terminal in tree.get_terminals():
        name = terminal.name
        if name in new_split:
            return name

def find_parent_clase(leaf_name,tree):
    for clade in tree.get_nonterminals():
        for child in clade.clades:
            if child.name == leaf_name:
                return clade

def split_and_add(tree,split_set,parent_clade):
    if len(parent_clade.clades) == split_set:
        return
    new_clade_1= Clade(branch_length=None)
    new_clade_2= Clade(branch_length=None)

    for i in parent_clade.clades:
        if i.name in split_set:
            new_clade_1.clades.append(Clade(name=i.name,branch_length=None))
        else:
            new_clade_2.clades.append(Clade(name=i.name,branch_length=None))
    
    parent_clade.clades = [new_clade_1,new_clade_2]

    return
    
with open(r"파일경로",'r') as f:
    species = f.readline().rstrip().split()
    chars=[]
    for i in f.readlines():
        chars.append(i.rstrip())
    root_clade = Clade(branch_length=None)
    root_clade.clades = [Clade(name=species[i],branch_length=None) for i in range(len(species))]

    tree = Tree(root=root_clade)
    print(tree.get_terminals())
    for table in chars:
        new_split=making_split(table,species)
        ter = find_terminal(new_split,tree)
        parent_clade = find_parent_clase(ter,tree)

        split_and_add(tree,new_split,parent_clade)
    
    Phylo.draw(tree)

그래서 처음에는 이렇게 복잡한 형식의 코드르 짰다. Bio.Phylo 모듈을 사용하여 terminal clade를 찾고 그 terminal clade의 부모 clade를 찾고 그 부모 clade를 이용하여 말단의 taxa를 추리고 그 taxa로 charator table의 taxa의 superset을 찾고 이후 나누어 clade를 또 만들어주고,,, 하는 형식의 코드를 짰다. 하지만 이렇게 하니 너무 복잡하였고 Bio.Phylo 모듈은 내가 원하는 형식으로 newick 형식을 출력해주지 않아 불편하고, 알고리즘적으로도 시간이 너무 오래 걸린 좋지 않은 누더기 코드이고 근본적으로 정답도 잘 나오지 않았다.

 그리하여 인터넷에서 정답을 몇 개 찾아본 결과, 가장 마음에 드는 코드를 찾게 되어, 이 코드로 공부하였다.https://github.com/danhalligan/rosalind.info/blob/main/rosalind/bioinformatics_stronghold/chbp.py github에서 찾은 코드이며, 내가 할려고 했던 방식이 top down 형식이라면 이 방식은 bottom up 형식으로 tuple을 이용하여 근본적으로 두 taxa를 하나로 합친다는 발상으로 코드를 짠 형식이다.

def pick_informative(characters):
    for i, ch in enumerate(characters):
        if sum(ch) == 2:#2개민 있음
            return i, 1
        if sum(ch) == len(ch) - 2:#2개만 있음
            return i, 0


def drop_uninformative(characters):
    return [x for x in characters if 1 < sum(x) < len(x) - 1]


def flatten(x):
    if isinstance(x, (list, tuple)):
        return "(" + ",".join(flatten(e) for e in x) + ")"
    else:
        return str(x)


def chbp(names, characters):
    while characters:
        i, v = pick_informative(characters)#일단 두 개 있는거 부터 시작
        ind = [i for i, x in enumerate(characters[i]) if x == v]#2개 짜리의 인덱스를 2개 특정
        names[ind[0]] = tuple(names[i] for i in ind)#2개를 합침
        del names[ind[1]]#합친거 없앰
        for ch in characters:
            del ch[ind[1]]#캐릭터 테이블에서도 없앰
        characters = drop_uninformative(characters)#이제 쓴거 버림
    return tuple(names)

def main(file):
    names, *characters = open(file).read().splitlines()
    names = names.split()
    characters = [[int(x) for x in list(ch)] for ch in characters]
    res = chbp(names, characters)#res는 튜플로
    wf=open(r"파일경로",'w')
    print(flatten(res), ";", sep="",file=wf)

main(r"파일경로")

tuple형식의 names를 이용하여 charator table과 names를 직접 고쳐가며, 두 개씩 괄호로 묶어서 만들어내는 형식이다. 그리고 그렇게하여 나온 tuple은 DFS식의 함수를 사용하여 문자열이 나올때까지, 계속 수정하여, newick 형식의 tree를 출력하였다. 필자는 개인적으로 Bio.Phylo module을 맹신하였으나, 어쩌면 직접 다루는 것이 알고리즘적으로 더 편할 수도 있겠다는 생각을 하게 되었다.

'문제해결(PS) > ROSALIND' 카테고리의 다른 글

Counting Quartets  (0) 2024.12.17
Encoding Suffix Trees  (0) 2024.11.23
Using the Spectrum Graph to Infer Peptides  (0) 2024.11.23
Quartets  (0) 2024.11.19
Matching a Spectrum to a Protein  (1) 2024.11.18