Entering edit mode
7.8 years ago
sheenams
▴
50
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?
