bedtoolsでbamファイルのread coverageを求める方法は以前にも書いた。
現在のbedtoolsのバージョン(2.27.1)ではbedファイルではなくbamファイルから直接read coverageを出すことが可能になっている
$ bedtools coverage -a referense.bed -b sample.bam > sample_cov.bed
というフォーマットで問題なく解析できるのだが。
ただしこの操作は非常にメモリを食うらしく、4GB超のbamファイルを解析しようとすると64GBのメモリを使い切ってスワップも15GBくらい使って、どうにもならなくなって落ちる。
ところが、
$ bedtools bamtobed -i sample.bam > sample.bed
と一旦bedに変換してから、先程のコマンドにつなげると、bedへの変換にはほとんどメモリを消費せず、coverageの計算も30GB程度のメモリで十分に収まることがわかった。
というわけで、bam直の解析はやめといたほうが良い。
