Intensity-based Z-score (シグナル値を考慮した Z-score )

マイクロアレイにおいて、シグナル値の低いものは、信頼性(データの再現性)がよくない場合があります。これは、MAプロットや散布図の形状からも推測されます。(同一条件で比較したサンプルにおいても、シグナル値の低い部分に変動が見られがち。)

例えば、シグナル値が、1000から10000になったら、ratio=10ですが、10から100になっても同じ ratio=10 です。ratio だけで、発現変動遺伝子を判定すると、特にシグナル値の低い部分に、ratio が高い遺伝子が見つかりがちですが、再現性を考えた場合、少し不安です。

一方、 Z-score で判定した場合も、MA プロットから分かるように、シグナル値の高い部分は、Z-score が高くなりにくい傾向があります。

この2点の問題を解消するために、  Intensity-based Z-score (Intensity-dependent Z-score)* というものがあります。(以前から提案されているもので、珍しい手法ではありません。)

Intensity-based Z-score

Intensity-based Z-score では、文字通り、Intensity (=シグナル値の強さ) を考慮した Z-score です。マイクロアレイデータをシグナル値の大きさごとに一定の区間に分割して、その区間ごとに logFC の平均値と標準偏差を求めます。そして、それぞれの区間ごとに、Z-score を算出します。

そうすることで、シグナル値の低い部分は、大きい標準偏差 (SD) を用い、一方、シグナル値の高い部分は、小さい SD を用いて Z-score を判定することになります。したがって、シグナル値の低い部分は、通常の Z-score より判定条件が厳しくなり、シグナル値の高い部分は、判定条件が緩くなります。

MAプロットで確認すると、違いがよくわかります。ここでは、便宜的に、シグナル値の大きさによって8区間に分割して、Z-score を算出しています(図下)。(8つという区間に必然性はありません。もう少し細く区切ってもよいです。)また、Intensity-based Z-score を使用した図は、事前に使用した2サンプルのフラグが、共にA(absent: 未検出)のプローブを除外しています。

Z-score で色付けした MA プロット。
Z-score で色付けした MA プロット。
Intensity-based Z-score で色付けしたMAプロット。(フラグでカットしている。)
Intensity-based Z-score で色付けしたMAプロット。(フラグでカットしている。)

図形的には、変動している遺伝子を外側から数パーセント取るという形になります。単純に ratio で判定するよりも慎重な方法といえるでしょう。しかしながら、これも万能な手法というわけではありません。可能であれば、繰り返しサンプル (replicates) が、3サンプル以上 (n=3以上) あって、検定などの手法が使える方が望ましいことには変わりません。

参考

* Quackenbush. Microarray data normalization and transformation. Nat Genet (2002) vol. 32 Suppl pp. 496-501.

 

Z-score の利点と欠点

MAプロットで確認すると、 ratio と Z-score の関係が、分かりやすいと思います。解析例1のデータ (sample2/control2) においては、ratio=2 のラインより、2SDのラインが外側にあります。この場合は、 Z-score のほうが厳しい判定と言えます。(追記:この例の場合、2SD=約1.344 (logFCで)、ratio > 2.54 が Z-score > 2 のラインです。)

MAプロットに示した Z-score の関係。
MAプロットに示した Z-score の関係。

また、 logFC を Z-score 化して、発現変動遺伝子の判定を行う場合、利点と欠点があります。

Z-score を判定に用いた場合の利点

単純に ratio で判定する場合に比べて、下記のような利点があります。

  • 判定の基準で迷わない。
  • 発現変動が大きいデータの場合、カットオフに使える。
  • 発現変動が小さいデータでも、発現変動遺伝子を判定できる。

Z-score > 2 または、 Z-score < -2 ということに、統計学的な意味があるので、判定の基準をいくらぐらいにしようかということで迷う必要はありません。(2 以下にあまり意味がありません。)

また、発現変動が大きいデータの場合、ratio > 2, ratio < 0.5 などで判定すると、4000個、5000個の遺伝子が、発現変動ありとなることもあります。そのような大きいデータの場合、 Z-score で判定することで、発現変動遺伝子の数を1000個から2000個程度に減らせます。(発現変動遺伝子が多すぎても困ることがあります。例えば、DAVID で1度で解析できる遺伝子の数は、3000個です。)

逆に、サンプルによっては、たかだか1.5倍程度しか変動しないような、小さい発現変動のデータもあります。そのような場合、Z-score だけで判定すれば、(ratioの大きさに関係なく)変動している遺伝子の上位5%程度を取ることになるので、必ず、変動遺伝子を得られます。

Z-score を判定に用いた場合の欠点

利点の裏返しです。

  • 判定が厳しすぎるときがある。
  • シグナル値が高い部分で、|Z-score| > 2 に比較的なりにくい。

例えば、研究対象としている遺伝子があり、 ratio だと 2 倍を超えているのに、 Z-score だと、 1.7 というようなケースです。Z-score  が2以下だからといって、「変動していない」という証明になるわけではありませんが、印象が良くないかもしれません。

また、散布図やMAプロットの先が細いということからも見えることですが、シグナル値の高い部分は、変動が小さく見えがちです。シグナルの高い部分にも、シグナルの低い部分と同程度の標準偏差を求めるならば、高い部分は、|Z-score| > 2 になりにくいといえます。

しかしながら、そもそも、変動が大きすぎたり、変動が小さい場合は、まず、サンプル数(n数)を増やすことを考えるべきです。n=3以上にして、limma, SAM, t-test などの検定を行ったほうが、p-value も算出できますし、より再現性の高い遺伝子を変動ありと判定できます。

 

