Maqの使い方まとめ
自分用まとめです
レファレンスのシークエンスをとってくる
Humanのレファレンス配列をUCSCからとることにする
UCSC Genome Browser: DownloadsのHuman Genome (hg19, GRCh37)
Data set by chromosomeにある
染色体毎の配列(chr1.fa.gzなど)をダウンロード
解凍して、一つのファイルにまとめておく
cat chr*.fa > allchr.fa
シークエンサーのデータをとってくる
ウチには次世代シークエンサーは無いので(^^;)、NCBIのSRA(Short Read Archive)からとってくる
今回の目的に合いそうなものとして
Study: Targeted capture and massively parallel sequencing of human exomes JS0001 (SRP000910)
のデータを使うことにする
そのうちの一つ
Human HapMap individual NA19240
のデータ(SRX005930)のSRR017992をダウンロードしてみた。解凍してみたら10GB近くあった(@_@)
ところで、
SRAのFastqデータをMaqで使うには、ちょっと加工してやらなければならないらしい
タイトル行にあるスペースを削除すればいいそうな GCOEのテキストに習ってsedで加工してやる
sed 's/ //g' < SRR017992.fastq > SRR017992.fastq_ns
MAQでマッピング
いろいろ細かいことは後でマニュアルを読んで調べることとして、まずはeary runでマッピングしてみた
$ maq.pl earyruy -1 60 -d SRR017992 allchr.fa SRR017992.fasta_ns
-1 60 のオプションは、Fastq配列の頭から60ベースのみをマッピングに使うためにつけた*1
後は待つこと半日〜1日
ウチのマシン(2x2.26 GHz Quad Core Xeon, 8GB)で36,394,026spotsのデータをマッピングするのに、約24時間かかった
マルチコアに最適化していないと作者もいっているとおりで、シングルコアしか使っていなかった。
MaqViewrでデータを眺める
まずはインデックスをつけてやる必要がある
$ maqindex -i -c consensus.cns all.map
そしたらMaqviewでデータを眺めてみる
$ maqview -c consensus.cns all.map