前回の DROMPA3: その2 parse2wig - Palmsonntagmorgen で登場した評価指標である Library complexity (PCR bias) の解説です。
PCR biasとは
クロマチン免疫沈降法（ChIP）で得られたリードをゲノムにマップすると、以下のようなマップリード分布が得られるでしょう。
この図ではゲノム上（下部の黒い線）にリードがマップされている様子を示しています。ゲノム上にはリピート領域やCopy number variation (CNV)などがあるためリード分布は完全に一様にはなりませんが、それでもある程度ランダムにゲノム上に張り付きます。 リードがたくさんマップされている（濃縮している、enrichしている、とも言います）領域は、目的のタンパクが結合する領域かもしれません。
それでは、リードが以下のような分布をしていたらどうでしょうか。
リード数は最初の例と同じです。このような分布をしていた場合、リードが張り付いている領域をタンパク結合領域だと考えてもよいのでしょうか。
結論から言えば、この例のような場合、リードがたくさんマップされていても濃縮領域であるとは考えません。
ヒトのChIP-seqサンプルのリード数は通常多くとも3000万リード程度ですが、ヒトゲノムは約30億塩基対ありますので、マップリードが一様に分布しているとすると、100塩基対あたり1リードがマップされる計算になります。順鎖、逆鎖に5割ずつマップされていることを考慮すれば、鎖あたりでは200塩基対あたり1リードになります。
その程度の確率の中で、上図のように全く同じポジションに多量のリードが張り付くというのは偶然の産物とは考えられませんので、これらはPCR増幅の産物であると考えられます。
このような偏りのある（PCR biasの大きい）データをそのまま解析に用いると偽陽性（false positives）の原因になりますので、このような同一箇所にマップされたリードは1つを残して全て以下のように除去されます。
PCR biasをフィルタした結果、この例ではこの領域に4リードしか残りませんでした。PCR biasとしてフィルタされるリードのことを「冗長なリード(redundant reads)」、フィルタ後に残るリードを「非冗長なリード(nonredundant reads)」と呼びます。このように、見かけのマップリード数が多くても、実際に解析に用いられる非冗長なリード数が少ない場合、解析に必要なリード数を満たさない可能性がありますので注意が必要です。
なお、このPCR biasのフィルタリングはChIP-seqだけでなく、オープンクロマチン解析やHi-C解析などでも利用されています。RNA-seqではリードが同一のポジションにマップされることはあり得るため、通常用いられません。
Library complexityとは
読まれたサンプル内に存在するPCR biasの割合を定量的に評価するため、ENCODE projectで"Library complexity"という指標が提案されました*1。Library complexityはnonredundant fraction (NRF)、すなわちマップリード数に対する非冗長なリードの割合によって計測されます。
しかしこれだけだとまだ問題があります。同じサンプルでもリード数が増えるほど偶然同じ場所にマップされるリードの数が増えるため、同一サンプルでもマップリード数によってNRFの値が変わってしまいます。そこでENCODEでは1000万マップリードあたりの非冗長なリードの数をスコアに用いることを提唱しています。
parse2wigのスタッツファイルで出力されるLibrary complexityはこの値を示します。 ヒトのサンプルでは、この値が0.8以上あることが望ましいとされています。
PCR bias の要因
私の経験上、PCR biasが大きい、すなわちLibrary complexityが低くなる原因の多くは、PCR増幅をかける前の初期DNA量が少ないことによります。
たとえば抗体のS/N比が悪く、タンパク結合DNAがほとんど回収できていない、あるいは細胞数が少なすぎるような場合に、PCR増幅を多くかけることによってシーケンシングのためのDNA量を確保しようとしますが、そのようなサンプルは得てしてLibrary complexityが低くなります。そうすると、解析に使えるリード数は結局足りないということになりがちです。更に、Library complexityの非常に低いサンプルはリードのGC含量に強い偏りがあるなど、他の品質評価基準を満たさないことも多いため、リードを読み足すよりはサンプルを作り直す方が望ましいです。
一方、DNAの断片化に超音波破砕ではなくDNaseを用いた場合、制限酵素サイトによってDNA断片の端が切り揃うのでLibrary complexityが低くなるのでは、という疑問が湧く方もいるかもしれません。私の解析したサンプルの中では、体感ではほとんど違いがないか、あるとしても2~3%程度です。ので、断片化手法の違いはあまり気にする必要はありません。
Library complexityについての考察
私のこれまでの経験では、クオリティが良いとされるサンプルは確かにだいたいLibrary complexityの値が0.8以上ありました。が、0.8を切っているサンプルが直ちに使えないという訳ではありません。0.6程度だったとしても、フィルタ後に残った非冗長なリードの絶対数が十分あれば（最低1000万くらい）、問題なく解析できることもあります。
また、RNA pol2のように非常に抗体のS/N比が高いサンプルでは、リードの濃縮領域に多量にマップされた「本物の」リードがある程度フィルタされてしまうため、Library complexityが低くなります。真のマップリードの一部を除去してしまうことの是非には議論がありますが、ピーク抽出の感度への影響は「ほとんど気にしなくてよいレベル」であると報告されています*2。
また、酵母のような小さいゲノムの場合、最初に計算した「塩基あたりにマップされるリードの数」が1以上となります。この場合、1カ所につき1リードのみ許すという基準は厳しすぎるものになりますので、基準を多少緩める必要があります。
まとめ
まとめると、Library complexityのスコアと非冗長なリードの絶対数を確認するのは重要である一方、その値が問題なければ、実際の解析においてPCR biasを除去する必要性自体はそこまでない（かもしれない）ということになります。公開されているChIP-seq解析ツールの大半は（DROMPAも含め）PCR biasをデフォルトで除去しますので、敢えてオフにする必要は特にありませんが、S. cerevisiaeのrDNA 領域のようなhighly repetitiveな領域を観察したい場合には、このフィルタリングをオフにする必要があります。
