How to count the number of alignments in a bam file without using samtools?

Question

[steve]

Preferably, using either standard GNU tools, or perhaps something in the Python standard library. Any ideas?

[Carlo Yague]

if the bam is uncompressed, you could juste use wc -l which will return the number of lines (= the number of alignments assuming there is no unaligned entry in the file). But samtools is much better.

edit: I got it wrong, see comments below.

[steve]

dont think I have ever seen a case where people actually use uncompressed bam file, at that point you might as well be using sam format

[Carlo Yague]

It is used when one need to pipe a bam file in another process (saves compression time, as said below by Jorge).

[swbarnes2]

I think by definition, a bam is compressed.

[Carlo Yague]

not with samtools view -bu I think

[Jorge Amigo]

Even if uncompressed, which is recommended for piping purposes in order to save time compressing and decompressing output and input respectively, a bam file is still binary, so wc -l will not work.

[Carlo Yague]

You are both right. I tested it to make sure:

samtools view -bu input.bam | wc -l
459910
samtools flagstat  input.bam | head -n 1
27053366 + 0 in total (QC-passed reads + QC-failed reads)

[Jorge Amigo]

Just a couple of comments:

-u implies -b, therefore samtools view -u input.bam is enough to get an uncompressed bam file.

Also, if you only need to know the number of reads, generating all flagstat metrics is not that efficient, and samtools view -c input.bam would be sufficient. It won't be that faster, since most of the time is spent in decompressing and reading the bam file, but it neither will be slower plus it's simpler to write down.

Answer 1

[xiaoguang]

Bam file is a binary file stroing alignment information, so we must need htslib.h api to interprete that. if you use GNU, you need htslib.h. if you use python, pysam must be the best choice for you.