ノイズはただのゴミか?

同じサンプルをマイクロアレイで2回測定して比べることを想像してください。 二つの測定値で散布図を描くと、Fig4Fig5のようになりました。 あなたなら、どちらが質の高いデータだと思いますか?

Fig4のほうがばらついているように見えます。 相関係数を計算すると、Fig4では0.919、Fig5では0.988となりました。 多くの人がFig5のほうが良いデータだと思うのではないでしょうか。 実はこの二つは、もともと同じ生データをMAS5(Fig4)とRMA(Fig5)というアルゴリズムを使ってシグナル値を計算したものです。 二つのアルゴリズムの特徴についてはMAS5とRMAについての解説をご覧いただくとして、ここでは「再現性」について考えます。

たしかに高い相関係数は再現性が高いことを表しますが、「相関係数が高いから再現性が高く質のいいデータである」とは言えないのです。

あるサンプルにおいてすべての遺伝子が発現しているとは生物学的に考えられません。 また、すべての測定システムはきちんと測定できるダイナミックレンジがあり、その限界を超えた領域では正しく測定できません。 マイクロアレイで測定したシグナル値の中には、発現していない遺伝子と、測定できないほど発現量が少ない遺伝子のものも必ず含まれています。 このような測定値は、値がとても低い領域で、繰返しサンプルにおいて再現性がないように見えることが予想されます。

Fig4では予想どおりノイズ領域が見られるのですが、Fig5では見られません。 しかしノイズが無いということは実験生物学の立場からあり得ないことなので、Fig5ではノイズが隠れていると考えられます。

データを解釈する立場からは、ノイズが隠れていることはとても厄介です。 なぜならノイズがノイズとして見えるから、シグナルとノイズの境界がわかり、測定値が信頼できるシグナル領域の遺伝子群を簡単に抽出できます(Fig4の黒)。 しかし、ノイズが隠れているFig5では、これが簡単にできないのです。 ただし、ヒストグラムを注意深く見ると、左側(値の低い領域)に高い山が見えます。 Fig4で抽出したシグナル領域の遺伝子群(黒)は、Fig5においてもよく似た分布をしていますから、やはりこの山がノイズ領域と考えられます。

「高い相関係数は再現性が高いことを表す」のは、シグナル領域において成立することです。 ノイズ領域では通用しません。 従って、相関係数はシグナル領域の遺伝子だけで計算するべきものであって、すべての測定値を使うのは適切ではありません。 これを混同するから、相関係数を高く見せるためにノイズを隠すようなアルゴリズムが作られます。 また、ダイナミックレンジとはシグナル領域を指すのであって、すべての値の分布域のことではありません。 これを混同するから、むやみにノイズ領域の広いデータが作られます。 ノイズがノイズとして見えることには価値があるのです。

それでは以上をふまえて、マイクロアレイやRNA-Seqのダイナミックレンジを見ていきましょう。

Paper fig04 thumb

GSE11670もHG-U133 Plus2.0アレイを使った測定データです。 非常に質の高いテクニカルレプリケートの例です。 散布図の左下では大きなばらつきが見られますので、シグナル値が100より低い領域では測定値の信頼性が低いと言えます。 ヒストグラムを見ると、MAS5のシグナル値の分布は一山型に見えますが、実は本質的に二山形です。 シグナル領域とノイズ領域の重なる領域が広いので、ここにピークが見えているのです。 つまり、シグナル値が⒑から100の間の測定値は、信頼性は低いですが、完全にノイズではありません。

Paper fig05 thumb

Fig4と同じGSE11670の生データをRMAで再解析したものです。 散布図の左下に大きなばらつきが見えません。 しかし、だからと言ってこのデータにノイズがなく、すべての測定値がシグナルだというわけではありません。 ヒストグラムの3および4行目をご覧ください。 非常に高いピークが左端に見られます。 つまり、ノイズ領域にあった値を狭い領域に押し込めたにすぎないのです。 従って、シグナル値が10から100の領域は、MAS5の時よりノイズ成分が濃くなっているため信頼性が劣ります。

