지난 일주일 동안 실험동물의 장에서 분리된 것으로 알려진 어느 미생물의 일루미나 시퀀싱 결과물을 조립하느라 무척 애를 먹었다. 16S rRNA gene sequence가 여러 조각이 난 채로 조립이 되고, contig의 수마저 1400개 정도에 육박하는 아주 골치 아픈 상황이었다. 두 가지 문제는 서로 다른 곳에 원인이 있었는데, 문제 발생의 원인이 실험 단계(균주를 분리 및 배양하고 시퀀싱용 genomic DNA를 추출하는 사람)에 있었는지 혹은 시퀀싱 업체에서 있었는지를 빠르고 정확하게 파악하기가 쉽지는 않다.
꽤 오래전에 같은 연구소에 근무하는 선배 과학자의 점액세균 샘플 자료를 분석한 일이 있다. 조립을 하고 나니 상상을 초월하는 수의 짧은 contig가 나왔고, BLAST를 해 보니 사람과 소 등 미생물 배양 실험 환경에서는 유입되기 어려운 생명체의 genomic fragment가 확인되었다. 이것이 미생물을 배양하고 DNA를 뽑던 사람에 의해 유입되었을 가능성은 극히 낮다. 다른 고객이 의뢰한 DNA 샘플을 처리했던 시퀀싱 업체의 장비 세척 불량이라는 강한 의심을 갖는 것도 당연하지 않겠는가?
k-mer frequency를 기반으로 하여 read를 전처리하면 그 원인이야 어떻든 간에 낮은 빈도로 존재하는 오염 유래 read를 제거하게 되어 양호한 조립 결과를 얻을 수가 있다. 혹은 read subsampling(~150X)을 하면 오염 read의 빈도가 현저히 떨어지므로 마찬가지로 좋은 결과를 얻을 수 있기도 하다. 이번 프로젝트에서는 sequencing depth가 600X를 넘는 수준이었기 때문에 아무 생각 없이 subsampling을 진행했더라면 결과에 별로 신경을 쓸 필요가 없었을 것이다(실제로 시도해 보았으나 약간의 개선은 있었지만 만족할 수준은 아니었다). 그러나 오염의 구체적인 원인을 알아야 미생물을 다루는 실험자에게 개선 사항을 전달할 수 있기 때문에 이 문제를 분석하는데 집착을 하지 아니할 수가 없었다.
시퀀싱 업체에서 발생한 오염이라고 가정해 보자. (1) 바로 직전에 러닝을 했던 라이브러리 샘플이 관로 내에 남아서 read contamination을 유발한 것인가, 혹은 (2) 다른 바코드를 붙여서 내 샘플과 동시에 러닝을 했는데 바코드 해독이 불충분해서 내 샘플의 데이터에 섞인 것인가? 아마 일루미나 장비를 직접 가동하는 기술자라면 (1)이나 (2) 중 어느 것이 원인일지 잘 알고 있을 것 같다. 내 생각에는 (1)이 주요한 원인이 되지 않을까 싶다. 아주 단순한 수준의 샘플 간 cross-contamination이 일어날 수도 있을 것이다.
6년 전에 이 문제에 관한 논문을 bioRxiv(doi) 에 실은 일이 있었다. bioRxiv라는 프리프린터 서버에 논문을 게재한 경험에 대한 기록도 내 블로그에 남아 있다. 이번 장내 미생물 조립 불량 사태를 겪으면서 즉각적으로 오염이 있을 것이라는 생각은 했지만, 단순히 colony isolation이 불충분해서 그랬을 것이라고 생각했었다. 그러나 오염원은 아주 엉뚱한 것이었다. 지난 3년 동안 인류를 괴롭혔던 바로 그 병원체의 게놈이었다!
(업데이트: 오전에는 그렇게 생각했었는지 오후에는 생각이 달라졌다. smalt가 너무 sensitive해서 빚어진 오해라고 여겨진다. 왜냐하면 smalt로는 의심스런 오염원 DNA 염기서열에 대한 read mapping이 이루어졌지만 bowtie2나 bbmap으로는 그렇지 않았기 때문이다. 오염원이 정확히 뭔지 확신하기가 어렵다...)
Dataset에서 오염에 해당하는 read를 사전 제거하기 위해 smalt를 사용하였다. Smalt는 영국의 Sanger Institute에서 개발한 read mapper로서 SourceForge에서 제공하는 최신 버전은 무려 8년 전에 나온 v0.7.6이다. 이 mapper program은 BWA나 bowtie2만큼 인기가 있지는 않지만 사용법이 단순하고(매뉴얼) 박테리아 게놈에는 매우 잘 작동하므로 나를 포함한 일부 매니아가 사용한다.
[Stack Exchange] Why is SMALT is better for microbial genomics than other mappers?
미생물 유전체의 조립을 마치고 나서 특정 reference에 대하여 매칭하는 contig가 무엇인지를 찾을 때에는 주로 cross_match를 사용한다. 이것 역시 대단히 고전적인 프로그램이다. Sanger dideoxy chain termination chemistry에 의하여 유전체 프로젝트를 진행하던 수십 년 전에 개발된 'word-nucleated banded Smith-Waterman comparison'의 구현으로서, 속도는 느리지만 결과 요약줄이 매우 직관적이라서 나는 이를 여전히 즐겨 사용한다. cross_match의 결과를 시각화하는 도구인 XMatchView라는 것도 있다. Cross_match는 오픈소스가 아니므로 이를 이용하려면 라이선스 계약을 맺어야 한다. 내가 (주)제노텍에 근무하던 시절, Phred/Phrap(+cross_match)/Consed의 라이선스를 구입했던 것이 아마 2000년이었을 것이다. 당시 1천만원 정도의 비용이 들었고, 정부출연연구소로 이직을 하면서 다행스럽게도 비용이 들지 않는 라이선스 계약을 맺을 수 있었다. 지금 쓰는 바이너리는 언제 컴파일한 것인지는 모르겠지만 버전 1.080812이다. Phred, Phrap과 Consed는 3GS 시대로 넘어오면서 이제는 사용하지 않게 되었으나, 일종의 범용 도구라고 할 수 있는 cross_match는 여전히 내 컴퓨터에서 작동하고 있다. MUMmer v3도 그런 성격의 고전적인 프로그램 아니겠는가?
Muscle와 usearch의 개발자 Robert Edgar는 2020년에 URMAP이라는 고속 read mapper를 개발하였다고 한다. 개발자에 의하면 URMAP은 정확도를 희생하지 않으면서도 BWA나 bowtie2보다 'an order of magnitude faster'라고 하니 조만간 써 봐야 되겠다. R. Edgar는 소속이 없이 독립적으로 일하는 사람이라서 bioRxiv 논문(링크)의 첫 쪽에도 'Independent Investigator'라고 소개를 해 놓았다. 참 멋있지 않은가? 이 논문은 2020년 6월에 PeerJ에 출판된 상태이다.
학술지에 논문을 발표하는 저자가 꼭 소속을 갖출 필요가 없음을 잘 보여주는 R. Edgar의 사례. 'Unaffiliated', 멋지지 않은가? 출처 링크. |
이번 작업을 하면서 오염에 해당하는 read를 제거한 뒤 pair를 맞추기 위해 Brian Bushnell의 repair.sh 스크립트(링크 @JGI)를 자주 사용하였다. B. Bushnell의 read mapper인 BBMap도 가끔씩 사용한다. 그가 만든 프로그램 모음인 BBTools를 종종 사용하면서 종종 감탄을 하게 된다. NGS 자료 처리를 위한 이렇게 다양한 유틸리티 모음을, 그것도 완성도가 아주 높은 결과물을 개인이 개발해 내다니!
파견 근무 3개월 차에 접어들면서 연구 업무(미생물 유전체의 NGS 분석)에 대한 감을 잃을까 걱정을 하고 있었는데, 이번에 이렇게 일을 몰아치기로 하면서 아직 나의 칼날이 무뎌지지는 않았다는 안도감을 되찾았다. 연구소에서는 쓰지 못하던 Oracle VirtualBox를 다시 자유롭게 사용하면서 사용 편의성이 예전보다 정말 좋아졌다는 느낌도 갖게 되었다. 연구소에서 VirtualBox 웹사이트 접근이 차단된 사정은 여기에 적어 놓았다. 하지만 파견 근무지의 전산 환경에서는 Windows Subsystem Linux가 설치되지 않는 기이한 현상이 벌어지고 있다. Microsoft Store도 설치되어 있지 않았고, 이를 어떻게 가능하게 만들지도 알기 어렵다. 보안 정책 때문일까? 모든 것이 다 원하는 대로 돌아갈 것을 기대하면 곤란할 것이다.
2022년 11월 8일 업데이트
Technology dictates algorithms: recent developments in read alignment. Genome Biology volume 22, Article number: 249 (2021)
댓글 없음:
댓글 쓰기