DAVIDを使ってマイクロアレイデータを解析する 動画うp

お久しぶりです…。
こちらでの報告がだいぶ遅れてしまいましたが、とにかく前回勤務でDAVID の動画をアップいたしました。こちら。http://togotv.dbcls.jp/20120927.html#p01
今回、例として使用する遺伝子リストについて、これぞというものがなかなか見つからなかったため、結局前回動画と同じ実験のデータから GEO2R で取得したものを用いましたが、遺伝子数が 250 とやや少ないことや、発現が増加したものと減少したものが混じっていることなどのため、機能アノテーションクラスタでは最上位のクラスタでも Enrichment Score が 2.16 と小さく、やや例としては相応しくなかったかもしれません。その辺の修正も含めて、続編を作ることになりました。DAVID による解析は DBCLS の講習会でもしばしば用いられ、その際は DAVID で選択可能な多くのアノテーションデータベースのうち、Gene Ontology に限って使用して説明しているそうです。次回動画では、Gene Ontology の説明もしながら、より実践的な使い方について説明していけたらよいと思います。

遅れてすみませんでした。後回しはダメですね…

DAVIDを使ってマイクロアレイデータを解析する 調査5〜8

4回分、特にまとまりがないですが。

Functional Annotation Clustering/Chart で term をクリックすると多くの場合その term についてのより詳しい情報が掲載されたページ、普通はそのアノテーションの本家のサイトの該当ページに移動する。ものによっては本家のサイトのトップページに飛ばされるだけだったり、そもそもリンクが貼られていなかったり。
KEGG pathway など一部のパスウェイデータベース由来の term をクリックすると、リスト中の遺伝子をそのパスウェイ上にマッピングした画像を表示できる。こんな↓

demolist1 でKEGG の MAPK Signaling Pathway を見たものの一部。星付きがリストの遺伝子。

Background について
デフォルトでは、対応する生物の遺伝子のうち、選択されているアノテーションが一つ以上ついているものの全体。ゲノムワイドな解析の場合はこれを用いれば良い。また、用いているアレイに対応したリストも用意されている。これらでは気に入らない場合は自分でアップロードして使うこともできる。

調べることはだいたい終わったので、動画内で例として用いる遺伝子リストについて検討していましたが、なかなかしっくりくるのに出会えず結局前作 GEO2R で使ったこれ Effects of Cigarette Smoke on the Human Oral Mucosal Transcriptome を今回も使うことにしました。
喫煙者女性・非喫煙者女性間で有意な発現差がある遺伝子上位 250 個のリストを GEO2R で得、DAVID に投げました。Functional Annotation Clustering の結果、"oxidoreductase" を含むクラスタや "xenobiotic metabolic process" を含むクラスタが上位にできます。Enrichment Score は 2.16 と高くありませんがまあいいでしょう。

以下動画の筋書き。あまり細かいことはやらない方針にしました。

Start Analysis をクリック
リストをアップロード
Gene ID の識別子を選択(AFFYMETRIX_3PRIME_IVT_ID)
List Type を選択(Gene List)
 Background の説明
Submit List
List Manager の "Rename" "Combine" などを説明
Show Gene List をクリック
 適当な Gene Name をクリックして Gene Report を表示
 RG をクリックして 関連遺伝子を表示、カッパ係数について説明
Functional Annotation Tool へ移動
 アノテーションの選択。デフォルトのまま変更はしない
 Functional Annotation Clustering を選択
  Enrichment Score を説明し、上位のクラスタの特徴を指摘
  タームを選択して元ページヘ
  BioCarta か KEGG のタームを選択しパスウェイマップを表示
  RT をクリックし関連タームを表示
  グラフのバーをクリックしそのタームをアノテーションに持つ遺伝子を一覧
  2D view を表示
 Functional Annotation Chart を選択、説明
 Functional Annotation Table を選択、説明
Gene Functional Classification へ移動、説明

DAVIDを使ってマイクロアレイデータを解析する 調査2〜4

