Wet系バイオインフォマティシャンの災難

機械に弱いバイオ研究者に降りかかる災難を綴っていきます。 Twitter: @1wantphd

【mapping】BWAでマッピング・出力したファイルはsam形式ではない(場合がある)

【要約】

・ロングリードのマッピングbwa memではsam形式のファイルが出力される

・ショートリードのマッピングbwa alnでは"sai形式"のファイルが出力される

そのため、samファイルに書き換えるためにはbwa samseをする必要がある

 

以下詳細

【BWAを使い始めた理由】

・Bowtie2が遅すぎる

・Hisat2の爆速マッピングに慣れすぎてBowtie2のスピードは止まって見える

 

【BWAとは】

ロングリードを高速でマッピングするために開発されたツール

メモリの消費量が数GBで済む

もちろんショートリードにも対応している

 

【いざ実践】

いつもお世話になっているhttps://bi.biopapyrus.jp/様を参考にマッピングに挑戦

https://bi.biopapyrus.jp/rnaseq/mapping/bwa/

 

【インデックスの作成】

bwa index -p [index_file名] UCSC_hg38.fa

hg38のヒトゲノムでインデックスを作成

Mac book proなら30分ほどで終わる

 

マッピング:問題発生】

ショートリード(70 bp以下)なのでalnのオプションをつける

bwa aln -t 3 [index_file名] DRRXXXXX.fastq > DRRXXXXX.sam

 

出力されたファイルを確認すると…

head DRRXXXXX.sam 
?92?92?9\6s}\6s}??A???A?????\?U\?U?????????m??m??3A?3AA?XA?X??Q??Q?BR?BR?A?X?A?X?@?X?@?XA?XA?XA?XA?X??z???z??)?)N?1@?zT?ZzT?Z?R?_?R?_???]??]{T?Z~T?ZO_O_??v??v?w?]?w?]??]??]?w?]?w?]Fw$

明らかにバイナリ形式でsamファイルではない

(sam ファイルは人の読める形式のファイルであるhttps://bi.biopapyrus.jp/format/sam.html

 

Githubのbwa readmeを確認すると(https://github.com/lh3/bwa/blob/master/README.md)

イルミナの70 bp以下のシングルエンドリードをbwa memを使ってマッピングした場合、saiファイルとして出力される…だと…?(なんてこった…)

saiファイルはbwa samseを使ってsamに変換できるので

以下のコマンドで目的とするsamファイルを得ることができました。

bwa aln -t 3 [index_file名] DRRXXXXX.fastq > DRRXXXXX.sai
bwa samse [index_file名] DRRXXXXX.sai DRRXXXXX.fastq > DRRXXXXX.sam

※注意:各ファイル名間のスペースは”1つ”でなければならない

 

 【感想】

10Gb(約4000万リード)のマッピングをsamファイルにするまでの所要時間は

30分から1時間程度と十二分な速さを見せつけてくれました。BWAに感謝。

同じデータをbowtie2でマッピングしようとしたら20時間経っても終わりませんでした。(これはこれで何かおかしい? マウスのマッピングなら2,3時間でいつも終わるのに)

マッピングの精度とかについてはよくわからないですが素晴らしいMapperであることに間違いなさそうですので是非使ってみてください。

 

【教訓】

ツールを使う際は開発元の説明をちゃんと読むべし