How to calculate bam coverage apart from Qualimap?

Question

[marongiu.luigi]

Hello,

I have been using qualimap to generate stats of the bam files I got, including coverage. Recently, I am getting crashes because there is not enough memory, even if I use the memory parameter --java-mem-size=20000M.

Is there another way to calculate the mean coverage of a bam file? I have tried with bedtools genomecov -bga FILE.bam > FILE.coverage.bed but I got a usage error (ERROR: Unrecognized parameter: the input file and ERROR: Need both a BED and a genome file). I am trying with samtools depth -a file.bam | awk '{c++; if($3>0) total+=1}END{print (total/c)*100}'. Is this one good or is there a fastest way?

Thanks

[Pierre Lindenbaum]

mosdepth .

Answer 1

[ATpoint]

The tool does not produce (the expected) output because you are not using the correct syntax/options.

genomecov expects the input via -ibam, not just the plain file name.

Just type bedtools genomecov and the helps pops up.

For alternatives please use the search function, start here:

Tools To Calculate Average Coverage For A Bam File?

[marongiu.luigi]

bedtools genomecov -ibam works, thanks, but it gives coverage for each position. can I get a total average coverage?

[ATpoint]

What do you mean by total coverage? Something like this, simply based on the % of mapped reads including the read length? How To Calculate Coverage

[marongiu.luigi]

Yes, Qualimap gives me, for instance, Coverage mean: 28.4187 plus the std dev and other stats. I only need the coverage mean. But the file I have is too big (300 Gb) and eats all the memory causing crashing. The methods listed here and in the other post are not good. For instance, samtools depth *bamfile* | awk '{sum+=$3; sumsq+=$3*$3} END { print "Average = ",sum/NR; print "Stdev = ",sqrt(sumsq/NR - (sum/NR)**2)}' return an empty line. I need a work around...