3回分です。

DAVID で使える解析ツールのまとめ。

Functional Annotation Tool
・Functional Annotation Clustering
リスト中の遺伝子が持っているアノテーションクラスタリングしたものを表示。リストの遺伝子の機能を概観できる。クラスタリング手法は下記。
・Functional Annotation Chart
リスト中で出現する回数に有意性が高い(大雑把に言えば、リスト中の多くの遺伝子が持っている)アノテーションを順に表示。
・Functional Annotation Table
リストの各遺伝子の持っているアノテーションを遺伝子ごとに一覧

Gene Functional Classification Tool
アノテーションをもとにリストの遺伝子をクラスタリング。"RG" (Related Genes) をクリックすることで、リスト中での有無にかかわらずそのクラスタの機能に関連する遺伝子を表示できる。

Gene ID Conversion Tool
リストの遺伝子 ID を、別の ID に変換 (Affymetrix ID から DAVID の ID に変換、など)

Gene Name Batch Viewer
リストの遺伝子 ID を遺伝子名に変換

NIAID Pathogen Annotation Browser
病原体の遺伝子検索。病原体のリストからいくつか選んでキーワードを入れると、キーワードに関連した遺伝子のリストを表示する。このリストを DAVID の他の解析に用いることができる。

アノテーションのリストでの出現回数の有意性を示すのに使われている P-value は modified Fischer exact p-value です。
http://david.abcc.ncifcrf.gov/helps/functional_annotation.html#fisher
全遺伝子の中であるアノテーションを持っているものの割合をバックグラウンドに、リスト中でのそのアノテーションの出現回数の有意性を計算します。
クラスタには、その有意性を表す指標として enrichment score が定義されています。所属する各アノテーションの p-value の、-log の幾何平均です。

クラスタリングについて。
2 遺伝子間または 2 アノテーション間の距離はκ係数というもので定義しています。
http://david.abcc.ncifcrf.gov/helps/linear_search.html#kappa
http://www.med.osaka-u.ac.jp/pub/kid/clinicaljournalclub12.html
持っているアノテーションが似ている遺伝子同士は近い
持たれている遺伝子が似ているアノテーション同士は近い
というのが大雑把なところ。

クラスタリングアルゴリズムは独自のものを使っており、
・各属性が複数のクラスタに所属しうる
クラスタの個数を自動決定
などの特徴があります。
 詳細についてはこちら。
http://david.abcc.ncifcrf.gov/helps/functional_classification.html#clustering
各点間のκ係数を計算する
κ係数が閾値(例えば 0.35)以上の点が、閾値(例えば 2 個)以上個存在する点を initial seed とする
ただし、クラスタ内の seed でない点同士の距離が閾値に満たないものが半数を超える場合は除外する
閾値(例えば 50%) 以上のメンバーを共有するクラスタをマージする

50% 以上のメンバーを共有するってどういうことでしょう。誰にとっての 50% なのか。

A∩B / A >= 0.5 かつ A∩B / B >= 0.5

なのかそれとも

A∩B / A∪B >= 0.5

なのか。
A={a, b, c, d}, B={c, d, e, f} に対して、前者ならA, B はマージされ後者ならされない。

また、いずれにしても順序に依存があります。
前者なら
{a, b, c}
{b, c, d}
{c, d, e}
{d, e, f}
から最終的にできるクラスタが、{a, b, c, d, e, f} だけになる場合と {a, b, c, d, e} と {d, e, f} などのようになる場合があります。

後者なら
{a, b, c}
{b, c, d}
{c, d, e}
から {a, b, c, d} {c, d, e} もしくは {a, b, c} {b, c, d, e} ができます。

ヘルプにはこの辺のことが書かれていないようです。うーん。

DAVIDを使ってマイクロアレイデータを解析する 調査1

更新がだいぶ遅れました、すみません。

新しいお仕事として、マイクロアレイのデータ解析ツール DAVID の動画のアップデートをいただきました。旧作はこちらhttp://togotv.dbcls.jp/movie/090925david_f.html