マイクロアレイのダイナミックレンジの比較

ダイナミックレンジを4桁とか5桁とかいう値を以て比較するのは、100キロメートルは80マイルより長いと主張するようなものです。 マイクロアレイの測定値は、それぞれの機械と数値化アルゴリズムに依存し、同じ単位の数値ではありません。 従って、測定値の分布の幅を比較してもしかたがありません。 ダイナミックレンジの比較は、シグナル領域にある遺伝子の数でするのです。 言い換えれば、何個の遺伝子がちゃんと測定できるかを問うのです。

Affymetrix 3' IVT GeneChip

2002年から2004年にかけて大量に使われたマイクロアレイで、代名詞的な製品でした。 Fig4およびFig5のデータはHuman Genome U133 Plus 2.0 Arrayのもので、このタイプのマイクロアレイの典型なデータでは、シグナル領域(黒い領域)にだいたい25,000個くらいのプローブセットがあって、遺伝子に対して複数のプローブセットが設計されていることも多いので、それを差し引いておよそ12,000個の遺伝子を測定できます。

理論的にはノイズ成分とシグナル成分を分けて考えることができますが、実際のデータではきれいに分かれていないので、境界は線ではなく帯またはスペクトルのようなものになります。 ノイズ領域とシグナル領域の重なりが大きければ、ヒストグラムで見ると境界領域は山型になります。 MAS5 (Fig4) では、黒いシグナル領域のすぐ左に山型が見えます。 ここがノイズ成分とシグナル成分が混じっている境界領域です。 そして、その左側の裾野はノイズだけの領域です。 信頼性は落ちますが、境界領域のデータも使えなくはありません。 しかし、RMAのデータ (Fig5) ではそのような分離はできません。

GSE11670のデータをSSA形式でダウンロードして、Subio Platformで詳しくご覧ください。

Illumina BeadChip

GeneChipより後発で製品化されたもので、非常に高い再現性を謳っていました。 HumanHT-12 V4.0 expression beadchipのデータを散布図(Fig6)を見るとばらつきが無いように見えますが、これは上述のRMAの効果にすぎません。 シグナル領域にあるプローブはおよそ15,000個で、遺伝子数にして11,000個くらいです。 Affymetrixの3' IVT GeneChipと比べると、ダイナミックレンジは狭くなっています。 また、RMAによりデータの解釈が難しくなるのですが、それ以外の数値化手法を選択できません。

GSE25315のデータをSSA形式でダウンロードして、Subio Platformで詳しくご覧ください。

Affymetrix Gene ST Array

Paper fig06 thumb

GSE25315は、Illumina BeadChipシステムで測定されたデータです。 RMAで数値化されているため、散布図で見ると再現性が非常に高いように見えます。 しかし、ヒストグラムで見ると、ノイズ領域に大量の遺伝子があることが一目でわかります。 散布図だけでなく、ヒストグラムと合わせてデータを見ることが大事です。

Affymetrixによる3' IVT GeneChipの後発なのですが、このデータの特徴は、IlluminaのBeadChipと同様に数値化ソフトウェアによってバイアスとノイズを隠していることと、シグナル成分とノイズ成分の境目が非常に分かりにくいということです。 Fig7でいうと、シグナル値がだいたい6以上の領域では、発現量が下がるにしたがってばらつきが大きくなっています。 これは自然な形です。 しかしそこから下は逆に収束していきます。 これはソフトウェアによって作られた不自然な形です。 このような意味でシグナル領域をとると、プローブセットで12,000~15,000個、遺伝子数で9,000~11,000個が測定可能のようです。

ちなみに、ノイズ領域をもっときちんと知りたいと思うなら、その細胞で発現していない遺伝子を検索してみてください。 そのような遺伝子が観測される領域がノイズ領域です。 ただし、ノイズ成分とシグナル成分が混在する領域が広いので、実際の解釈はとても難しいです。

