DAVID – 遺伝子発現解析（マイクロアレイ解析, RNA-seq）

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

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

16hr の発現変動遺伝子

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

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

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

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

40hr の発現変動遺伝子

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

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

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

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

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

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

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

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

発現変動遺伝子

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

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

DAVIDによる機能解析

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

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” を検索して見ると、テーブルの下の方に見つかりました。

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

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

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

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

マイクロアレイ解析のフローチャート2: Next step

発現変動遺伝子の抽出後、まず、必要なことは、増加、減少した遺伝子（発現変動遺伝子）が、どのような遺伝子であるかを分析することです。この方法は、大きく分けて、次の2通りの方法のいずれかを用います。

A. （生物学的な）機能で見る。

B. 発現変動のパターンで見る。

A. （生物学的な）機能で見る

何系の遺伝子が多いのかを見ます。（膜タンパク？転写因子？EMT関連遺伝子？薬剤耐性遺伝子？脂質代謝？）これに該当する解析が、GO解析（DAVIDなどを使用）です。

B. 発現変動のパターンで見る

特定の注目する遺伝子に似た変動パターンの遺伝子を探します。また、どのサンプルの変動パターンが似ているかを見ます。（MDM4といっしょに増減した遺伝子はどれ？健常者に近いサンプルはどれ？）これに該当する解析が、クラスタリング(k-means, SOM)、ヒートマップです。

このA., B. の2通りの解析は、独立したものではなく、また、それだけでは終わりません。GO解析の結果をさらにヒートマップで表示して、発現変動変動のパターンを確認するケース(A. –> B.)や、特定の発現変動パターンの遺伝子を選択して、その機能をGO解析で確認するケース(B. –> A.) が通常です。それぞれの解析を単独で行っても効果的ではありません。（セルイノベーターの解析サービスでは、初めからこれらの解析サービスを含めて提供しています。）

解析はさらに続きます。

さらに続く解析。A. から B. や、B. の後に A. を行って遺伝子群を絞り込む。

最終的にマイクロアレイ解析の結果として、「特定の（生物学的な）機能を持ち、かつ、特定の発現変動パターンを示した遺伝子群」が得られます。（数十個から数百個の遺伝子）

DAVID 操作ガイド3

8. 解析結果１：アノテーションの解析結果を見る例です。ここでは、”Functional
Annotation Clustering” を見る方法を解説します。中程の ”Functional Annotation Clustering” のボタンをクリックします。

新しいウィンドウに解析結果が表示されます。Functional Annotation Clustering では、同じような機能を持った遺伝子群を１つのクラスターとして考え、スコアの高いクラスターの順に表示されます。１つのクラスターには、metabolic process, anion transport などの Gene Ontology (GO) が複数含まれています。アップロードした遺伝子リスト中、それぞれの GO をアノテーションに持つ遺伝子は、GO の隣の青いバーをクリックすることで見ることができます(Gene Report)。 P_Value は、0.05 (5.0E-2)以下が統計的に有意と判断される目安です。

青いバーをクリックして表示される Gene Report は、右上の “Download file” をクリックすることで、ダウンロードできます。

ダウンロードした Gene Report は、Excel で読み込むことができます。開くときに、選択対象を ”すべてのファイル” とします（“すべての読み込み可能なファイル”ではなく）。

9. 解析結果２：パスウェイを表示させる例です。 Annotation Summary Results の画面（7. の画面）の “Pathways (3 selected)” を選択します。

同じウィンドウの中にアノテーション情報として登録されているパスウェイデータベースの一覧が表示されます。パスウェイデータベースの横の数字は、アップロードした遺伝子リストのうち、パスウェイデータベースにヒットした割合と個数です。 “Chart” ボタンをクリックすると、 “Functional Annotation Chart” のウィンドウが表示されます。また、Chart ボタンの隣の ”青いバー” をクリックすると、 ”Functional Annotation Table” のウィンドウが表示されます。さらに、それぞれのウィンドウにおいて、各パスウェイ名をクリックすると、パスウェイの画像を表示できます。

パスウェイデータベース（ここでは KEGG ）に登録されている ”Apoptosis” や ”Cell Cycle” といった個々のパスウェイの中から、アップロードした遺伝子リストに含まれる遺伝子が載っている（マップされる）パスウェイを探すには、 ”Functional Annotation Chart” を参照します。スコアに関係なく、マップされるパスウェイを知りたい場合は、”Functional Annotation Table” を参照します。

Functional Annotation Chart ウィンドウには、アップロードした遺伝子リスト中の遺伝子を、統計的に有意な割合で含むパスウェイだけが表示されます。 “Term” 中の各パスウェイ名をクリックするとパスウェイが表示されます。

Functional Annotation Table ウィンドウには、アップロードした遺伝子リストのうち、パスウェイデータベースにヒットした遺伝子の一覧が表示されます。各パスウェイ名をクリックするとパスウェイが表示されます。パスウェイの数が多い場合は、ブラウザの検索機能を利用すると便利です。

マップされたパスウェイの例：下図は、”Cell Cycle” のパスウェイにマップされた場合の例です。アップロードした遺伝子リストに含まれていた遺伝子は、星印をつけて表示されます。各遺伝子をクリックすると、遺伝子について詳細な情報を見ることができます。

Reference
Huang et al. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nature Protocols (2009) vol. 4 (1) pp. 44-57.

DAVID 操作ガイド2

1. DAVIDにアクセスします。
http://david.abcc.ncifcrf.gov/

2. ショートカットメニューから、 “Functional Annotation” を選択します。

3. Functional Annotation Tool が表示されます。左側のメニューの “Upload” タブを選択します。遺伝子リストをアップロードするフォームが表示されます。

4. Excel で、解析結果を開き、ProbeSetID の列をコピーします。（下図は Affymetrix の場合です。それ以外は、GenbankAccesstion などを用います。）

5. Step 1: Enter Gene List の “A: Paste a list” の枠内に、コピーしたリストを貼付けます（この場合、B: のファイルを選択する必要はありません）。

6. Step 2: Select Identifier が、 ”AFFY_ID” になっていることを確認します。（ 4. で、GenbankAccession を選択した場合は、 Select Identifier に、”GENBANK_ACCESSION” を選択します。）さらに、 Step 3: List Type が、”GeneList” になっていることを確認します。最後に、 Step4: Submit List の “Submit List” のボタンをクリックします。（数十秒から数分程度時間がかかります。）

7. 遺伝子リストのアップロードが完了すると、次のような画面になります。対象となる生物種名を候補の中から選択し、”Select” ボタンをクリックします。右側に選択した生物種が反映された Annotation Summary Results が表示されます。