2010년 9월 15일 수요일

Illumina-Solexa read 다루기

Sanger Institute에서 제안한 FASTQ format은 서열과 quality score를 하나의 파일에 축약해서 표연하는 방식으로서 high-throughout short read data를 다루기에 적합하며, 많은 de novo assembler 또는 reference mapper 프로그램이 입력 파일의 형태로서 요구하고 있다. 원래 quality score(Phred score)는 십진수 0~93으로 로 표현되지만, 이를 FASTQ format에서는 인쇄 가능한 아스키 문자로 대체하고 있다. 98, 99와 같은 높은 값의 quality score는 사용자가 편집한 값을 표현하는 용도로 쓰인다.

ASCII 32는 공백이므로, 인쇄 가능한 형태의 문자로 표현되는 ASCII 33~126이 FASTQ에서의 quality score로 쓰인다. 그러면 Phred score를 FASTQ format으로 표현하려면 어떻게 해야 하는가? 93보다 큰 값은 93으로, 이보다 작은 값은 그대로 둔 뒤 33을 더한 다음, 이 ASCII 값에 해당하는 문자로 변환하면 된다. 요약한다면 0~93까지의 Phred score를 33~126 범위의 숫자로 변환하여, 문자로 전환한다는 뜻이다.

이를 간단한 Perl code로 표현해 보자. Phred quality를 $Q, Sanger FASTQ의 값을 $q라 하면,

$q = chr(($Q<=93? $Q : 93) + 33); $Q = ord($q) - 33; Solexa/Illumina platform에서 생산되는 값은 조금 달라서, 0~62 범위의 값을 만들어 낸다. 이 경우에는 다음의 식을 적용하라. $Q = 10 * log(1 + 10 ** (ord($sq) - 64) / 10.0)) / log(10); http://maq.sourceforge.net/fastq.shtml에서 참조.

ord() 함수는 인수로 주어진 문자의 ASCII 값을 반환한다. Nucleic Acids Research 38(6) 1767-71에 Sanger FASTQ file format과 Solexa/Illumina FASTQ variant에 대한 상세한 내용이 나와있다.

maq를 이용해서 Illumina-Solexa 결과를 Sanger type의 FASTQ로 바꿀 수 있다.

maq sol2sanger s_1_sequence.txt s_1_sequence.fastq

s_1_sequence.txt가 Solexa read sequene file이다.

[2012년 7월 11일에 추가 작성]

quality score를 ASCII 값으로 전환하기 위해 더해주는 기본 값이 Sanger standard FASTQ(64)와 Illumina FASTAQ(33)가 서로 달랐었지만, CASAVA(Illumina pipeline) 1.6부터인가는 64로 통일되었다. 따라서 2011년 봄 이후로 생산되는 일루미나 데이터에 대해서는 quality encoding의 호환성에 대해 고민할 필요가 없어졌다. 매우 다행스런 일이다.

댓글 없음: