This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Bam-readcounts min-base-quality not giving expected output

I would like to limit the base-quality-score to phred 13 (similar to default of samtools mpileup (1.3.1)), using bam-readcount (0.8.0). Despite various command modifications, I'm not getting the expected output.

Default:

bam-readcount -f human_g1k_v37.hotspot.fasta final.bam -l F5DNA.txt --max-warnings 0

1       169519049       C       8085    =:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00  A:21:60.00:20.05:0.00:15:6:0.20:0.01:33.00:15:0.85:124.00:0.62  C:39:60.00:17.72:1.54:24:15:0.18:0.01:21.10:24:0.85:123.92:0.54 G:8:60.00:11.38:0.00:1:7:0.07:0.02:37.50:1:0.85:124.00:0.13     T:8014:60.00:26.69:0.04:4033:3981:0.15:0.01:35.90:4033:0.85:123.97:0.44 N:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00  -C:3:60.00:0.00:0.00:2:1:0.18:0.01:9.33:2:0.86:124.00:0.58

min-base-quality 13:

bam-readcount -f human_g1k_v37.hotspot.fasta final.bam -l F5DNA.txt --min-base-quality 13 --max-warnings 0

1       169519049       C       8067    =:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00  A:21:60.00:20.05:0.00:15:6:0.20:0.01:33.00:15:0.85:124.00:0.62  C:27:60.00:20.70:2.22:14:13:0.16:0.01:21.04:14:0.85:123.96:0.45 G:2:60.00:24.50:0.00:1:1:0.15:0.01:24.50:1:0.85:124.00:0.44     T:8014:60.00:26.69:0.04:4033:3981:0.15:0.01:35.90:4033:0.85:123.97:0.44 N:0:0.00:0.00:0.00:0:0:0.00:0.00:0.00:0:0.00:0.00:0.00  -C:3:60.00:0.00:0.00:2:1:0.18:0.01:9.33:2:0.86:124.00:0.58

Mpileup: (-Q0, not showing the variants)

 samtools mpileup -f human_g1k_v37.hotspot.fasta -d 1000000 -A -l F5DNA.txt -Q0 final.bam 
1       169519049       C       8085

Mpileup: (default of -Q13, not showing the variants

samtools mpileup -f human_g1k_v37.hotspot.fasta -d 1000000 -A -l F5DNA.txt -Q13 final.bam 
1       169519049       C       4072

Looking at IGV, the readcounts at this position are 8085, which matches both mpileup and bam-readcount output when minimum base quality is 0. Any idea why bam-readcounts does not appear to filter the same when the minimum base quality is set to 13?

bam-readcount mpileup

1 answer

answered here: https://github.com/genome/bam-readcount/issues/37

Log in to answer this question.