東大医科研ヒトゲノム解析センター 中 井 謙 太 knakai@ims.u-tokyo.ac.jp バイオインフォマティクス 遺伝子発見 東大医科研ヒトゲノム解析センター 中 井 謙 太 knakai@ims.u-tokyo.ac.jp
いまこそ DNA 言語学を!
遺伝子発見 (gene finding) ゲノムDNA配列中に存在する遺伝子の位置を同定する 狭義では、タンパク質コード領域(coding exons)を同定 開始コドンから終止コドンまで 広義では、RNA遺伝子やプロモータの同定なども含む
主なアプローチ 配列比較に基づく方法(実用的) 遺伝子に関する一般的知識に基づく方法(眉唾?) 遺伝子などの重要領域は進化的に保存される 近縁生物種のゲノム配列 cDNA、EST などのマッピング 参考:RepeatMasker (既知の繰り返し配列 DB と照合) 遺伝子に関する一般的知識に基づく方法(眉唾?) Ab initio 法とよばれる(第一原理に基づく)
配列比較に基づく方法 意外と難しい!(高等真核生物) 計算時間の問題 非常に長いイントロン(非常に短いエキソン) シークエンシングのエラー 同じスコアを与える複数のアラインメント 特に両端のあたりは保存されていない 計算時間の問題 EST の数は巨大 ゲノムは絶えず更新される
アラインメントと進化過程 二つの配列の距離 最適アラインメント(並置) GSSLRWA--PGP GS++ W P P 両者をつなぐ進化過程の変異数(のようなもの)で近似 置換、挿入、欠失 挿入・欠失は隣り合う位置に連続して二回起こるよりも一度長いものが起こる確率が高い 最適アラインメント(並置) 進化的に対応する位置を並べる 通常はアミノ酸配列レベルで比較した方が正確だが、偽遺伝子の判定は塩基レベルで行う GSSLRWA--PGP GS++ W P P GSTINWGGTPLP 配列1 配列2 置換 挿入・欠失 (ギャップ)
BLAST の原理 近似を導入して高速にアラインメントを求める 類似部分配列を検出して、その周囲を調べる 検索配列 LAALLNKCKTPQGQRLVNQWIK … PGG 18 PEG 15 PRG 14 PKG 14 PNG 13 PDG 13 PMG 13 PQA 12 .. 検索配列 325 DB 配列 290 SLAALLNKCKTPQGQRLVNQWIKQ .. +LA++L+ TP G R** +W+ TLASVLDCTVTPMGSRMLKRWLHM .. 検索配列 近縁配列
ギャップ付き BLAST 最初のBLASTはギャップを扱えなかった 新BLASTはギャップを扱い、しかも検索速度を高めた 2ヒット法 ギャップの部分はちぎれた 複数の類似領域の存在を確率的に評価した 新BLASTはギャップを扱い、しかも検索速度を高めた 2ヒット法 DPの計算範囲も限定
Ab initio 法 一番単純には ORF 検索 Search-by-Signal Search-by-Content 正統的な?方法だが、現在は単独では全く無力 ただし、原核生物の翻訳開始点を決めるのには有効 最近は、偽遺伝子を判定する手がかりとして脚光を浴びる Search-by-Content 現在の主流 コード領域、非コード領域の非ランダム性を利用
ORF 検索とその問題点 A G C T A G C T S A L Open Reading Frame ある読み枠で二つの終止コドンの間 開始コドンから拾う場合も ランダムなら 3 / 64 の確率で終止コドンにぶつかるので、100コドンも続けて終止コドンがないと偶然の可能性小 イントロンが入って、読み枠は変わってしまうため、この方法で遺伝子を探すのは難しい A G C T A G C T S A L
シグナルの表現法 TACGAT TATAAT Consensus (1) TATAAT GATACT Consensus (2) 一番単純なのは、コンセンサス表示 (1) で不一致を許さないと、2個ヒットで、4000bp に一回偶然のヒット (1) で1個不一致を許せば、3個ヒットで、200bp に一回 (1) で2個不一致を許せば、全部ヒットで、30bp に一回 (2) で不一致を許さないと、4個ヒットで、500bp に一回 (2) で1個不一致を許せば、全部ヒットで、30bp に一回 TACGAT TATAAT GATACT TATGAT TATGTT Consensus (1) TATAAT Consensus (2) TATRNT
重み行列法 (WM, PSSM) A C T A T A A A C T A T A A A C T A T A A Score = –41 g t -28 18 1 12 -15 -31 -12 -10 -18 -50 -11 -7 17 -17 10 -10 Score = –41 A C T A T A A -28 18 1 12 -15 -31 -12 -10 -18 -50 -11 -7 17 -17 10 -10 Score = 57 A C T A T A A -28 18 1 12 -15 -31 -12 -10 -18 -50 -11 -7 17 -17 10 -10 Score = –32 A C T A T A A
原核生物における関連シグナル プロモーター ターミネーター リボゾーム結合部位 非AUGスタートコドン UUUUU いろいろなシグマ因子依存性がある 主要なシグマ因子でも長いゲノム中での予測は難しい ターミネーター ρ因子非依存性のものはそれなりに リボゾーム結合部位 生物種によって差;そもそも正解があいまい 非AUGスタートコドン 既知の例が多いと統計的に処理できる UUUUU
大腸菌 (E. coli) Promoters Major sigma factor: s70 mRNA –1 +1 Transcriptional start sites –10 box (TATA box) mRNA –35 box spacer TTGACA 16〜19bp TATAAT 5〜9bp –1 +1
真核生物における関連シグナル プロモーター CpG アイランド スプライス・シグナル 既知の転写因子結合部位が集中しているところや、塩基組成に特徴のあるところを探すなど まだほとんど当たらない? 発現条件の読み取り問題:次回のテーマ CpG アイランド スプライス・シグナル
CpG アイランド 5’-CG-3’ の C はメチル化されやすい メチル化された遺伝子は不活化される ハウスキーピング遺伝子などではメチル化されていない
RNA Splicing Mature mRNA in Tissue A mRNAprecursor Mature mRNA in Tissue B
RNA Splice Sites (exon) A64 G73 | G100 T100 A62 A68 G84 T63 ... GT-AG Rule (exon) A64 G73 | G100 T100 A62 A68 G84 T63 ... ... 12Py N C65 A100 G100 | N (exon) Py: C or T
多すぎる偽物
統計情報に基づく遺伝子発見 コーディング・ポテンシャル 一般にシグナル情報に基づくより高い信頼度 タンパク質のコード領域はランダム配列ではない コドンの自由度とも関連した3文字周期性 たとえば、 6文字単語の使用頻度統計に基づき統計的に判定 一般にシグナル情報に基づくより高い信頼度 イントロンがあると、数段難しくなる(特にコード領域の両端は短いため)
コーディング・ポテンシャル 各フレームで6文字単語の出現頻度を、コード領域・非コード領域で集計して、未知配列の頻度とのずれを調べる 非コード領域に出現しにくい単語が固まって分布する領域は情報を多くもつ(情報量)
隠れマルコフモデルの考え方 ATATGGCAG… Hidden Markov model (HMM) 確率的な状態遷移を使って現象をモデル化 音声認識などの分野で成功 解くべき問題 何らかの規則性をもった記号列の確率的モデルを作る ATATGGCAG… C A G T C T A T A A G G A C
重み行列と HMM TATAAT TATAAA TATACA TTTAAT .... TaTAa(t/a) Start Mj End G: 0.0 C: 0.0 A: 0.6 T: 0.2 G: 0.0 C: 0.2 A: 0.5 T: 0.5 G: 0.0 C: 0.0
マルコフ連鎖 A G T C 現在の状態からある確率で次の状態に遷移(一次) 過去のn個の状態列から次の状態が確率的に決まる それ以前の履歴は考慮しないモデル 過去のn個の状態列から次の状態が確率的に決まる 過去1週間の天気から明日の天気を予測 初代GeneMarkはこの理論でコード領域をモデル化した A G T C
隠れマルコフモデル (HMM) . . . State 1 2 N P(R) 0.2 0.2 . . . 0.4 確率が二重に入る 各々決まった確率で記号を出力する状態間の遷移 どの状態にいるかは隠れており、出力記号だけを観測可能 . . . State 1 2 N P(R) 0.2 0.2 . . . 0.4 P(B) 0.8 0.2 . . . 0.4 P(Y) 0.0 0.6 . . . 0.2 Time 1 2 3 . . . T State q3 q1 q2 . . . q2 Color R B Y . . . R
モデルの定式化 l = (A, B, p) N 個の状態、M個の出力記号 A = {aij}, aij= Pr(qj at t + 1 | qi at t) B = {bj(k)}, bj(k) = Pr(vk at t | qj at t) p = {pi}, pi = Pr(qi at t = 1) 観測記号列 O = O1, O2, …, OT
3つの基本問題 モデルによる評価 遷移状態の推定 学習によるモデルの最適化 任意の記号列に対して与えられたモデルのスコアを計算 前向き・後ろ向きアルゴリズム 遷移状態の推定 任意の記号列に対するベストパス(遷移状態列)を計算 Viterbi アルゴリズムなど 学習によるモデルの最適化 与えられた記号列群に基づきパラメータを最適化 Baum-Welch アルゴリズムなど
HMM による遺伝子発見 遺伝子構造をモデル化 実際の学習データで最適化 モデルに一番良く合う領域を検索 Baum-Welch 法 Viterbi 法
原核生物の遺伝子発見 翻訳開始点の精密なモデルを利用 性能の正確な評価は難しい GeneHacker Plus (Yada ら) True negative はわかっても、False positive はわからない そもそも実験的に翻訳開始点が決まっている例は少ない GeneHacker Plus (Yada ら) 通常型 コード領域 SD 配列 スペーサ 開始コドン 終止コドン Start End 水平遺伝型 コード領域
真核生物の遺伝子発見 ... ... GENSCANにおける真核生物遺伝子モデル (1997) エキソン・イントロン長統計 同一状態持続確率も考慮 フレームの無矛盾性 (+) strand (-) strand E0 I0 5' UTR ... Einit Promoter E1 I1 Single Exon Gene Intergenic region ... Eterm Poly(A) Signal E2 I2 3' UTR Mirror image
プログラムの性能評価 Sensitivity (Recall) 予測正解数 真の正解数 Specificity (Precision) 予測正解数 真の正解数 Specificity (Precision) 予測した数 データの準備 そこに遺伝子があるとわかっても、遺伝子がないことは保証できない(人工データという手もあるが)
高等真核生物における予測精度 RepeatMasked human Chromosome 22 sequence Program Exon level Gene level # Sn Sp Annotation 3660 522 DIGIT 4513 0.695 0.563 654 0.123 0.098 FGENESH 5734 0.708 0.452 866 0.115 0.069 GENSCAN 6588 0.707 0.393 803 0.067 0.044 HMMgene 6840 0.629 0.336 1508 0.090 0.031
遺伝子発見の現代的意義 現在はとにかく、決定した配列から遺伝子を探したい 配列決定と代表的な生物の遺伝子機能解析が進めば、ほとんど「予測」プログラムの出番はなくなる? ゲノム情報をどこまで解釈できたかを定量的に示す例題としての重要性は残る 将来的には、統計的な方法よりも、細胞内での遺伝子認識をシミュレートする予測法が重要になるだろう
参考書 D. W. Mount, “Bioinformatics: Sequence and Genome Analysis”, Cold Spring Harbor Lab. Press (2001)(岡崎・坊農監訳、バイオインフォマティクス:ゲノム配列から機能解析へ、メディカル・サイエンス・インターナショナル (2002) R. Durbin 他(阿久津、浅井、矢田訳)バイオインフォマティクス:確率モデルによる遺伝子配列解析、医学出版 (2001) 美宅成樹・榊佳之編 応用生命科学シリーズ 9:バイオインフォマティクス、東京化学同人(2003) JST 技術者向け教材「配列解析の応用コース」(7月正式公開) http://www.sequence.info/JSTWebLearningPlaza/course/course.html
余談:Profile HMM マルチプルアラインメントのモデル化に使われる Pfamデータベースなど 各位置でマッチ状態と挿入/欠失状態 開始 Mj 終了 Ij Dj