機能ゲノム学（第6回）

Transcription

1 RNAseqによる 定 量 的 解 析 とqPCR マイクロアレイなど との 比 較 東 京 大 学 大 学 院 農 学 生 命 科 学 研 究 科 アグリバイオインフォマティクス 教 育 研 究 ユニット 門 田 幸 二 (かどた こうじ) 1

2 自 己 紹 介 1995 年 3 月 高 知 工 業 高 等 専 門 学 校 工 業 化 学 科 卒 業 1997 年 3 月 東 京 農 工 大 学 工 学 部 物 質 生 物 工 学 科 卒 業 1999 年 3 月 東 京 農 工 大 学 大 学 院 工 学 研 究 科 物 質 生 物 工 学 専 攻 修 士 課 程 修 了 2002 年 3 月 東 京 大 学 大 学 院 農 学 生 命 科 学 研 究 科 応 用 生 命 工 学 専 攻 博 士 課 程 修 了 学 位 論 文 : cdnaマイクロアレイを 用 いた 遺 伝 子 発 現 解 析 手 法 の 開 発 ( 指 導 教 官 : 清 水 謙 多 郎 教 授 ) 2002/4/1~ 産 総 研 生 命 情 報 科 学 研 究 センター 産 総 研 特 別 研 究 員 2003/11/1~ 放 医 研 先 端 遺 伝 子 発 現 研 究 センター 研 究 員 2005/2/16~ 東 京 大 学 大 学 院 農 学 生 命 科 学 研 究 科 特 任 助 手 2007/4/1~ 現 在 東 京 大 学 大 学 院 農 学 生 命 科 学 研 究 科 特 任 助 教 アグリバイオインフォマティクス プログラム 2

3 Contents イントロダクション( 発 現 レベルの 数 値 化 ( 定 量 化 )) マイクロアレイ RNA-seq(ゲノム 配 列 既 知 のモデル 生 物 の 場 合 ) 前 処 理 ( 定 量 化 や 正 規 化 ) RPKM NAC FVKM など 他 のプラットフォーム(qPCRやマイクロアレイ)との 比 較 発 現 量 レベル(intra-sample) サンプル 間 比 較 レベル(inter-sample) 非 モデル 生 物 のRNA-seq 解 析 戦 略 de novo transcriptome assembly 発 現 変 動 コンティグ 同 定 3

4 トランスクリプトームとは ある 状 態 のあるサンプル( 例 : 目 )のあるゲノムの 領 域 ヒト 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 遺 伝 子 全 体 (ゲノム) どの 染 色 体 上 のどの 領 域 にどの 遺 伝 子 が あるかは 調 べる 個 体 ( 例 :ヒト)が 同 じなら 不 変 ( 目 だろうが 心 臓 だろうが ) 転 写 物 全 体 (トランスクリプトーム) 遺 伝 子 1は 沢 山 転 写 されている( 発 現 している) 遺 伝 子 4はごくわずかしか 転 写 されてない 4

5 トランスクリプトームとは ある 状 態 のあるサンプル( 例 : 目 )のあるゲノムの 領 域 光 刺 激 ヒト 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 遺 伝 子 全 体 (ゲノム) どの 染 色 体 上 のどの 領 域 にどの 遺 伝 子 が あるかは 調 べる 個 体 ( 例 :ヒト)が 同 じなら 不 変 ( 目 だろうが 心 臓 だろうが ) 転 写 物 全 体 (トランスクリプトーム) 遺 伝 子 2は 光 刺 激 に 応 答 して 発 現 亢 進 遺 伝 子 4も 光 刺 激 に 応 答 して 発 現 亢 進 5

6 トランスクリプトーム 情 報 を 得 る 手 段 光 刺 激 前 (T1)の 目 のトランスクリプトーム 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 これがいわゆる 遺 伝 子 発 現 行 列 光 刺 激 後 (T2)の 目 のトランスクリプトーム 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 マイクロアレイ ( 電 気 泳 動 に 基 づく 方 法 ) 配 列 決 定 に 基 づく 方 法 6

7 トランスクリプトーム 取 得 (マイクロアレイ) よく 研 究 されている 生 き 物 は 多 数 の 遺 伝 子 (の 配 列 情 報 )がわかっている 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 光 刺 激 前 (T1)の 目 の トランスクリプトーム Image courtesy of Affymetrix 蛍 光 標 識 ハイブリダイゼーション ( 二 本 鎖 形 成 ) わかっている 遺 伝 子 (の 配 列 の 相 補 鎖 )を 搭 載 した チップ メーカーによって 搭 載 されている 遺 伝 子 の 種 類 が 異 なる 搭 載 されていない 遺 伝 子 ( 未 知 遺 伝 子 含 む 例 : 遺 伝 子 4)の 発 現 情 報 は 測 定 不 可 7

