How Can I Get The Field Of Ac Key-Value Pair In The Info Column In Vcf File?

Question

[cli2]

Dear all,

I am new to DNA sequence analysis, especially variant calling. I had some sequence data and used bwa to perform the alignment. After I had sam file from bwa, I use samtools and bcftools to perform variant calling. However, when I more the final vcf file, I don't see the AC key-value pairs in the vcf file. I cannot calculate the frequency of alleles at each locus. Would anyone have an idea about what I should do?

My purpose is to find the allele frequencies at each locus. If there are other methods to do this, please give me a help. Thank you.

In fact, when I view the sorted bam file, I can see the alignment and the reference. It looks fine.

---

I tried the following commands in my workflow:

samtools mpileup -g -f ref.fa short_reads_aligned.sorted.bam > short_reads_variants.bcf

bcftools view -c -v short_reads_variants.bcf > short_reads_variants.vcf

---

I also tried the following

samtools mpileup -uf ref.fa short_reads_aligned.sorted.bam | bcftools view -bvcg - > var.raw.bcf bcftools view var.raw.bcf > var.flt.vcf

Thanks in advance C

[cli2]

Dear pd3,

I don't know how to use this "vcf-annotate" command. Could you please provide more details? Thank you.

Answer 1

[pd3]

One of many ways to do it is

vcf-annotate --fill-AC-AN

from http://vcftools.sourceforge.net/perl_module.html#vcf-annotate.

Newer versions of samtools+bcftools (develop branch on github) do output AC, AN natively.

Answer 2

[cli2]

I finally use GATK and obtain the AC field in my vcf file. I think the problem is due to that the mpileup does not output the allele counts for each locus in the sequence although mpileup used allele counts in the process of generating bcf and vcf files. So the vcf-annotation can add AC field in the vcf but it is not the correct one.