bcftools
bcftools copied to clipboard
samtools
Heterozygous SNP sites calling
Hello,
I am currently using the mpileup | bcftools call pipeline to call multiple SNPs in diploid organisms. If I understood this correctly, the call command exploits the genotype likelihoods generated by mpileup (based on a summary of mapped reads on the reference genome) to make the actual calling. I was wondering if there is a way to tweak bcftools to assess the zigosity of SNPs based on a customised AD ratio, such as that all SNPs where ad1/ad2 or ad2/ad1 is >0.2 are called heterozygous. Where ad1 = REF depth and ad2 = ALT depth. Here is an example: #CHROM POS ID REF ALT QUAL FILTER INFO FORMAT sample_sorted.bam chr_1 131328013 . A G 950 . 0 GT:PL:AD:GQ 1/1:168,106,0:779,3193:99
In the example above AD = 779 , 3193 so the fraction between 779 (ad1) and 3193 (ad2) is ca 0.244. In this case I would like to tell bcftools to call this SNP site heterozygous simply based on my criterion, instead it calls it homozygous. Is there a way around this? Thank you.
Can you please show also the raw mpileup line? The calling is informed also by the QS tag.
Hi, sorry for the delay in replying.
I am using the older version samtools mpileup | bcftools call . Would that generate the QS tag? I cannot see it in my VCF.
You'll need a much newer version of bcftools (ideally the latest) and use bcftools mpileup rather than samtools mpileup. By the way, it always helps to show the full command.