8 マイクロアレイデータ 遺 伝 子 発 現 行 列 光 刺 激 前 (T1)の 目 のトランスクリプトーム 蛍 光 標 識 光 刺 激 後 (T2)の 目 の トランスクリプトーム ハイブリダイゼーション ( 二 本 鎖 形 成 ) 専 用 の 検 出 器 で 各 遺 伝 子 に 対 応 する 領 域 の 蛍 光 シグナ ル 強 度 を 測 定 ハイブリダイゼーション と シグナル 検 出 正 規 化 8

9 RNA-seqデータ 遺 伝 子 発 現 行 列 RNA-seq 光 刺 激 前 (T1)の 目 のトランスクリプトーム ゲノム 配 列 にマッピング -イメージ 塩 基 程 度 から なる 配 列 が 沢 山 ある - 実 際 - 数 百 万 個 の 配 列 が あり どの 遺 伝 子 に 対 応 するか 不 明 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 定 量 化 ( 例 : 生 の リード 数 をカウント) 正 規 化 ( 短 い) 配 列 を 読 んだものという 意 味 で(ショート)リードなどと 呼 ばれる 9

10 前 処 理 ( 定 量 化 や 正 規 化 ) 基 本 的 な 考 え サンプル 間 の 総 リード 数 の 違 いをいかに 補 正 するか 配 列 長 由 来 の 偏 り( 長 いほど 沢 山 sequenceされる)をいかに 補 正 するか ( 長 さの 異 なる 複 数 のisoformsが 存 在 する 場 合 にその 遺 伝 子 の 配 列 長 をいかに 定 義 するか) RPKM (Mortazavi et al., Nat Methods, 2008; ERANGEの 論 文 ) Reads per kilobase of exon per million mapped reads NAC (Griffith et al., Nat Methods, 2010; ALEXA-seqの 論 文 ) Normalized average coverage FPKM (Trapnell et al., Nat Biotechnol., 2010; Cufflinksの 論 文 ) Fragments per kilobase of transcript per million mapped fragments FVKM (Lee et al., Nucleic Acids Res., 2010; NEUMAの 論 文 ) Fragments per virtual kilobase per million mapped reads 本 質 的 に 同 じ Multiple isoforms 10

11 RPKM (Mortazavi et al., 2008) Reads per kilobase of exon per million mapped reads 対 応 reads ( 重 複 を 含 む 塩 基 配 列 数 )がマップされており そのうち 744 readsがa1bgという 遺 伝 子 のエ クソン 上 にマップされていて この 遺 伝 子 をRPKMという 単 位 で 定 量 化 すると82.9となる どうやって 計 算 してる? 11

12 RPKM (Mortazavi et al., 2008) RPM 正 規 化 (マイクロアレイなどと 同 じところ) Reads per million mapped reads サンプルごとにマップされた 総 リード( 塩 基 配 列 ) 数 が 異 なる 各 遺 伝 子 のマップされたリード 数 を 総 read 数 が100 万 (one million)だった 場 合 に 補 正 raw counts:all reads= RPM : 1,000,000 A1BGの 場 合 は 744 : 5,087,097 = RPM : 1,000,000 1,000,000 1,000,000 RPM raw counts all reads 5,087,097 RPKM 正 規 化 (RNA-seq 特 有 ) Reads per kilobase of exon per million mapped reads 遺 伝 子 の 配 列 長 が 長 いほど 配 列 決 定 (sequence)される 確 率 が 上 昇 各 遺 伝 子 の 配 列 長 を 1000 塩 基 (one kilobase)の 長 さだった 場 合 に 補 正 RPKM 1,000,000 1,000 raw counts all reads gene length 1,000,000,000 A1BG ,764 5,087,097 1,000,000,000 raw counts gene length all reads 12

13 NAC (Griffith et al., 2010) Normalized average coverage 1リードがx 塩 基 の 長 さとして 考 える 長 さ 補 正 ある 遺 伝 子 のaverage coverage (AC)は その 遺 伝 子 上 にマップされた 総 塩 基 数 を その 遺 伝 子 の 長 さ で 割 ったものなので AC 総 リード 数 補 正 raw counts x 744 x gene length 1,764 サンプルごとにマップされたリードの 総 塩 基 数 が 異 なるので マップされたリードの 総 塩 基 数 が10,000,000,000 塩 基 だった 場 合 に 補 正 10,000,000,000 NAC AC all reads x NAC 遺 伝 子 3 NACとRPKMは 本 質 的 に 同 じ だが NACのほうがより 厳 密 10,000,000,000 raw counts 10 RPKM gene length all reads RPKM? 13

