This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Allele frequency of vcf file

I know it is a stupid question, but I'm a beginner of snp calling and really confused.

I run samtools and bcftools get vcf files like this

#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  accepted_hits.bam
1       565286  .       C       T       222     .       DP=27;VDB=2.078631e-01;AF1=1;AC1=2;DP4=0,0,14,12;MQ=40;FQ=-105  GT:PL:GQ        1/1:255,78,0:99

However, I not only want AF1 (Max-likelihood estimate of the first ALT allele frequency), but also want AF (allele frequency). How can I get a result containing both information?

I run samtools and bcftools like this:

samtools mpileup -ugf /data/Genomes/Homo_sapiens/Homo_sapiens_assembly19.fasta accepted_hits.bam | bcftools view -bvcg - > SRR1294493.bcf
bcftools view SRR1294493.bcf | vcfutils.pl varFilter -D 44 > SRR1295163.vcf

Thanks a lot.

snp
• 10,731 views • link
updated by steve steve • written by juncheng juncheng

seems you are missing RO and AO tags, reference observation count and alternate observation count.

• 0 views • link
written by Adrian Pelin Adrian Pelin

Yes, it is. Do you know why is like this?

• 0 views • link
written by juncheng juncheng

Not sure. Try to use FreeBayes, I know it reports more than enough information.

• 0 views • link
written by Adrian Pelin Adrian Pelin

4 answers

export BCFTOOLS_PLUGINS="$HOME/src/bcftools/plugins/"

bcftools +fill-AN-AC SRR1294493.bcf

• 1 views • link
written by Ryan Layer • Ryan Layer

can you explain what this is? Not clear what is going on in this answer. Thanks

• 0 views • link
written by steve steve

If by allele frequency you mean the fraction of reads supporting alternate allele, that information is provided in DP4 field.

DP4=0,0,14,12;

DP4 reports the number of reads covering the position with the reference allele mapped to forward and reverse strands, followed by alternate allele mapped to forward and reverse strands.

Allele frequency here would be 0 (0/26) for reference allele, and 1.0 (26/26) for alternate allele.

• 1 views • link
written by Noushin N Noushin N

right, but its more useful in most cases to have the pre-computed value saved and available for use

• 0 views • link
written by steve steve

If anyone knows anything, please help. I'm really struggle at the moment.

• 0 views • link
written by juncheng juncheng

Hi juncheng,

you can use bcftools plugins to fill in the missing AF tag to your output vcf file.

The following example could be useful to you :

export BCFTOOLS_PLUGINS="~/mytools/bcftools/bcftools-1.3.1/plugins/"

$samtools mpileup -t DP,AD,ADF,ADR,SP,INFO/AD,INFO/ADF,INFO/ADR -Bugf ~/reference/human_g1k_v37.fasta -l target_regions.bed na12878_AlignedReads.bam 2> mpileup_errorLog.log | \ bcftools call -vmO u -f GQ 2> bcftoolsCall_errorLog.log | \ bcftools plugin fill-tags -Oz -o output_variants.vcf -- -t AF 2> bcftoolsPlugin_errorLog

• 0 views • link
written by gsr9999 gsr9999

Log in to answer this question.

More posts like this