マイクロアレイデータの解析例 1.7 (変動パターン、ヒートマップ)

前回は、先に特定の機能の遺伝子をピックアップし、その後、ヒートマップで変動パターンをチェックする方法を紹介しました。これに対し、先に変動パターンをチェックするという方法もあります。

変動パターンで抽出したのち、ヒートマップを確認する

解析例1のデータは、各タイムポイントごとに、3つの比較を行いました。

  • 16hr: sample1/control1
  • 24hr: sample2/control2
  • 40hr: sample3/control3

この3つの比較のいずれかの比較において、変動している遺伝子増加:ratio>2 かつ Z-score > 2、減少:ratio < 0.5 かつ Z-score < -2)の遺伝子を抽出して、ヒートマップを作成しました。ヒートマップの色付けなどは、これまでと同様です(logFC の中央値からの距離)。図には、2678個の変動遺伝子が含まれています。縮小しているため、画像のサイズは前回のアポトーシスのヒートマップと同じですが、こちらが多くの遺伝子が含まれています。

いずれかの比較で変動した遺伝子のヒートマップ。変動パターンを確認できる。
いずれかの比較で変動した遺伝子のヒートマップ。変動パターンを確認できる。

ヒートマップを見ることで、変動パターンを把握できます。例えば、16hrと24hrで共通に増加した遺伝子があるのか、24hrと40hrで共通に減少した遺伝子があるのか、あるとしたら何個くらいか、すべてヒートマップから見て取れます(数値は、おおよそではありますが)。

「共通の変動遺伝子が何個ある」というのをチェックするために、「ベン図」が用いられることもありますが、比較の組み合わせが、2つのときはよくても、今回のように3つ以上の場合は、複雑になりますので、ヒートマップで確認することをお勧めします。(また、マイクロアレイデータの場合、データの見せ方にもよりますが、変動遺伝子の個数には、あまり意味がない場合が多いです。)

特徴的な変動パターンから機能解析 (DAVID)

実験条件から、いくつかの変動パターンが予想されるかもしれませんが、もし、特徴的なパターンがあるのなら、ヒートマップで見えるはずです。このデータの場合、24hrだけで増加、40hrだけで増加という部分が特徴的なようです。(=なぜか、16hrだけ増加は少ない?)

すべての変動パターンが見られる。
すべての変動パターンが見られる。

ヒートマップの結果、気になる変動パターンが見つかれば、その部分だけを取り出して、再び、DAVIDなどで確認します。そうすれば、24hrだけで増加する遺伝子に、何系の遺伝子が多いのか、判断できます。しかしながら、その結果、生物学的な意味があるのかは、個人の仮説次第です。(ヒートマップから結論が得られるわけではありません。)

参考:

  • ヒートマップの拡大図
  • 発現変動遺伝子のエクセルファイル。エクセルのフィルター機能を使えば、特定のパターンの部分だけ取り出せます。
 

マイクロアレイデータの解析例 1.6 (機能で抽出、ヒートマップ)

解析例1の続きです。あらかじめ、アポトーシス関連の遺伝子に変動が見られることが予想されるのであれば、 アノテーションをもとにそれらの遺伝子をピックアップして(抽出して)、ヒートマップを書いてみても良いでしょう。

アポトーシス関連遺伝子を抽出してヒートマップを作成

アポトーシス関連遺伝子を抽出して、ヒートマップを作成した例を示します。(ヒートマップの作成方法は、リンク先を参照してください。)

アポトーシス関連遺伝子のヒートマップ(縮小図)。
アポトーシス関連遺伝子のヒートマップ(縮小図)。

ここでは、概略を見るために、発現変動の有無に関わらず、アポトーシス関連のアノテーションを持つ遺伝子を全て抽出しています。また、ヒートマップの色付けは、各遺伝子ごとにログ変換されたシグナル値の中央値からの距離を色付けしています。結果、色の意味としては、その他のサンプルより低いものは緑に、中程度のものは黒に、その他のサンプルより高いものは赤に見えています。(横での比較(サンプルどうし)のみに意味があります。縦での比較(遺伝子どうし)に意味はありません。)

遺伝子名を確認するには、拡大図を参照してください。