14 おおざっぱにはこんな 感 じ 複 数 アイソフォーム 対 策 元 の 遺 伝 子 ( 補 正 後 )のgene length 値 をいかに 見 積 もるか? FPKM (Trapnell et al., Nat Biotechnol., 2010; Cufflinksの 論 文 ) 複 数 のisoformsの 長 さと 発 現 量 をもとに 発 現 量 で 重 みをつけた 平 均 値 を 採 用 gene1 長 さ 発 現 量 isoform1 69bp 20 isoform2 65bp 7 isoform3 60bp 5 補 正 後 のgene length *65 5* bp raw counts gene 定 数 length all reads 14

15 おおざっぱにはこんな 感 じ 複 数 アイソフォーム 対 策 元 の 遺 伝 子 ( 補 正 後 )のgene length 値 をいかに 見 積 もるか? FVKM (Lee et al., Nucleic Acids Res., 2010; NEUMAの 論 文 ) 共 通 部 分 のみを 利 用 して 他 の 遺 伝 子 にもマップされるものやisoform-specificなものは 使 わない gene1 isoform1 isoform2 isoform3 raw count ( 原 著 論 文 ではgNIR) = 3 raw counts gene 定 数 length all reads 15

16 おおざっぱにはこんな 感 じ 複 数 アイソフォーム 対 策 元 の 遺 伝 子 ( 補 正 後 )のgene length 値 をいかに 見 積 もるか? FVKM (Lee et al., Nucleic Acids Res., 2010; NEUMAの 論 文 ) 共 通 部 分 のみを 利 用 ( 他 の 遺 伝 子 にもマップされるものやisoform-specificなものは 使 わない) gene1 全 ての 可 能 なx bpの オリゴマー gene2 gene3 x-mers EUMA= 12 EUMA= 22 EUMA= 31 1,000,000,000 FVKM 3 12 all reads raw counts gene 定 数 length all reads 16

17 (gene length - x +1) 通 り おおざっぱにはこんな 感 じ 複 数 アイソフォーム 対 策 元 の 遺 伝 子 ( 補 正 後 )のgene length 値 をいかに 見 積 もるか? virtual length (Sultan et al., Science, 2008) 全 エクソンの 領 域 を 利 用 gene1 他 の 遺 伝 子 上 にはなくユニークに ヒットするx-merの 数 の 期 待 値 (theoretical total number of unique x-mers)を virtual length と 定 義 raw countsのほうも100%マッチでユ ニークにマップされるリード 数 のみ をカウント 定 数 raw counts gene length all reads 17

18 exon array 他 のプラットフォームとの 比 較 (vs. microarray) 発 現 量 レベル(intra-sample) 2,434 genes log 2 (NAC) Mortazavi et al., Nat Methods, 2008のFig. 3c Griffith et al., Nat Methods, 2010のSuppl. Fig. 9a(A) 18

19 exon array 他 のプラットフォームとの 比 較 (vs. microarray) サンプル 間 比 較 レベル(inter-sample) 217 genes 2,434 genes Roche 454 Mane et al., BMC Genomics, 2009のSuppl. Fig.の 下 半 分 log 2 (NAC) Griffith et al., Nat Methods, 2010のSuppl. Fig. 9b(A) 19

20 他 のプラットフォームとの 比 較 (vs. qpcr) 発 現 量 レベル(intra-sample) FVKM FPKM 27 genes RPKM RPKM Lee et al., Nucleic Acids Res., 2010のFig. 2 20

21 他 のプラットフォームとの 比 較 (vs. qpcr) サンプル 間 比 較 レベル(inter-sample) Griffith et al., Nat Methods, 2010のFig. 2 21

