samtools depth print out all positions I am trying to use samtools depth (v1.4) with the -a option and a bed file listing the human chromosomes chr1-chr22, chrX, chrY, and chrM to print out the coverage at every position:

cat GRCh38.karyo.bed | awk '{print $3}' | datamash sum 1
3088286401


I would like to know how to run samtools depth so that it produces 3,088,286,401 entries when run against a GRCh38 bam file:

samtools depth -b $bedfile -a $inputfile


I tried it for a few bam files that were aligned the same way, and I get differing number of entries:

3087003274
3087005666
3087007158
3087009435
3087009439
3087009621
3087009818
3087010065
3087010408
3087010477
3087010481
3087012115
3087013147
3087013186
3087013500
3087149616


Is there a special flag in samtools depth so that it reports all entries from the bed file?

If samtools depth is not the best tool for this, what would be the equivalent with sambamba depth base?

sambamba depth base --min-coverage=0 --regions $bedfile $inputfile


Any other options?
You might try using bedtools genomecov instead. If you provide the -d option, it reports the coverage at every position in the BAM file.

bedtools genomecov -d -ibam $inputfile &gt; "${inputfile}.genomecov"


You can also provide a BED file if you just want to calculate in the target  region.