の順に高くなるため、左から、の順に色付けされていれば、いずれのタイムポイントにおいても減少する傾向にあると言えます。

逆に、左から、の順に色付けされていれば、いずれのタイムポイントにおいても増加する傾向にあると言えます。

縮小図を見た印象としては、アポトーシス関連遺伝子全体としては、それほど、タイムポイントによる違いが見えている部分はないようです。(DAVIDでも20個程度しか見つからないようです。)また、増加減少の傾向が見えている部分も変動幅は大きくないように見えます。さらに、コントロールのサンプル内でもタイムポイントによって増減しているように見える部分があるのも気になります。

アポトーシス関連遺伝子に見られる変動パターン。
アポトーシス関連遺伝子に見られる変動パターン。

色だけでは、正確な変動は分からないため、最終的には、シグナル値、ratio, Z-score の値を確認します。ヒートマップは、あくまで、おおまかな傾向を掴むものと考えてください。

 

マイクロアレイデータの解析例 1.5 (機能解析の続き)

解析例1のデータのほかのタイムポイントの発現変動遺伝子も、DAVIDで確認してみましょう。

16hr の発現変動遺伝子

前回の 24hr の発現変動遺伝子と同様に、 ratio と Z-score によって判定された発現変動遺伝子を DAVID で確認してみます。

16hr の発現変動遺伝子を DAVID で解析した結果。
16hr の発現変動遺伝子を DAVID で解析した結果。

24hr のときとは異なり、かなり膜系のタンパク、GPCRなどのレセプターの変動が多い印象です。

16hr の比較におけるケラチンのアノテーションクラスター。
16hr の比較におけるケラチンのアノテーションクラスター。

また、ケラチンを含むアノテーションクラスターも上位にあります。Enrichment Score = 1.88  で有意といえます。(1.3を超えているので。)

16hr の比較におけるアポトーシスのアノテーションクラスター。
16hr の比較におけるアポトーシスのアノテーションクラスター。

一方、アポトーシスを含むアノテーションクラスターは、24hr のときと同様に、それほど多くありません。

40hr の発現変動遺伝子

同様に 40hr のタイムポイントにおける発現変動遺伝子も確認してみました。

40hr の発現変動遺伝子を DAVID で解析した結果。
40hr の発現変動遺伝子を DAVID で解析した結果。

トップは、 T-cell activation です。また、16hr と同様に膜系、レセプターが上位にきていますが、Enrichment Score は下がっています。(他の機能に分散している?)

40hr の比較におけるケラチンのアノテーションクラスター。
40hr の比較におけるケラチンのアノテーションクラスター。

ケラチンの変動はあまり見られなくなっています。16hr, 24hr, 40hr とだんだん変化がなくなるということでしょうか。

40hr の比較におけるアポトーシスのアノテーションクラスター。
40hr の比較におけるアポトーシスのアノテーションクラスター。

アポトーシス関連遺伝子の変動は、やはり少ないようです。

全体を通して、個人的には、アポトーシスの影響は小さいように感じました。また、増殖系への影響も少なそうです。膜タンパクの影響もありそうですが、24hr でいったん少なくなるのが解釈に困るような気がします。

今回の解析の内容は、あくまで主観的なものです。論文の内容を確認していませんし、フェノタイプに変化があるなら、必ずしも、DAVIDのスコアが高い必要はないと考えます。

 

マイクロアレイデータの解析例 1.4 (発現変動遺伝子、機能解析)

正規化後のシグナル値から、ratio と Z-score を算出して、発現変動を判定します。その結果、発現変動している遺伝子のリスト(例えば、エクセルファイルのようなテーブル)が得られます。ここでは、 sample2/control2 = 24hr どうしを比較したデータを使用しています。

発現変動遺伝子のテーブル。
発現変動遺伝子のテーブル。

発現変動遺伝子

この例では、下のMAプロットに示されたような部分が変動していると判定されています。判定条件としては、 ratio > 2 かつ Z-score > 2 (または ratio < 0.5 かつ Z-score < -2)という少し厳しめの条件にしています。増加しているもの、減少しているもの合わせて、1032個を変動していると判定しています。

発現変動遺伝子を色づけしたMAプロット。
発現変動遺伝子を色づけしたMAプロット。

