Allele frequency of vcf file
4
1
Entering edit mode
8.0 years ago
juncheng ▴ 200

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 • 7.0k views
ADD COMMENT
0
Entering edit mode

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

ADD REPLY
0
Entering edit mode

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

ADD REPLY
1
Entering edit mode

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

ADD REPLY
2
Entering edit mode
5.9 years ago
Ryan Layer ▴ 60

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

bcftools +fill-AN-AC SRR1294493.bcf

ADD COMMENT
0
Entering edit mode

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

ADD REPLY
2
Entering edit mode
5.9 years ago
Noushin N ▴ 600

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.

ADD COMMENT
0
Entering edit mode

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

ADD REPLY
0
Entering edit mode
8.0 years ago
juncheng ▴ 200

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

ADD COMMENT
0
Entering edit mode
5.3 years ago
gsr9999 ▴ 230

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

ADD COMMENT

Login before adding your answer.

Traffic: 1835 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