Affymetrixの製品は、Plus 2.0シリーズ以降、新しい製品ほどダイナミックレンジが狭くなっていたのですが、最新のClariom S Arrayに至っては大幅に縮小しました。 かつてのテクノロジーリーダーの面影は、完全に見えなくなりました。

GSE22288のデータをSSA形式でダウンロードして、Subio Platformで詳しくご覧ください。

Agilent Whole Genome 4x44k Microarray

Paper fig07 thumb

GSE22288は、Affymetrixが3' IVT Array(HG-U133 Plus2.0など)の後継として出したGene ST Arrayによる測定データです。 3' IVT Arrayでは、エクソンよりユニークな配列の多い3' UTRにプローブが設計されていたのに対し、Gene ST Arrayのプローブはエクソンに設計されています。 新しい技術が古い技術よりも良いとはかぎりません。 歴史的に見ても成熟した技術のほうが優れている例は少なくありません。 また、3' IVT GeneChipのデータはインターネットから大量に入手できるというメリットも考慮に値します。

Agilentの、これより前の世代のマイクロアレイは質の良いものではなかったのですが、この製品以降は他の追随を許さない圧倒的なダイナミックレンジの広さを誇っており、その地位は2004年の発売以来一度も揺らいでいません。 Fig8のように、シグナル領域にある数はプローブ数でおよそ34,000個、遺伝子数で16,000個以上になります。 ダイナミックレンジが広いということは、中~低発現領域の遺伝子の発現量を測定できるようになったということですが、生物学的にはさらに大きな意味があります。

上述のとおりノイズ成分とは、発現していない遺伝子と、発現が低すぎて正確に測定できない遺伝子が不可分に混ざり合ったものです。 発現していない遺伝子を分離するには、測定感度が上げて後者の成分を減らすしかありません。 つまり、4,000~5,000個の低発現領域の遺伝子を測定できるようになるだけでなく、発現していないと考えられる遺伝子もわかるようになるということです。 発現のON/OFFスイッチの切り替えを見るためには高い感度が必要なのです。

GSE36082のデータをSSA形式でダウンロードして、Subio Platformで詳しくご覧ください。

Paper fig08 thumb

GSE36082は、マイクロアレイの歴史において大きな進歩を遂げたAgilent Whole Genome 4x44kによって測定された実験データです。 この後継であるSurePrintマイクロアレイも高い品質を維持しています。 ヒストグラムを見ると、左端に急峻なピークが見えます。 これはネガティブコントロースのプローブと、それと同等のシグナル値を示すプローブによって形成されるものです。 優れている理由は、ノイズ領域とシグナル領域の重複が少なく、境目がくっきりしていることです。 これにより、低くても発現している遺伝子と、発現していない遺伝子を区別して抽出できる可能性があるのです。

RNA-Seqのダイナミックレンジ

マイクロアレイでのノイズ成分は、発現していない遺伝子と、発現が低すぎて正確に測定できない遺伝子でしたが、RNA-Seqでは少なくとも前者の成分は原理的に無いと考えられますが、後者のノイズ成分は免れないと考えられます。 ただしマイクロアレイのノイズとは異なる特徴がありますので、詳しく見ていきましょう。

Countデータ

Fig9はHiSeq2000によるRNA-Seqのデータで、値はcountです。 値が20より下で急にばらつきが大きくなっており、ノイズ領域があることがわかります。 シグナル領域にはおよそ12,000個の遺伝子があります。 つまりダイナミックレンジはAffymetrixの3' IVT GeneChipとほとんど変わらないということです。

しかし、さらに詳しくRNA-Seqのノイズ領域を見ていると、データの解釈を難しくするある特徴に気が付きます。 それは、片方のサンプルでは値があって、もう片方のサンプルでは値が0という遺伝子があるということです。 PCRによる増幅がかかったりかからなかったり、シーケンサーが読むのはごく一部のRNAなので読まれたり読まれなかったり、そして読んだ配列がマッピングしやすいところだったりそうでなかったり、などの偶然によって値があったりなかったりするのです。 これはつまり、値が0だからとって発現していないとは断言できないことを表します。

RNA-Seqでもまた、ノイズ成分は発現していない遺伝子と、発現が低すぎて正確に測定できない遺伝子が不可分に混ざり合ったものとなります。 これを克服することは原理的には可能でも、コスト面から実用化は不可能でしょう。 RNA-Seqが感度においてアジレントのマイクロアレイに匹敵することは、これからもずっとないだろうと予想しています。

補足ですが、繰返しサンプルの平均値を計算するときに、0を入れるかどうかを検討する必要がありますので、データの解釈もマイクロアレイより難しくなります。

GSE49110のデータをSSA形式でダウンロードして、Subio Platformで詳しくご覧ください。

FPKMデータ

Paper fig09 thumb

GSE49110は、Illumina HiSeq2000を使ったRNA-Seqの実験データです。 RNA-Seqは原理的に発現していない遺伝子によるノイズを含まないと考えられますが、散布図の左下の領域にはノイズ領域があります。 シグナル領域にある遺伝子はおよそ12,000個です。 つまり、中~低発現の遺伝子の発現量は測定できません。 これでHiSeq2000のデータです。 Genomic Analyzerの時代のRNA-Seqのデータがどれほど品質の低いものだったか想像に難くないでしょう。

Fig10もHiSeq2000によるRNA-Seqのデータですが、値はFPKMです。 Countのデータでははっきりと見えていたノイズ領域が、なんとなく見えなくなっています。 FPKMは、100万リードあたり・エクソン1 kbpあたりのフラグメントの数ということです。 「100万リードあたり」はサンプル間のリード総数の差を均すノーマライズに相当します。 「エクソン1 kbpあたり」にというのは、大きい遺伝子ほどたくさんマッピングされるリードがあって、小さい遺伝子はマッピングされるリードが少ないことが予想されるので、遺伝子の大きさによらず発現量を見積もるために行うものです。

これだけ聞くと、なるほどそういうものかと思うでしょう。 しかしよく考えてみてください。 長さで均すということは、さきほど見たcountの点が、大きい遺伝子は左下方向に、小さな遺伝子は右上方向に移動するということになります。 これが、ノイズ領域がはっきりと見えなくなる理由です。 FPKMにはノイズがないのではなく、見えにくくなっていると知るべきです。 エクソンの長さの平均値は1kbよりも長いので、多くの遺伝子は左下方向に移動し、値の範囲は2桁ほど下に延びますが、これは長さで均したことの副作用にすぎません。 これを以て「RNA-Seqはダイナミックレンジが広い」とか、「低発現領域の感度が高い」と言うのはまったくナンセンスです。

もし、同一サンプルにおいてどちらの遺伝子の発現が多いのかを知りたいのであれば、長さで均すことは必要でしょう。 しかし、ある遺伝子が別のサンプルで増えるか減るかを知りたい場合は、長さで均す必要はありません。 シグナル領域とノイズ領域を区別するのが簡単なのは、countのデータです。 FPKMを使って解析をする場合でも、まずはcountのデータを使ってノイズ領域の遺伝子を除く作業を行ったほうが簡単です。 また、これも論文で見かけるのですが、たとえば「FPKMが10以上は信頼できる」というロジックには注意が必要です。 なぜなら、countの値がノイズ領域にあっても、非常に小さい遺伝子はFPKMの値が大きくなるからです。

以上の理由でcountよりもシグナル領域は明確でないのですが、FPKMでもシグナル領域にある遺伝子の数はおよそ12,000個です。

GSE53567のデータをSSA形式でダウンロードして、Subio Platformで詳しくご覧ください。

Paper fig10 thumb

GSE53567も、Illumina HiSeq2000を使ったRNA-Seqの実験データです。 ただし、数値の出し方が異なっています。 Fig9はcountの値ですが、こちらはFPKMの値です。 FPKMの場合、遺伝子の大きさを均す処理が入っているため、大きな遺伝子は左下方向に、小さな遺伝子は右上方向に動きます。 この結果、シグナル領域とノイズ領域がcountのとき(Fig9)よりも曖昧になるのですが、当然ながらノイズがなくなるわけではありません。

質の高い発現データとは何か

これまでさまざまな種類のマイクロアレイ、そしてRNA-Seqのデータのダイナミックレンジを見てきました。 もしあなたが、発現の高い上位10,000個~12,000個の遺伝子を解析したいのであれば、解釈のしやすさには違いがあるものの、おおむねどのプラットフォームを選んでもかまわないです。 しかし、それよりも低い発現量の4,000~5,000個の遺伝子を見たい場合や、発現していない遺伝子は何かを知りたい場合は、アジレントのマイクロアレイ以外の選択肢はありません。

2004年後半にアジレントの革命的なマイクロアレイが出ていたのですが、時代はRNA-Seqをもてはやすようになり、10年以上も無駄に停滞しているという状況は驚きではないでしょうか。 しかし、 これがトランスクリプトミクスの現状です。 オミクスデータの解析についてオープンディスカッションが十分されてこなかった状況が招いた膨大な損失です。 今後このような無駄を避けるために、私たちは次の4つの観点からデータの品質を検討することをおすすめします。

  • ノイズがノイズとして明確に見えている。
  • シグナル領域にあって、きちんと測定できる遺伝子の数が多い。
  • シグナル領域とノイズ領域の重なりが少ない。
  • 発現差の検出感度が高い。

最後の点についてはこれまで述べていませんでしたが、再現性をよくするためにサンプル間のばらつきを抑えることと、発現差を感度良く検出することは、互いにトレードオフの関係にあります。 再現性がとても良いように見える測定装置あるいはアルゴリズムを使った場合は、小さな発現差を検出できない可能性が高いでしょう。 逆に、研究対象の発現差が小さいと予想される場合は、できるだけ高品質な生データをとることと、生データをあまり加工しないアルゴリズムを選択することが重要です。

最後に実験生物学の観点から、一般に正確性と精度はトレードオフの関係にあるので、このバランスを考慮することが大事です(Fig12)。 きれい過ぎるデータは、かえって注意が必要かもしれません。 バイアスのかかったデータは生物学的発見の妨げとなるどころか、誤った結論を導き出しかねません。 詳しくは実験誤差に関する理論の解説をご覧ください。

Paper fig12 thumb

正確性と精度はトレードオフの関係にあるため、同時に高めることは一般的にはできません。 正確性の高いデータであれば、たとえ精度が低くても、繰り返しの測定値の平均値から真の値をうまく推定できると考えられます。 しかし、精度が異常に高く見えるようなデータの場合、むしろバイアスがかかっていないか注意が必要です。 バイアスがかかっているとしたら、その原因はたいてい数値化と正規化のアルゴリズムです。 オミクスデータを見るときは、このことを念頭においたほうがいいでしょう。

私たちの提案するソリューション、Subio Platform

オミクスデータをどう解析すればいいのか、知っている人は一人もいません。 解析の専門家などいないのです。 おそらく決定的なアイデアを持ち込んでくるのは、統計学者やバイオインフォマティシャンではなく、実験と現象に精通している実験生物学者でしょう。 彼らが積極的に加わって、すべての参加者が等しくアマチュアとしてアイデアを出し合って検討しなければ、オミクスを活用した新しい生命科学はこれ以上先に進むことができないというところに来ています。

私たちがこの問題にどのように向き合い、なぜSubio Platformを提案しているかご覧ください。

Back to Top