Z-score (発現変動遺伝子を判定するもう1つの方法)

ratio (logFC) 以外に発現変動遺伝子を判定する方法として用いられるものに、 Z-score があります。あまり聞きなれない用語かもしれませんが、偏差値というとどうでしょうか?

ある値が、その群の平均値から、標準偏差 (SD) の何個ぶん離れているかを求めたものが、 Z-score です。Z-score の考え方自体は、特殊なものではありません。Z-検定 (Z-test) という使われ方もあります。(Z-score 化されるのは、logFC に限りません。例えば、 cBioPortal では、シグナル値が Z-score 化されています。)

ここで用いるのは、 logFC を Z-score 化したものです。ある logFC  が、 logFC の平均値から、標準偏差の何個ぶん離れているか計算します。

ratio (logFC) のヒストグラム

ratio (logFC) のヒストグラムを示します。データは、解析例1のものを使っています。

ratio (logFC) のヒストグラム。
ratio (logFC) のヒストグラム。

ほとんどの遺伝子の logFC が0(つまりratio=1)、変動していないことが分かるかと思います。ratio で判定するなら、 ratio > 2 または、 ratio < 0.5 を満たす部分に含まれる遺伝子を発現変動遺伝子と判定していることになります。(ヒストグラムの右側と左側)。

標準偏差と Z-score

次に、同じヒストグラムに標準偏差 (SD) と Z-score の関係を示してみます。

標準偏差とZ-scoreを示したヒストグラム。
標準偏差とZ-scoreを示したヒストグラム。

解析例1のデータ (sample2/control2) の場合、計算してみると、logFC の平均値は、ほぼ0になります。SDは、約 0.672 です。また、分布関数の考え方から、平均値から標準偏差2個分離れていれば、稀な値 (p-value < 0.05) ということになります。

したがって、 Z-score で判定する場合は、下記のようになります。

  • 増加した遺伝子: Z-score > 2
  • 減少した遺伝子: Z-score < -2

絶対値の記号を用いれば、両方を合わせて、  |Z-score| > 2 とも書けます。

追記:2SDだと、2 x 0.672 = 1.344 なので、 logFC > 1.344 であれば、 Z-score > 2 です。ratio で言い換えると、ratio > 2.54 であれば、Z-score > 2 となります。

参考

https://www.khanacademy.org/math/probability/statistics-inferential/normal_distribution/v/ck12-org-normal-distribution-problems-z-score

 

ヒートマップの色づけ (6): 中央値からの距離とシグナル値の関係

log変換した後、さらに、中央値(または平均値でも可)の距離に変換することで、コントロール群と実験群の差が見えやすくなります。このとき、もとのシグナル値との関係はどうなっているでしょうか?

中央値からの距離に変換したヒートマップ。
中央値からの距離に変換したヒートマップ。

シグナル値が異なる複数の遺伝子を色づけ

上記は、2つの遺伝子のシグナル値を、中央値からの距離に変換して色付けした場合の例です。中央値は、それぞれの遺伝子ごとに算出します。つまり、 Gene A の中央値は、4.5 であり、 Gene B の中央値は、9.5 です。それぞれの中央値を使って、距離を求めます。

したがって、 Gene A の WT1 は、 4 – 4.5 = -0.5 となり、 Gene B の WT1 は、9 – 9.5 = -0.5 となります。同様にすべてのサンプルについて算出し、色づけします。凡例にあるように、色づけの基準は、1つで複数の遺伝子に対応できます。

結果として、どちらも WT は -0.5 、 KO は、 +0.5 となり、色づけは、 Gene A も Gene B も同じパターンになります。

もとのシグナル値の高低に注意

このように、中央値からの距離に変換した場合、コントロール群と実験群の差は分かりやすくなりますが、 各遺伝子のシグナル値の高低は分からなくなります。

よって、ヒートマップからは、一見、差があるように見えても、マイクロアレイデータからシグナル値を確認すると、低い部分の変化であることもあります。シグナル値が低い場合は、ノイズである可能性もあるので、重要な遺伝子であれば、ヒートマップの色だけで判断せずに、もとの数値をチェックした方がよいでしょう。

 

ヒートマップの色づけ (5): 中央値からの距離に変換して色付け

「コントロールで低いところは、実験群で高いところは」のように色づけするには、 log変換した後、さらに、中央値からの距離に変換することが必要です。

中央値からの距離

中央値 (median) は、文字通り、真ん中の値です。例の場合、6サンプルのシグナル値は、9から10の間に分布しています。よって、その真ん中である 9.5 が中央値となります。(奇数サンプルの場合は、平均を取る必要はありません。5サンプルなら、3番目に大きい値が中央値です。)

中央値からの距離は、単純に「差」を取ればよいです。例えば、WT1のシグナル値は、9なので、9-9.5 = -0.5 となります。同様にKO1のシグナル値は、10-9.5 = 0.5 です。

中央値からの距離で色づけ。
中央値からの距離で色づけ。

中央値からの距離で色づけ

中央値からの距離に変換したら、色づけします。「どのくらいの差を何色にするか」は作者が決めます。ここでは、-1であれば、、+1であれば、としました。中央値からの距離なので、必然的に中央値に該当するものが、真ん中=となります。どのような色づけにしているかは、必ず凡例を付けましょう。

MeV の使い方で紹介している色づけ方法は、この作業を意味しています。

ここでは、一例として、中央値を使いました。これ以外に、平均値を使う場合もあります。また、標準偏差を求めて、平均値からSD何個分離れているか、という値に変換する場合もあります。