ZGA pipeline에 두 쌍의 paired end read FASTQ 파일을 투입한 결과가 좀 이상하다

서로 다른 종에 속하는 미생물 균주로부터 유래한 일루미나 sequencing read를 한꺼번에 조립하면 결과는 엉망이 된다. 특별히 metagenomic assembler를 쓰지 않은 상황이라면 이런 불량한 결과가 나오는 것은 너무나 당연하다. 가끔은 일부러 테스트를 위해 이런 조립을 해 보기도 한다. 

매우 유사하지만 ANI 기준으로 분명히 서로 다른 종에 속하는 두 균주에서 생산한 일루미나 sequencing read가 있다. 전부 오래전에 NCBI SRA에 등록해 놓았는데 그 분량이 너무나 많은 관계로 다운로드 후 각각 150x 분량으로 서브샘플링을 하였다. 

• 균주 A: A_1.fastq, A_2.fastq
• 균주 B: B_1.fastq, B_2.fastq

SPAdes로 조립을 하려면 다음과 같이 입력해야 한다.

spades.py --threads 16 -o outdir \
--pe-1 1 A_1.fastq --pe-2 1 A_2.fastq \
--pe-1 2 B_1.fastq --pe-2 2 B_2.fastq

혼동하기 쉽게 생겼다. 예전 방식(deprecated syntax)으로는 '--pe1-1 A_1.fastq --pe1-2 A_2.fastq --pe2-1 B_1.fastq --pe2-2 B_2.fastq'의 형태로 입력하란다. Library ID number에 해당하는 숫자(빨간색)의 위치가 미묘하게 다르다. 

이렇게 조립을 하면 당연히 결과는 엉망이 된다. Contig의 수는 많고, 길이는 매우 짧으며, total length도 크게 늘어난다. 그러면 ZGA pipeline을 쓸 때에는 어떻게 명령어를 조합해야 하는가? Assembly engine으로는 SPAdes를 사용한다. 

ZGA 웹사이트의 설명을 보자.

zga -1 Lib1.R1.fq.gz Lib2.R1.fq -2 Lib1.R2.fq Lib2.R2.fq combination of reads from two sequencing libraries

SPAdes보다는 명령어 라인 구성이 간결해서 좋다. 라이브러리의 순서를 지켜서 forward read file은 -1 뒤에, reverse read file은 -2 뒤에 나열하면 된다. 그런데 결과가 너무나 좋다! A 유래 샘플만 사용한 것과 거의 다르지 않은 것이다! 이건 말이 되지 않는다. 그래서 다음과 같이 두 파일을 합쳐서 하나로 만든 뒤 이를 이용하여 ZGA를 실행하였다.

cat A_1.fastq B_1.fastq > C_1.fastq
cat A_2.fastq B_2.fastq > C_2.fastq
zga -o output -1 C_1.fastq -2 C_2.fastq

심하게 조각난 assembly가 산출되었다. 이렇게 나오는 것이 정상이다. ZGA 소스 코드에 약간의 버그가 있는 것 같다. 그것이 아니라면 나의 착각 또는 실수?

SPAdes version 3.14.0부터 도입된 '--isolate' 옵션을 쓰는 이유나 철저하게 이해하도록 하자. metaSPAdes를 돌리는 것이 아닌데도 '--isolate' 옵션을 따로 준다는 것이 잘 받아들여지지는 않지만 말이다.

Special --isolate option for assembly of standard datasets with good coverage (>100x). If you have high-coverage data for bacterial/viral isolate or multi-cell organism, we highly recommend to use --isolate option.