生物情報科学科 情報基礎実験 3年生冬学期水曜日 (森下研究室担当)
本演習では次世代シーケンサーデータを解析する基本的なプログラムの作成を目指します。
ホームページ (質問があるたびに大切なことは追記してゆく形で更新します)
http://mlab.cb.k.u-tokyo.ac.jp/~moris/lecture/2009-winter-experiments/schedule.htm
5つの課題を3つに減らして、各々の〆切を延長しました。
(0) C/C++ および Java のプログラミング環境をインストール
l 資料 http://www.xerial.org/wiki/lecture/2009/Setup
l 使用言語は C/C++, Java どれでも構いません。
(1) Suffix Array の実装
l
Induced Sorting を利用して線形時間アルゴリズムを実装
できれば講義ノートを参考に自分で実装を考える。
参考論文:Nong,
Zhang, and Chan. Two
Efficient Algorithms for Linear Suffix Array Construction.
Appendix に Cプログラムあるので、自作が難しいようでしたら、このプログラムを解読。提出用のプログラムにはできるだけコメントを入れてほしい。
なおこのプログラムの動きを例題を使って解説した資料はココ。
l チェック用のアルゴリズム(たとえば Larsson-Sadakane’s O(n log n)-time algorithm)を実装して、同一の答えが出ることを確認する。C++ で書いた Larsson-Sadakane’s O(n log n)-time algorithmの プログラム Larsson_Sadakane.cpp ( Visual C++ 等のC++ プログラムを認識するツールで開くこと)
l テスト用に、自分でランダムな配列を生成し、長さを 10^4, 10^5, 10^6, 10^7, 10^8 と伸ばしたとき、計算時間が長さに対してほぼ線形に増えてゆくことを確認する。主記憶の大きさを便宜大きくする。
l
実データとして、以下のページから最新のヒトゲノムデータ(2009年2月版 hg19はhuman genomeの19版目を意味する)を入手する。
http://hgdownload.cse.ucsc.edu/goldenPath/hg19/chromosomes/
たとえば
chr1.fa.gz は、1番染色体(chr1) の塩基配列データを fasta 形式(最初の行に “>” からはじまるコメントあり)で格納したファイルで gzip で圧縮してある。 chr *** random.fa.gz および churn**** の名前をもつデータは染色体上での位置が確定していないデータなので、今回は使わないこと。
l
ゲノム中の N および大文字と小文字の扱い方 (2009/10/28 追加)
ゲノム配列中には大文字の A, C, G, T (コード領域)および小文字の a, c, g, t (コード領域以外)が現れるが、今回は区別の必要が無いので、すべて大文字に変換してよい。その他にゲノム配列には、A, C, G, T 以外に N という文字列が現れる。N は塩基を確定できなかったことを意味しており、テロメアやセントロメアのように繰返し配列が連続していて解読することが諦められた領域、それ以外の領域でも繰返し配列の影響で曖昧性の無い解読ができなかった領域、リードエラーで正確な解読ができなかった部分を示している。N のコードの仕方はいくつかある。A, C, G, T に 1, 2, 3, 4 の番号を割り当てて N には 5 を割り当てる方法が自然である。他には N に A, C, G, T を割り当てることで文字を 2 bit 表現して記憶領域を圧縮する方式もあり、主記憶や二次記憶装置が高価で潤沢に使えなかった時代に愛用された。N に A,C,G,Tを割り当てても、その位置に問合せ配列が誤って写像されることを避けるようにできれば悪影響が少ないので、それで我慢することにする。たとえば、ゲノム中には poly-A (A が連続して現れる配列)が多い性質を利用して N に A を割り当てる方式がある。poly-Aを含む問合せ配列は、元からある多数の poly-A に写像されてしまい、位置を確定できない。したがって N に A を割り当てても、問合せ配列の位置を確定できない状況には変化が無いと言える。もうひとつ愛用されているのは N に A,C,G,T をランダムに割り当てる考え方である。問合せ配列がある程度の長さ(長さ 20 以上)の時には、ランダムに生成された場所に誤って当たる確率は非常に低くなるため、間違いは起こりにくい。自分の好きな方式で N を扱ってください。
l 染色体もしくはヒトゲノム全体で suffix array を構築。まずは 21,22番染色体(5000万塩基程度)を対象に実験。オプションとして、ヒトゲノム全体を処理するときには、大きな主記憶が必要なので、workstation (メモリ16GB) を利用してください。suffix array が正しく構築できたことを、induced sorting と Larsson-Sadakane 双方で同じ結果が得られたということで確認できたということをレポートしてくれれば十分です。 ヒトゲノムでの suffix array はその後の演習でも使うので大事です。うまくゆかない場合は、相談してください。(2009/10/28 追加)
l レポート〆切: 10/29(木)(演習期間 3週間)
(2) 問合せ処理とBurrows-Wheeler Transform (BWT) の実装
l suffix array から longest common prefix を計算するプログラムを実装し、長さ k の塩基配列が出現する位置を出力するプログラムを作成する。ヒトゲノム中に出現する回数が非常に多い長い配列を列挙してみる(たとえば k > 100)。
l suffix array から Burrows-Wheeler Transform を構築し、対象とするゲノム配列中に出現する問合せ配列の位置を枚挙するプログラムを作成する。長さn の配列(ヒトゲノム等)に対して、長さ k (= 10〜50) の問合せ配列をランダムに生成して、問合せ時間の平均時間を求め、k (および n) に関して線形比例するかどうか確認する。
l Burrows-Wheeler Transform のデータ構造である Occ の主記憶容量を抑えるために、L 個毎に Occ の行を保存して、on-the-fly で他の行を補完する方法を実装する。L の値を 2, 4, 8, 16, 32, 64, 128, 256, 512, 1024 ぐらいまで変化させたとき、記憶領域の変化と、長さ k (= 20 等) のランダム配列に対する平均計算時間がどのように変化するか観察する。
l
下記の論文を参考にして、ミスマッチおよびギャップを許して問合せ配列をゲノム上にアラインメントするプログラムを作成する。より多くのミスマッチおよびギャップを許容するに従って、速度が変化する様子を分析する。
参考論文:
Heng Li and Richard Durbin. Fast
and accurate short read alignment with Burrows–Wheeler transform.
Bioinformatics 2009 25(14):1754-1760.
なおこのプログラムの動きを例題を使って解説した資料はココ。
2010/01/13 再び更新(最後のスライドのアルゴリズムは何箇所か修正しないと動きません)
l レポート〆切 : 11/19 (木) 延長 2009/12/31 頃まで
(3) Chaining の実装と大規模ゲノムの比較
l 2つのゲノム中 G1, G2で対応する領域を探し出すために、G1 中の長さ k の塩基配列を、BWT を使ってG2中にアラインメントし、その結果を
Chaining するプログラムを作成。k の値や、ミスマッチおよびギャップの許容範囲を工夫。
Chaining Algorithm は生物情報ソフトウエア論の講義ノートおよび配布資料を参照のこと。
l 例としてヒト、マウス、ラットのX染色体の比較を行う。データは以下のサイトからダウンロード。
http://hgdownload.cse.ucsc.edu/downloads.html
l 課題: ヒト、マウス、ラットゲノムのX染色体を比較し、下記の参考論文中の図6a「X染色体の比較図」に示す図を描画する。
参考論文: Rat Genome Sequencing Project Consortium. Genome
sequence of the Brown Norway rat yields insights into mammalian evolution.
Nature 428, 493-521 (1 April 2004)
図6 X染色体の比較
l 担当: 中谷洋一郎助教
l レポート〆切: 延長 2010/1/31
(4) ゲノムブラウザーとデータベース (オプション課題)
l ゲノムブラウザーとデータベースシステムの利用
http://www.xerial.org/wiki/lecture/2009/Report4
l 担当: 斉藤太郎助教
l 演習+講義: 2010/1/13 (水) 演習資料を説明し配布
l レポート〆切: 延長 2010/2末日
下記の課題は4年生夏学期の演習に回します
転写開始点、クロマチン構造、DNAメチル化 (書きかけで将来更新)
l
次世代シークエンサーを利用した転写開始点、クロマチン構造、DNAメチル化情報の収集
参考文献: 蛋白質 核酸 酵素 2009年8月号 特集 次世代高速シークエンサーの応用と情報解析
l データソース TBA
l 「転写開始点は1点でなく分布する」 「ヌクレオソーム位置も一定していない」 などを確認