hi I have used following command to calculate depth of target regions in my bed file
samtools depth -ab target_interval.bed 12345XX.bam > sample_depth.txt
And I used following command to see how many reads are mapping to particular region
samtools view -c 12345XX.bam chr1:111212724-111212724
This command gave output as : 4
I checked the same region region in the sample_depth.txt. But the output there is 2
What is the difference between the above methods? Why does it give two different outputs?
1 answer
Samtools depth is using the mpileup algorithm to find overlapping data, along with all the nuances that involves. That means filtering by flags (unmapped data, secondary reads, duplicates, QC failure), limits of maximum depth, possibly some other things like removal of overlapping templates (I cannot recall if that is on or off by default), and definitions of what actually means aligned (CIGAR N operation isn't for example; I forget if D is). Note unmapped data here will mean placed (having CHR/POS) but not aligned.
View however is just an unfiltered count unless you also modify the filtering options, and even then it will be counting records spanning the region with N (and perhaps D) cigar ops.
Which is appropriate for you depends on the question you are asking.
Log in to answer this question.