シークエンサーは、フローセル上で塩基配列を伸長させるときに生じた光を補足し、その光の色と強さを機械学習モデルとに入力し、シークエンシングエラーを予測（計算）している。その計算過程はブラックボックスと言われることが多い。シークエンサーによって読み込まれた塩基配列とそのときのシークエンシングエラーが一通りの処理を経て、FASTQ ファイルに変換される。なお、FASTQ ファイルのなかで、シークエンシングエラーは、下記の式に基づいて、クオリティスコアに変換されて記載される。
\[ Q = -10\log_{10}p_{error} \]FASTQ ファイル中のクオリティスコアは、次のように、数値ではなく、文字で記述されている。各エントリーの 4 行目に書かれている文字列がクオリティスコアになる。
@SRR1170086.1 HWI-ST845:120525:D10G7ACXX:8:1101:9728:1999 length=50
GTTTTTAAAATGAGTTTGCAAATATTACTGTATTTTTNTCCCATGCTTTT
+SRR1170086.1 HWI-ST845:120525:D10G7ACXX:8:1101:9728:1999 length=50
11=DDFFFFHHGFIGGIIJJJJIJJJJJJIIIIIJJJ#1?DHIJJIJIII
@SRR1170086.2 HWI-ST845:120525:D10G7ACXX:8:2316:19897:100799 length=50
CGATGGAATAGATTTCTCCAAGTTAGTTGGAGGCTAGTCTTCTCTACATA
+SRR1170086.2 HWI-ST845:120525:D10G7ACXX:8:2316:19897:100799 length=50
???DADDA:++22C::4:::AA33C,,,2EFECB;ED9D4?4*0?D0B9Bこれら英数字で記述されているクオリティスコアを数値に変換するには、その英数字の ASCII コードを調べ、その ASCII コードから 33 を引く計算を行う。
\[ Q = (ASCII code) - 33 \]例えば、1 番目のリード（SRR1170086.1）の最初の塩基は G で、そのクオリティスコアは 1 である。1 の ASCII コードを調べると 49 であることがわかる。このとき、クオリティスコアは 49 - 33 = 16 となる。つまり、1 番目のリードの最初の塩基 G のクオリティスコアは 16 である。
51 3 71 G 91 [ 111 o
52 4 72 H 92 \ 112 p
33 ! 53 5 73 I 93 ] 113 q
34 " 54 6 74 J 94 ^ 114 r
35 # 55 7 75 K 95 _ 115 s
36 $ 56 8 76 L 96 ` 116 t
37 % 57 9 77 M 97 a 117 u
38 & 58 : 78 N 98 b 118 v
39 ' 59 ; 79 O 99 c 119 w
40 ( 60 < 80 P 100 d 120 x
41 ) 61 = 81 Q 101 e 121 y
42 * 62 > 82 R 102 f 122 z
43 + 63 ? 83 S 103 g 123 {
44 , 64 @ 84 T 104 h 124 |
45 - 65 A 85 U 105 i 125 }
46 . 66 B 86 V 106 j 126 ~
47 / 67 C 87 W 107 k
48 0 68 D 88 X 108 l
49 1 69 E 89 Y 109 m
50 2 70 F 90 Z 110 nクオリティスコアは、シークエンシングエラーから変換されているため、その逆変換を行えば、クオリティスコアからシークエンシングエラーを計算することができる。例えば、クオリティスコアが 16 のときのシークエンシングエラーは次の様に計算できる。
\[ p_{error} = 10 ^ {-\frac{Q}{10}} = 10^{-\frac{16}{10}} = 0.02511886 \]ここまでの話をまとめると、1 番目のリードの最初の塩基 G のクオリティスコアは 16 であり、そこから計算されるシークエンシングエラーは 0.025 である。つまり、この塩基 G は、0.025 の確率でシークエンシングエラーとなっている場合がある。
クオリティスコアとシークエンシングエラーの計算式に着目すれば、両者の間に次のような関係が見られる。クオリティスコアが 10, 20, 30, 40 のときのシークエンシングエラーを大まかに把握しておくとよい。
- クオリティスコアが 10 ならば、シークエンシングエラーが生じる確率は 0.1 であるから、読み取られた塩基の信頼度は 90.0% である。
- クオリティスコアが 20 ならば、シークエンシングエラーが生じる確率は 0.01 であるから、読み取られた塩基の信頼度は 99.0% である。
- クオリティスコアが 30 ならば、シークエンシングエラーが生じる確率は 0.001 であるから、読み取られた塩基の信頼度は 99.9% である。
RNA-Seq の場合、シークエンサーで読み取られたリードがリファレンスのどこにマッピングされているのかさえ明らかにすれば発現量を見積もることができるので、シークエンサーが出力するクオリティスコアをそのまま信用して解析を進めても差し支えない。一方で、GWAS の場合、解析個体のゲノム塩基を明らかにする必要があるため、リードのマッピング位置だけでなく、リード上の個々の塩基がどのぐらい信頼できるのかも重要になってくる。このため、GWAS の場合は、シークエンサーが出力するクオリティスコアをそのまま信用して解析するのではなく、それらのクオリティスコアを補正したものを使うのが一般的である。
現在では、FASTQ に使われているクオリティスコアはほとんど上述のように Q + 33 のものを利用している。古いシークエンサーが出力したデータを解析するときに、クオリティスコアが Q + 64 となっているかどうかを一度確認しておく必要がある。
SSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSS..................................................... ..........................XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX...................... ...............................IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII...................... .................................JJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJ...................... LLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLLL.................................................... !"#$%&'()*+,-./0123456789:;<=>?@ABCDEFGHIJKLMNOPQRSTUVWXYZ[\]^_`abcdefghijklmnopqrstuvwxyz{|}~ | | | | | | 33 59 64 73 104 126 0........................26...31.......40 -5....0........9.............................40 0........9.............................40 3.....9.............................40 0.2......................26...31........41 S - Sanger Phred + 33 X - Solexa Solexa + 64 I - Illumina 1.3+ Phred + 64 J - Illumina 1.5+ Phred + 64 L - Illumina 1.8+ Phred + 33
References
- FASTQ format. Wikipedia