22 前 処 理 は 重 要 ( 遺 伝 子 発 現 行 列 作 成 時 ) 発 現 量 補 正 の 基 本 形 定 数 raw counts gene length all reads 発 現 量 レベル(intra-sample)の(プラットフォーム 間 ) 比 較 all readsの 項 はなくてもよい サンプル 間 比 較 (inter-sample)の 場 合 基 本 形 ではまだ 不 十 分 Bullard et al., BMC Bioinformatics, 2010 RPKM 補 正 でもまだ 発 現 変 動 遺 伝 子 が 配 列 長 の 長 いものに 偏 る t 統 計 量 gene length で 若 干 緩 和 される Robinson and Oshlack, Genome Biol., 2010 サンプル 中 の RNA 組 成 の 違 い による 影 響 は 甚 大 付 加 的 な 正 規 化 係 数 (TMM)を 掛 けることで 影 響 が 緩 和 される 22

23 Robinson and Oshlack, Genome Biol., 2010 RNA 組 成 の 違 い のイメージ 仮 定 全 4 遺 伝 子 長 さが 同 じ(gene lengthの 項 を 無 視 できるので) 遺 伝 子 4だけが 発 現 変 動 遺 伝 子 raw counts gene 定 数 length all reads サンプルS1 (all reads = 30) 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 サンプルS1 (all reads = 30) 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 補 正 サンプルS2 (all reads = 15) 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 サンプルS2 (all reads = 30) 遺 伝 子 1 遺 伝 子 2 遺 伝 子 3 遺 伝 子 4 補 正 結 果 :S1で 高 発 現 が1 個, S2で 高 発 現 が3 個 23

24 Robinson and Oshlack, Genome Biol., 2010 M M-A plot (R-I plot) TMM = A 横 (and 縦 ) 軸 で 上 位 下 位 のx (and y)%をtrim 残 りのデータでMのMean(TMM)を 計 算 24

25 TMM 補 正 するしないで 得 られたDEGセット 中 の 割 合 TMM 補 正 なし(Marioni et al., Genome Res., 2008) サンプルS1(Liver):22% サンプルS2(Kidney):78% TMM 補 正 あり(Robinson and Oshlack, Genome Biol., 2010) サンプルS1(Liver):47% サンプルS2(Kidney):53% 基 本 形 で 発 現 量 補 正 追 加 補 正 その 後 の 解 析 raw counts gene 定 数 length all reads TMM median etc 発 現 変 動 遺 伝 子 検 出 分 類 クラスタリング etc 25

26 マイクロアレイからの 知 見 ( 発 現 変 動 遺 伝 子 ;DEG) Jeffery et al., BMC Bioinformatics, 2006 別 のランキング 法 を 用 いると 違 った 結 果 に ランキング 上 位 x 個 の 集 合 方 法 1 方 法 2 一 致 は8-21%!( 再 現 性 低 い ) Kadota et al., Algorithms Mol. Biol., 2008,2009 既 知 のDEGは 全 体 的 に 発 現 レベルが 高 い ランキング 法 は t-test 系 とFold Change 系 に 大 別 でき この 間 の 比 較 で 再 現 性 低 下 遺 伝 子 発 現 行 列 作 成 時 に 用 いる 前 処 理 法 (Affymetrixの 場 合 )の 違 いの 影 響 もある ランキング 法 と 前 処 理 法 の 組 合 せが 大 事 感 度 特 異 度 が 高 いランキング 法 :Rank products or WAD 再 現 性 :WAD( 前 処 理 法 によらず) Hu and Xu, BMC Genomics, 2010 感 度 特 異 度 :WAD > t-test > Fold change > Rank products 上 位 1,000 遺 伝 子 までで 評 価 前 処 理 法 として 何 を 使 ったか 不 明 ( 公 共 DBはMAS-preprocessed dataが 大 半 で Rank productsとの 相 性 悪 い) 26

27 非 モデル 生 物 のトランスクリプトーム 解 析 de novo genome assembly 用 プログラム Velvet (Zerbino and Birney, Genome Res., 2008) ABySS (Simpson et al., Genome Res., 2009) EULER-SR (Chaisson et al., Genome Res., 2009) etc de novo transcriptome assembly 用 プログラム( 特 にIllumina) Multiple-k (Surget-Groba and Montoya-Burgos, Genome Res., 2010) Trans-ABySS (Robertson et al., Nat Methods, 2010) Rnnotator (Martin et al., BMC Genomics, 2010) Oases (Schulz and Zerbino, unpublished) 2010 年 の 夏 以 降 からtranscriptome 用 のものが 続 々 登 場 27