例えば GEO2R を使って発現量に変化があった遺伝子を見つけてきたら、その遺伝子群がどのような機能を持つか、どのような経路を構成するかなどが次に知りたいことになってきます。DAVID (The Database for Annotation, Visualization and Integrated Discovery) では、遺伝子の ID リストを与えることで、Gene Ontology などの様々なアノテーションに基づいて遺伝子のクラスタリングなどを行うことができます。
基本的な使い方はとても簡単で、ID のリストを欄にペーストして、解析に用いるアノテーションを選択し、実行するだけです。

旧作動画からの大きなインターフェースの変更はなく、作りなおすだけならそんなに時間はかからないと思われます。だからといってあまりそのままでも面白くないので適当に付け加えます。

NCBI GEOの使い方5 GEO2Rの使い方 動画編集2〜4

 本日動画編集を完了し、アップいたしました。こちら
 前々回の作業中に、NCBI 本家によるチュートリアル動画の存在が判明(こちら)。後出しになってしまったのは若干悔しいですが、まあ気にしません。
 さて本日のアップ作業ですが、YouTube にアップしている最中に作業マシンが落ちるというたいへんたのしいできごとがありました。…ブルースクリーンは心臓に悪いですね。強制終了から起動しなおしたところ何事もなかったようにしていましたが。恐ろしい。
 さて、これで、私が昨年放牧を始めてからずっと続けてきた NCBI GEO シリーズが完結しました。長かったですねえ。で、今日チケットを色々見てみましたが、まだ次回作のネタは決まっておりません。次回決めます。

NCBI GEOの使い方5 GEO2Rの使い方 調査6・7, 動画編集1

動画の中で例として用いる実験データを吟味しておりました。今まではかなり適当に選んでいましたが、今回は感覚的に意味がわかりやすいものを選ぼうと思い、時間をかけました。その結果、
GSE17913 Effects of Cigarette Smoke on the Human Oral Mucosal Transcriptome
に決定。喫煙者と非喫煙者で口腔粘膜の遺伝子発現を比較した実験です。異物代謝に関わる遺伝子(ex. cytochrome P 450)や酸化ストレスに関わる遺伝子が、喫煙者では非喫煙者に対して高発現している、というようなわかりやすい結果が出ています。
これで動画を撮れる状態になりましたので、撮影・編集を開始しました。以下大まかな筋書きです。

GEO トップから series -> cigarette smoke で検索

GSE17913 を選択

never_smoker_female と smoker_female をそれぞれグループに

value distribution で発現分布を確認、export からテキストで表示

デフォルトの設定で top 250

各カラムの説明

CYP1B1 の発現量グラフを表示

select columns でカラムを変更(logFC を外し GO.Function を入れる)

save all result で結果をテキストで表示、保存

options の各項目の説明

profile graph で適当な遺伝子(NFE2L2; 酸化ストレスにより活性化する転写因子)の発現グラフを表示。遺伝子の ID は View data for GPL**** から探す

R script を表示

今日の放牧で、options の説明に入るあたりまで編集しました。長さは10分くらいになりそうです。

NCBI GEOの使い方5 GEO2Rの使い方 調査5

前回保留にした、多重検定の補正の話です。

 ある有意水準 α で検定を N 回行うと、本当は有意差がないのに誤って帰無仮説を棄却してしまうことが期待値 Nα で起こってしまいます。そこで、余計に帰無仮説を棄却しないように個々の検定の有意水準をいじってやる必要がある、というモチベーションです。
 有意水準をいじるのは p 値をいじるのと同等です。 GEO2R では調整した p 値を表示できます。
 以下では GEO2R (というか limma)で使える6つの補正手法について説明します。
 なお
 N 検定の回数
 α 有意水準
  元の p 値を昇順に並べたときに i 番目となる p 値
 つまり
  補正後の p 値
 とします。

