SAMtools
先日紹介したリードのマッピングの記事では、出力をSAM形式に指定していました。
rnakato.hatenablog.jp
SAM形式はファイルサイズが非常に大きく読み込みにも時間がかかるので、バイナリ形式のSAMであるBAM形式が下流解析でよく利用されています。
ここではSAM形式をSAMtoolsを使ってBAM形式に変換します。
SAMtoolsのインストール
Ubuntuであればaptでインストールできます（ややバージョンが古い）。
$ sudo apt install samtools
Anacondaを使っている人であればcondaでもインストールできます。
$ conda install -c bioconda samtools
どちらでもない人はSAMtoolsのWebサイトからダウンロードしましょう。
BAM作成
以下のコマンドを実行すると、SAM形式のファイルを読み込んでBAM形式に出力します。
$ samtools view -bS sample.sam > sample.bam
-b はBAM系形式で出力するオプションで、-Sは入力ファイルの形式を自動認識させるオプションです。-hをつけると、SAMファイルのヘッダを含めるようになります（デフォルトでは省略される）。
-@オプションをつけると使用CPU数を指定できます。
$ samtools view -bS sample.sam -@ 4 > sample.bam
BAMのソート
上記コマンドで得られたbamファイルはソートされていませんが、多くのツールではsorted BAMを要求されます。
以下のようにしてbam をソートしましょう。
$ samtools sort sample.bam > sample.sort.bam
これでソートされたbamファイルが生成されました。
ソート前のファイルは不必要なので削除しましょう。
注意点：古いバージョンのsamtools (0.1.19など)とはsortコマンドの指定方法が変わっていますので、
古いバージョンを想定したスクリプトを新しいバージョンのSAMtoolsで動かそうとするとエラーになります。
sortコマンドがエラーになる時はバージョンを意識してみましょう。
リダイレクトによるsorted bam生成
上記のコマンドはBAM作成 -> sorted BAM作成を別々に行っており、中間ファイルを削除する手間もあるのが面倒です。
そこで、以下のようにリダイレクトで生成されるBAMファイルを直接ソートしてみましょう。
$ samtools view -bS sample.sam -@ 4 | samtools sort > sample.sort.bam
| は、左側のコマンドから標準出力に出力される結果をそのまま右側のコマンドの入力にするという命令です。
samtools viewの出力を直接samtools sortに投げることで、中間ファイルを作成することなくsorted bamを作成することができます。
なお、以下のようにsortコマンドに直接SAMファイルを指定しても同じことができるようです。
$ samtools sort -@ 4 sample.sam > sample.sort.bam
注意点：samtoolsは内部で生成される中間ファイルをカレントディレクトリに生成しますので、
一度の多量のファイルを並列で変換すると、中間ファイルが衝突して結果がおかしくなる可能性があります。
BAMのindex作成
一部のツールはBAMの読み込みを高速化するためにBAMのindex（索引）ファイルを要求します。
以下のコマンドでindexを作成します。なお、indexコマンドはsort済のBAMファイルしか受け付けません。
$ samtools index sample.sort.bam -@ 4
こうすると、sample.sort.bam.bai という名前のindexファイルが生成されます。
BAMファイルを作成する時は、sortとindex作成を常に行うようにしておけば間違いないでしょう。
おまけ：マッピングツールの出力をそのままsorted bamに変換する
マッピングツールの出力をリダイレクトでsamtoolsに投げてやると、そもそもSAMファイルを作成・削除する手間が要りません。
$ bowtie2 -x index sample.fastq -p4 | samtools view -bS - | samtools sort > sample.sort.bam
for文でたくさんのサンプルを処理する時も１行で書けて楽ちんです。
$ for num in $(seq 1 20); do
$ bowtie2 -x index sample$num.fastq -p4 | samtools view -bS - | samtools sort > sample$num.sort.bam
$ done
