- Reworked IonHammer module
- Major improvements in metaSPAdes pipeline
- Support for strand-specific data and coverage-based isoform detection in rnaSPAdes
- exSPAnder repeat resolution is improved in terms of both performance and repeat resolution accuracy
- Better RAM consumption and running time in general
- Improved careful mode for isolate datasets
- Several other important bug-fixes
일루미나와 PacBio 등 다양한 시퀀싱 플랫폼을 지원할뿐만 아니라 RNA-seq과 메타게놈까지 그 적용 대상이 확대되고 있으며 지속적으로 개선이 되는 모습을 보면 참 부럽다는 생각이 든다. 흔히 IT 강국으로 불리는 우리나라에서는 왜 전세계적으로 사랑을 받을 killer app이 나오지 않는 것을까? 카카오톡이나 라인 정도?
내가 생명정보학 분야를 전부 다 파악하지는 못하지만 주로 관심을 갖는 genome assembly 영역만 보아도 전세계적으로 골고루 특징적인 도구를 꾸준히 내놓고 있다. phrap, cap3, Celera assembler, velvet, SPAdes, SOAPdenovo, IDBA, 그리고 최근의 HGAP, falcon, canu에 이르기까지. 중국과 대만도 이렇게 기여를 하고 있는데 우리나라는 없다. 참으로 안타까운 현실이다.
VelvetOptimiser 새소식
나도 모르던 사이에 VelvetOptimiser가 GitHub로 이사를 갔다. 현재 배포되는 위치는 https://github.com/tseemann/VelvetOptimiser이다. Linuxbrew를 사용해서 설치를 해 보니 원래 쓰던 것과 같은 2.2.5 버전이었다. 게다가 딸려온 velvet binary도 최대 kmer가 매우 낮은 기본 조건으로 컴파일된 상태였다. 어쩐지 테스트 실행을 할 때 설정한 kmer size와 무관하게 31로 시작을 하더라니... 지워버리고 말았다. 나중에 linuxbrew로 재설치를 하게 된다면 시스템에 이미 깔린 velvetg, velveth를 쓰도록 매만져야 하겠다. 가장 간단하게는 ~/.linuxbrew/bin/ 이하에 설치된 velvetg와 velveth를 지워버리면 된다.
새삼스럽게 왜 velvet인가? 요즘 미생물 유전체 시퀀싱 데이터를 관리하는 DB를 구상하는 중인데, library insert size를 계산하는 가장 간편한 방법이 무엇인지 찾고 있었다. De novo assembly report에는 결국 대부분 나오게 되는 수치이지만 개별 시퀀싱 샘플마다 전부 다른 방법을 쓴 관계로 이를 일일이 다시 찾아내려면 여간 성가신 일이 아니기 때문이다. 조립을 통해 만든 contig 위에 read를 다시 매핑한 뒤 getinsertsize.py를 실행하여 정확한 수치를 알아내는 것이 가능하나 reference를 다시 마련해야만 한다. 차라리 '가벼운' de novo assembler를 일괄적으로 한번 더 돌리는 것이 낫다. 가볍다는 것의 의미는 트리밍이나 오류 교정 등의 전처리를 되도록 하지 않는 것을 의미한다. 정확한 assembly를 얻고자함이 아니라서 문제가 될 것은 없다. 대신 velvet은 interleaved file을 입력물로 제공해야만 한다.
VelvetOptimiser를 실행할 때 -v(--verbose) 옵션을 주면 각 라이브러리에 대해서 추정된 insert length와 표준편차를 logfile에 기록한다. 실제 실행을 해 보니 -v 옵션을 주지 않아도 알아서 기록이 된다. library 크기와 관련된 수치가 왜 두 줄로 표시되는지는 잘 모르겠다.
Paired Library insert stats:
Paired-end library 1 has length: 609, sample standard deviation: 191
Paired-end library 1 has length: 609, sample standard deviation: 192
IDBA의 활용은 어떨까?
이번에는 심플한 assembler로서 IDBA를 선택해 보았다. IDBA는 interleaved fasta 파일을 입력물로 받아들인다. 따라서 fastq file의 쌍 또는 interleaved fastq file을 갖고 있다면 설치_디렉토리/bin/fq2fa를 사용하여 전환해야 한다. 'N'을 포함한 read를 제거하는 --filter 옵션은 필수는 아니다.
IDBA 역시 기본 조건으로 컴파일하면 다룰 수 있는 read length에 한계가 있어서(128 bp) MiSeq data를 쓰기에 곤란하다. 이러한 상황에서는 소스의 헤더 파일을 바꾸어서 다시 컴파일해야 된다(SEQanswers 250bp reads in idba_ud). 이 문제는 버전 1.1.1에서는 해결되었다.
실행 중에 화면으로 표시되는 결과는 고스란히 output_dir/log에도 기록된다. 설정된 kmer에 대한 조립을 마친 뒤 최적 결과를 선택하고 이어서 라이브러리 관련 수치가 나온다. 실행 속도는 velvet보다 조금 빠른 것처럼 느껴진다.
distance mean 583.363 sd 147.819
Velvet과 IDBA의 실행 속도 비교
'N' 제거용 필터를 거치지 않은 동일한 interleaved fasta file을 사용하여 velvet와 idba의 실행 속도를 비교하여 보았다. kmer의 크기는 91로 고정하고 thread의 수도 1로 하였다. 기본 thread의 수는 velvet의 경우 4, IDBA는 최대 가용 CPU이다.
$ VelvetOptimiser.pl -s 91 -e 91 -f '-shortPaired -fasta M217.fa' -t 1
...
Paired-end library 1 has length: 604, sample standard deviation: 185
...
real 35m53.620s
user 33m10.259s
sys 1m42.553s
$ idba -o idba_out -r M217.fa --mink 91 --maxk 91 --num_threads 1
...
distance mean 584.907 sd 147.223
...
real 17m26.175s
user 17m21.627s
sys 0m4.952s
IDBA는 입력파일을 fastq에서 fastq로 전환하는 수고를 들여야 하지만 실행 속도가 더 빠르다. 특히 thread 수 관련 옵션은 아무것도 지정하지 않으면 기본적으로 최대 가용한 것을 다 쓰게 되므로 속도가 더 빠르다. 그런데 왜 두 프로그램이 제공하는 수치가 다를까? IDBA는 paired read 사이의 distance를 계산하는 반면 velvet은 library size라서 그런가?
원래 (paired reads or insert) distance는 서로 마주보고 위치한 두 read의 내부 거리를 의미한다. 따라서 insert size = distance + 2 rd (rd: read length)가 되어야 한다. 그러나 위에서 얻은 두 수치는 이를 만족하지 못한다. 겨우 20 bp밖에 차이가 나지 않는다. 무엇을 믿어야 할까? mapping 이후의 distance 계산은 어쩌면 일부의 read에 대해서 샘플링을 하여 산출하는지도 모른다.
Simulated Illumina read를 사용하여 계산을 한다면?
Read simulation tool로 잘 알려진 ART(MountRainier-2016-06-05)를 사용하여 가짜 read를 만든 뒤 이를 velvet와 IDBA로 각각 조립해 보면 어느 프로그램이 좀 더 정확한 insert length 관련 수치를 계산해 내는지 알아낼 수 있을 것이다. 평균 fragment size가 600 bp이고 표준편차가 150 bp인 HiSeq 2500 모델의 paired read를 만들어서 작업을 해 보았다.
$ art_illumina -ss HS25 -i /data/project/00_references/NC_007492.fa -l 150 -f 200 -o simul -m 600 -s 150
$ fq2fa --merge simul1.fq simul2.fq simul.fa
$ time VelvetOptimiser-2.2.5/VelvetOptimiser.pl -s 91 -e 91 -f '-shortPaired -fasta simul.fa' -t 1
...
Paired-end library 1 has length: 599, sample standard deviation: 149
Paired-end library 1 has length: 600, sample standard deviation: 149
...
real 24m57.445s
user 16m30.927s
sys 4m21.573s
$ time /usr/local/apps/idba-1.1.1/bin/idba -o idba_out -r simul.fa --mink 91 --maxk 91 --num_threads 1
...
distance mean 599.669 sd 143.797
...
real 15m47.285s
user 13m30.796s
sys 1m21.300s
이번에는 프로그램이 산출한 수치가 거의 같다. 실행 속도는 IDBA가 더 빠르며, assembly 수치는 velvet이 더 우수하다(contig의 수가 적음). 따라서 library size 계산을 빨리 하려면 IDBA가 적당할 것이다.
댓글 없음:
댓글 쓰기