・Bonferroni の方法
 各検定の有意水準を α/N とします。これは各 p 値を N 倍するのと同じ。

 個々の帰無仮説が棄却される確率が α/N 以下ですから、1つでも帰無仮説が棄却される確率は次式により α 以下となります(途中の近似は α<
 「棄却すべき帰無仮説を棄却できなくても、棄却すべきでない帰無仮説を棄却してしまう確率を α で抑えたい」というようなスタンスです。
 全遺伝子の中から2群間で発現量に差があるものを探す場合だと、 N=30000 くらいになり得ます。 α=0.05 とすれば個々の検定の有意水準が 1.7*10^-6 と、非常に小さくなってしまいます。かなり厳しいです。

・Holm の方法
 Holm の方法について説明する前に、 step-up 方式と step-down 方式について説明します。
 Bonferroni の方法では各検定の有意水準を等しく N で割りましたが、以下で紹介する各方法では、 p 値を昇順に並べ、その順位に依存する値で調整します。その際に step-up 方式では p 値の大きい方から調べていき、調整後の p 値が有意水準より小さくなるものが見つかったら、以降の帰無仮説を全て棄却します。
step-down 方式ではその逆で p 値の小さい方から調べていきます。調整後の p 値が有意水準より大きくなるものが見つかったら、以降の帰無仮説を全て採択(あるいは判定保留) します。
 それで、 Holm の方法ですが、これは step-down 方式で各有意水準を N-i+1 で割るものです。つまり

 p 値は 1 を超えてはいけないので外側の min は明らかです。内側の max は step-down 方式の性質で、一つ順位が上の帰無仮説が採択されたら以下の帰無仮説も全て採択されることを表しています。次の Hochberg の方法の説明のあとで合わせて例を示します。
 p 値が大きくなるにつれて有意水準が緩くなりますが、間違えて棄却する確率は α で抑えられます。

・Hochberg の方法
 Holm の方法の step-up 版です。

 内側の min は step-up 方式の性質で、一つ順位が下の帰無仮説が採択されたら以上の順位の帰無仮説も全て棄却されることを表しています。
 では例です。帰無仮説が 5 つあり、 p 値が
 
 であるとします。この場合 Holm, Hochberg のそれぞれの方法による調整 p 値は以下の表のようになります。

i N-i+1 (Holm) (Hochberg)
1 0.011 5 0.055 0.055 0.044
2 0.012 4 0.048 0.055 0.044
3 0.02 3 0.06 0.06 0.044
4 0.022 2 0.044 0.06 0.044
5 0.07 1 0.07 0.07 0.07

 
・Benjamini & Hochberg の方法 (BH 法)
 棄却される仮説のうち、本当は差がないものの割合を FDR(False Discovery Rate) といいます。これを有意水準以下に抑える手法です。
 これまでと違い、「棄却すべきでない帰無仮説を棄却してもいいから、棄却すべき帰無仮説はできるだけ棄却したい」というスタンス。
 step-up で、各有意水準を i/N 倍します。

 有意水準 だとすると、本当は差がないのに有意水準を下回る帰無仮説の個数の期待値は です。この場合実際に棄却される帰無仮説の個数が i 個ですから、 となります。したがって、 とすれば FDR も α 以下で抑えられることになります。
 制限の度合いが適度で、マイクロアレイデータの補正においては最も広く使われています。 GEO2R のデフォルトでもこの方法が使われます。

・Benjamini & Yekutieli の方法
 BH 法を厳しくした方法です。
 
 とします。 k は N だけに依存する定数で、 1 以上の値をとります。この上で
 
 BH 法と比べて、全ての有意水準が 1/k 倍され、厳しくなっています。

・Hommel の方法
 少し複雑です。まず
 
 で与えられる j を求め、有意水準を α/j とするものです。
 Hochberg の方法より多くの帰無仮説を棄却することが知られていますが、実際にはあまり使われないようです。


参考
Bonferroni法、Holm法、False Discovery Rate 大阪大学腎臓内科
FDR を制御する多重比較法の性能評価
Multiple comparison procedures based on marginal p-values