FASTQファイル - 塩基配列を格納するためのフォーマット

お疲れ様です。9月です。寒いです。

本日はバイオインフォマティクスの基礎の基礎。
FASTQファイルについて見ていきましょう。

ショートリードを生成するシーケンサーだと、だいだい <~200 bpの塩基配列を取得できます。
大抵の場合、シーケンサから出力されたデータは、まずFASTQファイルに変換されます。

なぜでしょう。

世界基準だからですかね。

とりあえずその程度普及しているフォーマットということです。

これは知っとかナイト!です。

Phredクオリティスコア

現在多くのシーケンサでは、Phredクオリティスコアというスコアによってベースコールのクオリティを定義しています。
Phredは元々、ベースコールを行うプログラムだったらしく、
そこで使用されていたクオリティスコアがなかなか優秀であったため、現在でも使用されているとか何とか。

クオリティスコア(Q)の定義

Phredクオリティスコアでは、
任意の塩基のシーケンスクオリティスコア(Q)は、以下の式により定義されます。

Q = -10log_{10}(e)

ここで、eはベースコールが正しく行われない確率です。
つまり、Q値が高いほどエラーの確率が低く、正しいベースコールが行われていることを意味しています。

Q10、Q20、Q30の意味としては以下のようになります。

クオリティスコア ベースコールが正しく行われない確率 (e) 推定ベースコール精度
10(Q10) 10分の1 90%
20(Q20) 100分の1 99%
30(Q30) 1000分の1 99.9%

参考: https://jp.illumina.com/science/education/sequencing-quality-scores.html

FASTQファイル

ではではメインのFASTQの話に移りましょう。

FASTQフォーマットは、ハイスループットシーケンサから出力された塩基データを格納するファイル形式として登場した、
FASTAフォーマットの拡張版になります。

FASTQファイルの中身は、以下のようになっています。
これがリード1本に該当していて、これが何十万リード分含まれています。
つまり行数としては、 4 × (リード数) 行あるファイルということです。

見ていきましょう。

@HWI-ST999:102:D1N6AACXX:1:1101:1235:1936 1:N:0:
ATGTCTCCTGGACCCCTCTGTGCCCAAGCTCCTCATGCATCCTCCTCAGCAACTTGTCCTGTAGCTGAGGCTCACTGACTACCAGCTGCAG
+HWI-ST999:102:D1N6AACXX:1:1101:1235:1936 1:N:0:
1:DAADDDF\<B\<AGF=FGIEHCCD9DG=1E9?D>CF@HHG??B\<GEBGHCG;;CDB8==C@@>>GII@@5?A?@B>CEDCFCC:;?CCCAC

1行目は、なんとなく配列のIDぽいですね。
2行目は、なんとなくリードの塩基配列ぽいですね。
3行目は、1行目と同じぽいですね。 4行目は。。

そうです。これがPhredのクオリティスコアQです。

クオリティスコアQは数値じゃないん?

って思いますよね。ごもっともです。
でもクオリティースコアは2桁だったりするじゃないですか。
塩基と上下で見比べたいのに、縦が揃っていないと見にくいじゃないですか。

そんなわけで、ASCIIコードを使用しているのです。

#0-#32に対応するASCIIコードはプリントできない文字列なので、
#33〜#126の印刷可能な文字をFASTQファイルには使用します。

一方で、クオリティスコアQは0から定義したいので、
FASTQファイル内で使用するASCIIコードを以下のように定義します。

 Fastq ASCII code=Q+33=−10log_{10}e + 33

こうすると、スコア0には!、スコア1には"というように、スコアQとASCII文字が1対1に対応します。
このASCII文字が3行目に記されているということです!

ちなみに3行目は通常1行目と同じなのですが、optionalな情報やコメントなんかを載せるらしいです。
自分はそんなデータに遭遇したことはありませんのであんま気にしなくてよいと思います。