2020년 11월 17일 화요일

세균을 위한 최적의 SNP calling 방법

몇 메가염기쌍도 되지 않는 세균의 유전체 시퀀싱 데이터를 다루는 것은 흔히들 '간단한 일'로 인식하기 쉽다. 썩 틀린 말은 아니다. 그러나 여러 isolate의 변이체 분석을 하는 것은 인간 게놈을 다루는 것보다 더 까다로운 측면이 있다. 연구 대상 종에 따라 적용되는 biology가 전부 다르고, 같은 종에 속하는 isolate만을 다룬다 해도 최적의 reference genome을 선정하는 문제도 쉽게 풀리지는 않기 때문에 체계적인 벤치마킹은 인간 변이체에 대한 것보다도 그 사례가 더 적다고 보는 것이 맞을 것이다.

나는 습관적으로 CLC Genomics Workbench를 쓰고 있지만 이러저러한 사정으로 2018년에 마지막 업데이트를 하였다. 벤치마킹 사례에서 상용 소프트웨어에 대한 성능은 잘 다루어지지 않는다. 

Sequencing read로부터 phylogenetic analysis에 쓸 수 있는 alignment(core or SNP) 자료를 뽑고자 종종 snippy를 사용한다. Prokka, abricate, abricate 등 정말로 쓸만한 소프트웨어를 만들어 내는 Torsten Seemann(GitHub)을 존경한다!

그때 그때 기분에 따라서 snippy나 harvest suite를 쓰다가, 좀 더 새로운 시도를 해 보자는 마음을 먹었다. 이 부류의 종합 프로그램으로는 비교적 역사가 오래 된 NASP(Northern Arizona SNP Pipeline)이나, 아주 최근에 나온 PhaME 등을 검토하다가 일본 연구자들이 개발한 BactSNP라는 것을 만나게 되었다.

Evaluation of SNP calling methods for closely related bacterial isolate and a novel high-accuracy pipeline: BactSNP. Microbial Genomics 2019:5 [PubMed]

지금은 없어졌지만, 전주 영화의 거리에 [맛집 찾다가 내가 차린 집] 혹은 이와 비슷한 이름의 식당이 있었다. BactSNP가 개발된 배경도 그러하다. 기존의 도구로는 outbreak의 조사처럼 매우 가깝게 연관된 isolate를 비교할 때 false-positive SNP가 많이 나오는 것을 극복하기 어려워서 아예 파이프라인을 새로 만들었다는 것이다.

이 소프트웨어는 Platanus(PLATform for Assembling NUcleotide Sequence - 어떻게 여기에서 Platanus라는 acronym을 생각해 냈을까!) assembler를 개발한 팀에서 만든 것 같다. Platanus assembler 웹사이트에 BactSNP가 같이 소개되어 있기 때문이다. 중국의 SOAPdenovo, 러시아의 SPAdes, 홍콩(이제는 중국이구나!)의 idba, 일본과 중국의 합작품인 MEGAHIT... 이런 유명한 소프트웨어와  BactSNP가 어깨를 나란히 한다고 생각하니 문득 우리나라에서는 곧바로 떠오르는 생명정보학의 killer app이 없다는 현실이 서글프게 느껴졌다. 

내가 과문해서 그런가... 그렇다면 다행이고.

BactSNP에 관심을 갖는 것은 실행 명령어가 너무나 단순하기 때문이었다. 일루미나 read만을 사용한다면 snippy-multi에서 정의한 input.tab 파일을 그냥 써도 된다. 아시아서 만들어진 프로그램이 보통 그렇듯이 사용 설명은 그렇게 친절하지 않다.

SRA에서 다운로드한 raw data 몇 개를 넣고 테스트 러닝을 시작하였다. read의 트리밍을 하는가 싶더니 조립을 한다! 왜지? 처음부터 매핑을 하지 않고 조립이라니? BactSNP의 세부적인 방법은 Supplementary Material에 설명이 되어 있다. Simulated dataset의 생성과 다른 소프트웨어의 비교를 위한 명령행까지 아주 친절하게 잘 수록해 놓아서 공부하기에는 아주 좋다. GitHub 사이트의 다소 부실한(?) 설명이 이것으로 보상이 된다.

Supplementary Material의 설명에 의하면 각 샘플의 read를 조립(platanus 이용)하여 contig를 생성한 다음, nucmer를 사용하여 이를 reference 위에 1-to-1으로 위치시킨 뒤 pseudogenome을 만든다. 다음의 그림을 참조하면 pseudogenome이 어떻게 만들어지는지를 쉽게 이해할 수 있다. Reference의 coordinate를 뼈대로 하되 실제 염기서열은 sample read에서 조립된 contig인 것이 pseudogenome이다. Reference와 contig는 서로 차이가 나는 곳이 존재하므로, 세부적으로 고려할 곳이 생긴다. 예를 들어 reference 서열에 대한 insertion은 무시되고, Contig가 cover하지 않는 곳은 '-'으로 치환된다. Contig가 reference와 대하여 잘 align하다가 어느 위치부터 맞지 않게 된다면(contig 특이적인 insertion이 아닌 경우)? 그런 복잡한 사항은 supplementary material을 잘 읽어보라.


이어서 BWA-MEM으로 reference에 대한 read mapping을 실시한다. PCR 과정에서 만들어진 duplicated read는 MarkDuplicates로 제거하고, 신뢰하기 어려운 allele은 pseudogenome에서 가려진다. 결정된 SNP 정보는 snps_w_ref.tsv와 snps_wo_ref.tsv라는 tab-delimited file에 기록된다. w_ref와 wo_ref는 각각 reference 서열을 포함하거나('w' for with) 포함하지 않음('wo' for without)을 의미한다. 파이프라인 초기 단계에 작성된 pseudogenome에서 ambiguous base로 남아있던 위치는 snps_wo_ref.tsv에서 결정된 SNP로 대체된다. replaced의 의미는 바로 이것이었다. 계통발생학적 분석(예: gubbins와 tree inference 등)은 바로 이 replaced_pseudo_genomes_wo_ref.fa를 이용하면 된다. 

Method 페이퍼라서 그런지 전체 방법에 대한 간략한 설명과 결과만 본문에서 강조하고 상세한 설명은 전부 별도의 파일로 돌려 놓은 것이 아쉽다. 왜 요즘 논문은 인쇄한 몇 페이지의 본문을 읽는 것으로 해결이 되지 않을 정도로 점점 대형화하는지.. 저자도 많아지고 공동 주저자도 많아지고...


댓글 없음: