Top of Technical Information

マイクロアレイより確実に面倒な、RNA-Seqのデータ解析について

RNA-Seqでトランスクリプトームの測定を行うと、結果としてcountまたはFPKMのテーブルが得られます。これをエクセルで見てみると、0あるいは値が欠けているセルがたくさんあることがわかります。これらは「マッピングされたリードが無い」ことを表していますが、これが即ちその遺伝子が「発現していない」ことを表しているとは言えないのが、RNA-Seqのデータ解析の難しさの一つです。

その理由は、繰返し実験があると明らかです。例えば、2回の繰返しがあって、1回目のサンプルでは値があるのに、2回目では0ということがあります。一方のサンプルではcountの値がシグナル領域にあっても、もう片方では0となることもあります。これは、countが0でも「発現していない」とは言い切れず、「本当は発現しているにも関わらずシーケンサーにキャプチャーされなかった」可能性があることを示唆しています。これはRNA-Seqがデジタルな手法であることに起因する現象で、アナログなマイクロアレイでは起こらない問題です。そしてこのことは、RNA-Seqでは発現量の少ないRNAをすべて観測することは事実上不可能だということを示しています。(理論上は、数百回以上の繰り返し測定することで可能です。)

このような特徴を踏まえると、あるグループで平均値をとる場合、測定値があるものと0が混在していれば、0を含めた平均値を使うと、発現量がそのグループで不当に下がって見えることになります。0の影響を除いてから発現差解析をする必要があります。

ところで、一律で測定値に1や2を加えて0をなくす解析手法がありますが、これはまったく問題の解決にならず、問題を見えにくくするだけですのでやめたほうが良いでしょう。

Subio Platformでは、このような問題に配慮して解析をするために、RNA-Seqのデータ解析では、下記のNormalization Scenarioを提案します。プリセットのNormalize Scenarioで、「Profiling (GDC RNA-Seq counts)」を選択すると、自動的にNormalize Blockは配置されます。ただし、各ブロックのオプション設定は、データに合わせて調整する必要があります。どのようなRNA-Seqデータにも適合するわけではありませんのでご注意ください。

  1. Log Transformation
  2. Global Normalization
  3. Low Signal Cutoff
  4. Centering あるいは Ratio to Control Samples

ノーマライズシナリオの解説

最初にLog Transformationブロックを持ってくることで、countあるいはFPKMが0だった場合、これらは欠損値になります。こうしておくことで、その後グループ間で比較するときに、欠損値は自動的に排除されて計算されることになります。

次のGlobal Normalizationブロックについては、CountもFPKMも、すでにノーマライズされている場合がほとんどですので、無くてもいい場合は外してください。しかしヒストグラムを見て、正規化が必要と判断される場合は適用してください。適用する場合は、ヒストグラムを見ながらシグナル領域の山のピークが揃うように、基準値を75 th ~ 90th percentileあたりで調整してください。シグナル領域にある遺伝子が少ないほど、揃える基準は高いところになります。インプットのRNAがごく微量の場合は、99th percentileという極端な値になることもあります。

Low Signal Cutoffブロックは、ノイズ領域で発現変動が不当に大きく見えるのを抑制する目的で使います。ノイズ領域とシグナル領域の境界あたりからその半分の値までの範囲で閾値を調整するとよいでしょう。

最後のブロックは、比にするものです。コントロールサンプルが無い場合はCenteringブロックを選び、コントロールサンプルがある場合はRatio to Control Samplesブロックを選んでください。

発現差のある遺伝子の集め方の注意点

ところで、このような前処理で生成されたProcessed Signalは、分母または分子の値が欠損値の場合には計算不能のため、欠損値になります。そうすると、fold値やT検定のP値も計算不可能となります。しかし、たとえばグループAのサンプルではすべてcountが0だったのに対し、グループBではある程度信頼できるcountの値があるという場合は、グループBで発現が上昇していると言えます。このような遺伝子群は、フィルターツールを使って、グループAでは100%のサンプルでCh1RawSignalが0の遺伝子群と、グループBでは50%以上のサンプルでCh1RawSignalがたとえば10より大きい遺伝子群を集め、この二つをVenn Diagramツールを使って掛け合わせて抽出してください。

RNA-Seq のデータ解析手順

データインポートから解析まで一通り学びたい方は、オンライントレーニングをお申し込みください。

多変量解析

発現差解析は上記の前処理でいいのですが、クラスタリングやPCAには欠損値が邪魔になります。そこで、上記の前処理で発現差解析を行った後、欠損値の影響を避けたい場合は、Fill Missing Valuesブロックを追加してProcessed Signalを計算しなおしてから、クラスタリングやPCA Detailsを実行するほうがいいでしょう。おすすめの使い方は次の2通りです。

最後がCenteringになっている場合

Fill Missing Valuesブロックを 最後に 持ってきます。欠損値と入れ替える値は、-1(log2のスケールにおいて1/2に減少を表す)か-2にするといいでしょう。

最後がRatio to Control Samples になっている場合

F ill Missing ValuesブロックをRatio to Controls Samplesブロックの 前に 持ってきます。欠損値と入れ替える値は、Low Signal Cutoffで用いた閾値か、それよし少し低い値でいいでしょう。

ノーマライズ: Fill Missing Values ブロックの使い方」も合わせてご覧ください。