allele count in bam files

Question

[le.marcorin]

Hi! I'm trying to get the allele count for some snps in my bam files.

The variant calling was done by haplotypeCaller. The allele count in the output doesn't represent the real number of reads I have in my bam file, whitch is what I want.

I tried with bcftools mpileup, but the output only shows genotype liklyhoods. Also tried bam-readcount but it deosn't work with peired-end reads (I think).

I need to know the number of reads supports each allele in my vcf file.

Does anyone know how I coud do it?

[h.mon]

Could you post one vcf line from a SNP you are interested?

Answer 1

[h.mon]

I need to know the number of reads supports each allele in my vcf file.

If you are using HaplotypeCaller, the AD field is exactly what you want. See What is a VCF and how should I interpret it?, which explains it, and also briefly explains why AD may be different from what you expect.

AD is the unfiltered allele depth, i.e. the number of reads that support each of the reported alleles. All reads at the position (including reads that did not pass the variant caller’s filters) are included in this number, except reads that were considered uninformative. Reads are considered uninformative when they do not provide enough statistical evidence to support one allele over another.

See also Allele Depth (AD) is lower than expected, which enters into more details about AD quirks.