Calling variants on reads with MAPQ=0 on HaplotypeCaller or bcftools mpileup
3
0
Entering edit mode
3.1 years ago

I am working with about 500 samples of human exome data. used hg19 to align my reads and ran a standard best-practices GATK workflow. Later only to realise that a small 1Mb loci has not mapped properly due to the presence of a haplotype specific contig in the hg19 FASTA.

How do I force HaplotypeCaller or BCFTools to call a variant with mapping quality 0? I have disabled the read filters but that does not work. I just want the variants from one particular gene and since I have 500 samples, re-aligning them to a modified reference is not an option at the moment due to time constraints. I am aware of false positives occurring but we will validate any variants seen.

This is the command I've tried on BCFtools:

bcftools mpileup -q 0 -f ucsc.hg19.fasta file.bam -r chr17:43000000-45000000 | bcftools call -mv > file.vcf

this only outputs variants with MQ =/= 0

whole-exome-sequencing variants GATK hg19 • 1.4k views
ADD COMMENT
0
Entering edit mode
3.1 years ago

try

bcftools mpileup -q -1  -f ucsc.hg19.fasta file.bam -r chr17:43000000-45000000 | bcftools call -mv > file.vcf
ADD COMMENT
0
Entering edit mode

Thank you for responding! Unfortunately, this command returns the same result as the one I have posted above.

I've checked it on IGV and the variants are definitely present.

Can you think of any other solution to this problem, perhaps?

ADD REPLY
0
Entering edit mode
2.5 years ago
YT • 0

Hi,

-q 0 is the default value. Maybe try bcftools mpileup -A -f reference.fa Try with the -A option. Hope it works for you.

Best wishes, Yutang

ADD COMMENT
0
Entering edit mode
21 months ago
Gama313 ▴ 120

As far as I know, a MAPQ equal to 0 identifies reads mapping to multiple positions (i.e. pseudogenes). This means that you cannot call that mutation with enough confidence. If the region of interest is filled with mutations, it is definetly misaligned. Try blast sequences supporting your mutation and see where they come from. Let me know how it goes.

ADD COMMENT

Login before adding your answer.

Traffic: 1500 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6