28 de novo transcriptome assembly 目 的 :(short)readsのデータから 転 写 物 ごとのコンティグを 得 る アセンブリの 基 本 戦 略 1. ( 計 算 を 軽 くするため ユニークなリード 配 列 の 集 合 にしておく) 2. de novo genome assembly 用 プログラムを 複 数 のk 値 で 実 行 転 写 物 の 場 合 はcoverageが 多 様 である 転 写 物 が 高 (or 低 ) 発 現 のときはhigh (or low) coverageであることを 意 味 する kを 大 きくすると 高 発 現 転 写 物 がアセンブルされる 確 率 が 上 がる( 低 感 度 高 特 異 度 ) kを 小 さくすると 低 発 現 転 写 物 がアセンブルされる 確 率 が 上 がる(がキメラも 増 える; 高 感 度 低 特 異 度 ) Rnnotator:k=19, 21,, 33 Multiple-k: k=19, 21,, 29 Trans-ABySS: k=26, 27,, 49 いろいろ 試 して できるだけ 転 写 物 のcoverageを 上 げる ( 読 んだリードの 長 さLによってkの 探 索 範 囲 を 変 更 ) 28

29 35 bpのsingle-endでkを 考 える 各 リードを 全 ての 可 能 なk-mer (k < 35の 任 意 の 値 ; 例 えばk=25)に 分 割 して 有 向 グラフを 作 成 リード1 リード2 リード3 入 力 データ 35 bp read1 :TGCCGACATGCATCCAAGTAGGAATCCTTAGCTTA read1_1 read1_2 read1_3 read1_4 read1_5 read1_6 read1_7 read1_8 read1_9 read1_10 read1_11 有 向 グラフ の 作 成 read1_1 read1_2 read1_3 隣 接 するノード 間 は(k-1) bp のオーバーラップ 全 リードのグラフ 情 報 をもとに 同 一 ノードをマージしたグラフ(de Bruijn graph)を 作 成 し オイラーパス 問 題 として 解 く(=コンティグを 得 る) 29

30 de novo transcriptome assembly 目 的 :(short)readsのデータから 転 写 物 ごとのコンティグを 得 る アセンブリの 基 本 戦 略 3. それぞれのk 値 を 用 いて 独 立 してアセンブルを 行 った 結 果 から 長 いコンティ グ 中 に 短 いコンティグが100%マッチになるものはマージしていくことでnonredundant setにする k=25のときの ある 長 いコンティグ TGCCGACATGCATCCAAGTAGGAATCCTTAGCTTA k=19のときの ある 短 いコンティグ CGACATGCATCCAAGTAGGAATCCTTA マージ TGCCGACATGCATCCAAGTAGGAATCCTTAGCTTA 30

31 de novo transcriptome assembly 目 的 :(short)readsのデータから 転 写 物 ごとのコンティグを 得 る アセンブリの 基 本 戦 略 4. キメラコンティグを 分 割 Martin et al., BMC Genomics, 2010のFig. 3b コンティグに 再 びリードをマップさせてforward 側 と reverse 側 で 明 確 にcoverageが 異 なるところで 分 離 31

32 非 モデル 生 物 の 比 較 トランスクリプトーム 解 析 戦 略 1. 比 較 する 複 数 サンプル(samples A and B) 由 来 のリードを 一 つにま とめたセットを 用 意 2. de novo transcriptome assemblyプログラムを 実 行 し コンティグの セット(transcriptome sequence)を 得 る 3. Transcriptome sequenceに 各 サンプル 由 来 リードを(Bowtieなどを 用 いて)マップ 発 現 量 の 定 量 化 はNEUMA 的 な 考 え 方 でunique readsの 結 果 のみ 採 用 ( 正 規 化 は 二 つのサンプル 由 来 リードがマップされているコンティグの 発 現 レベ ルのみを 考 慮 し TMM 正 規 化 のような 考 え 方 を 採 用 ) 32

33 要 求 されること( 例 :Trans-ABySS) ABySS Trans-ABySS Pysam Blat Biopython Bowtie BWA Perl modules Cython pyrex Python Samtools Git ncurses zlib curl openssl expat 全 部 インストールす るまで 待 て! configure make make install 33

34 謝 辞 東 京 大 学 大 学 院 農 学 生 命 科 学 研 究 科 清 水 謙 多 郎 教 授 嶋 田 透 教 授 グラント 若 手 研 究 (B)(H21 年 度 - ): マイクロアレイ 解 析 の 再 現 性 感 度 特 異 度 を 飛 躍 的 に 向 上 させるデータ 解 析 手 法 の 開 発 ( 代 表 ) 新 学 術 領 域 研 究 ( 研 究 領 域 提 案 型 )(H22 年 度 -): 非 モデル 生 物 におけるゲノム 解 析 法 の 確 立 ( 分 担 ) 34