Skip to content

Wrong number of bases/reads reported by bedcov in very high coverage data #1950

Description

@ellisjj

I'm running samtools v1.18 on macOS/arm64

I have a dataset of very high coverage data, when running

samtools bedcov -c -g SECONDARY,QCFAIL,DUP $BED $BAM

it reports the number of reads in one interval as 1,300,694. However,

samtools view -c -F UNMAP $BAM $REGION

reports 6,277,273 which is significantly different. bedtools coverage also reports 6,277,273. Through trial and error, I discovered that

samtools bedcov -c -d 1000000 -g SECONDARY,QCFAIL,DUP $BED $BAM

reports 6,277,273. Different values for -d results in a number of reads reported between the 1.3M and 6.2M. My question is, is this a bug in bedcov when processing very high coverage data or a feature that's limiting the number of reads under some circumstances? It's not obvious to me at all why the -d option would change the total number of reads reported, and, if this isn't a bug, how can I ensure that the total number of reads is always reported (-d 1000000 worked for one particular BAM but might not universally work)?

Thanks

Metadata

Metadata

Assignees

Labels

No labels
No labels

Type

No type

Fields

No fields configured for issues without a type.

Projects

No projects

Milestone

No milestone

Relationships

None yet

Development

No branches or pull requests

Issue actions