まずは、変動している遺伝子の機能を調べてみましょう。機能を調べるには、データベースの DAVID が便利です。

DAVIDによる機能解析

DAVIDに変動している遺伝子の ID のリストをアップロードすると、そのリストに含まれる遺伝子にどのような機能の遺伝子が多いのか(膜タンパクが多いのか、転写因子が多いのか、アポトーシスに関連する遺伝子が多いのか)、確認することができます。操作方法は、リンク先を参照してください。

DAVID の Functional Annotation Clustering を実行すると、下記のようなテーブルが得られます。

DAVID の Functional Annotation Clustering の結果。
DAVID の Functional Annotation Clustering の結果。

左上の表示を確認すると、1032個の変動遺伝子のうち、558個が認識されたことがわかります。(ここでは、EntrezGeneID をアップロードしています。)変動しているもののうち、アノテーションがついていないもの、機能が未知のものが多かったようです。

テーブルの最初 (Annotation Cluster 1) から、signal peptide, glycoprotein 、つまり、シグナルを伝達する分子、糖タンパクの遺伝子が多かったことが分かります。だいたい、この機能のグループ(クラスター)は、どんなデータでも上位に来ることが多いです。これらの用語は、漠然とした分類のアノテーションなので、あまり、情報としての価値はないかもしれません。

アノテーションのクラスターとして、次は、ホルモンのようですが、遺伝子数はそれほど多くありません。

入手元のデータのサマリーによると、変異の結果、アポトーシスが活性化され、ケラチノサイトの分化が抑制とあります。

TRPV3 dysfunction may increase apoptotic activity, inhibit keratinocyte differentiation and disturb the intricate balance between proliferation and differentiation state of keratinocytes in the skin.

ブラウザの検索機能で、”apoptosis” を検索して見ると、テーブルの下の方に見つかりました。

apoptosis を含む Annotation Cluster.
apoptosis を含む Annotation Cluster.

確かに変動しているようですが、数は少ないようです。(1032個のうち、20個程度)Enrichment Score も 0.37 と低いです。

“keratinocyte” を検索すると、 “keratinocyte differentiation” も見つかりますが、これも数は少ないです。

"keratinocyte differentiation" を含む Annotation Cluster.
“keratinocyte differentiation” を含む Annotation Cluster.

24hr の比較では、あまり、変異体の影響が出ていないように思われます。ほかに 16hr, 40hr のタイムポイントもあるので確認してみる必要がありそうです。

 

マイクロアレイデータの解析例 1.3 (散布図、正規化の影響)

ボックスプロットを確認したら、次は散布図 (scatter plot) も確認してみましょう。

散布図

例として、 control2 と sample2 を比較した場合の散布図を示します。 ratio > 2 のプローブ(=遺伝子)を、 ratio < 0.5 の遺伝子をに色づけしています。正規化後データを用いています。

正規化後データの散布図。
正規化後データの散布図。

散布図の広がり方から、平均的なデータのように見えます。(がんのサンプルや、変動が大きいデータでは、もっと点が全体に散らばって見えます。)

正規化の影響

rawとして、正規化前のデータの散布図も作成しました。

raw データの散布図。シグナル値の分布に偏りが見られる。
raw データの散布図。シグナル値の分布に偏りが見られる。

よく見ると、rawデータと、正規化後のデータでは、(ratio で判定した場合)変動ありと判定される遺伝子に異なる部分があることがわかります。特にシグナル値の高い部分です。

散布図の左下から右上に引かれた赤線は、 y=x を意味しています。raw データの散布図は、集団が y=x より下に膨らんで見えています。そのため、rawデータをそのまま用いると、sample2で減少した遺伝子が多く見つかり、増加した遺伝子は少なく見つかることになります。

一方、正規化後データの散布図では、点の中心が y=x 上に載っていることが分かります。(=ほとんどの遺伝子が変動していない。偏りがない。)

このように散布図を確認すると、raw, 正規化前のデータに偏り(バイアス)がないか、また、正規化後のデータから偏りが解消されているのか確認できます。

スコアだけで判断して、結果を誤って解釈しないよう、散布図を必ずチェックするようにしましょう。