<p>It sounds like your main confusion is how to calculate AF. I'll answer that right away, but I'll also suggest a Python library that you can use to do this and much more.</p>
<p>You're correct, the AF is defined as the number of non-missing alternate alleles (AC) divided the number of non-missing alleles (AN). Calculating these are straightforward. Map each GT to its alternate alleles:</p>
<pre><code>0/0 =&gt; 0
0/1 =&gt; 1
1/1 =&gt; 2
./0 =&gt; 0
./1 =&gt; 1
./. =&gt; 0
</code></pre>
<p>Now, compute the sum. That is the AC. For the AN, map each GT to its number of non-missing alleles:</p>
<pre><code>0/0 =&gt; 2
0/1 =&gt; 2
1/1 =&gt; 2
./0 =&gt; 1
./1 =&gt; 1
./. =&gt; 0
</code></pre>
<p>Now, compute the sum. That is the AN. Now calculate AF as</p>
<pre><code>AC / AN
</code></pre>
<p>All of this is kind of annoying to do though! I suggest you use <a href="https://hail.is" rel="nofollow noreferrer">Hail</a> (disclaimer: I am a Hail developer).</p>
<pre class="lang-py prettyprint-override"><code>import hail as hl

# Hail is not yet aware of Plasmodium falciparum reference genome
mt = hl.import_vcf('sample.vcf', reference_genome=None)

# let's peek at our data
mt.show(n_rows=10, n_cols=5)

# Let's get some information about the GT specifically:
mt.GT.summarize()

# We can use Hail's built in variant_qc method to calculate AF for us:
mt = hl.variant_qc(mt)

# Let's take a look at that
mt.variant_qc.show()

# We can overwrite the old AF field and export back to a VCF file if necessary
mt = mt.annotate_rows(
    # AF[1] is the frequency of the first alternate allele
    # Your dataset appears to have some multi-allelics. The other
    # alternate allele frequencies are available in AF[2], ...
    info = mt.info.annotate(
        AF = mt.variant_qc.AF[1]
    )
)
hl.export_vcf(mt, 'sample-updated.vcf.bgz')
</code></pre>
<hr />
<p>This comment doesn't make sense to me:</p>
<blockquote>
<p>I learned that this calculation is not suitable for mithocondrial data.</p>
</blockquote>
<p>The definition of AF, AC, and AN doesn't change depending on the type of config (autosome, sex chromosome, mitochondrial). Perhaps the issue is ploidy encoding? A GT in a VCF can be:</p>
<ul>
<li>1-ploidy, e.g. <code>0</code> (reference), <code>1</code> (first alternate), <code>2</code> (second alternate)</li>
<li>2-ploidy, e.g. <code>0/0</code> (homozygous reference), <code>0/1</code> (heterozygous), <code>1/1</code> (homozygous first alternate)</li>
<li>etc.</li>
</ul>
<p>Sometimes people incorrectly encode a 1-ploidy config (for example the Y chromosome) as: <code>0/0</code> and <code>1/1</code> instead of <code>0</code> and <code>1</code>. That's not your fault though and it doesn't mean the definition of AF, AC, and AN change. It means someone gave you malformed data.</p>

CC BY-SA 4.0
2023-01-30T17:04:35.670Z
Thanks for the reply. You totally understood my question. Reagrding the part you said that does not make sense, what i meant was that cells can have multiple mts. From the info i gathered so far, due to this reason the simple equation of AF = AC / AN would not be appropriate for mt genome meaning using GT alone will not return the right result. The calculations used in softwares such as vcftools are only for non mt genome. I was trying to do it in Python anyway. I am glad you introduced this library. Does it work on Windows too?
eh329
CC BY-SA 4.0
2023-01-30T21:09:10.067Z
We don't regularly test against windows, but, in principle, it should work. Let me know if you run into bugs and we'll fix them!
Daniel King
CC BY-SA 4.0
2023-01-30T15:36:57.197Z
2023-01-30T15:36:57.197Z
Daniel King
<p>I have a VCF file with 200 samples (mitochondrial genome of Plasmodium falciparum). Here is a pic to take a look at:</p>
<p><a href="https://i.sstatic.net/M3wuG.jpg" rel="nofollow noreferrer"><img src="https://i.sstatic.net/M3wuG.jpg" alt="enter image description here" /></a></p>
<p>And a few relevant lines from the actual file:</p>
<pre><code>##INFO=&lt;ID=AC,Number=A,Type=Integer,Description=&quot;Allele count in genotypes, for each ALT allele, in the same order as listed&quot;&gt;
##INFO=&lt;ID=AF,Number=A,Type=Float,Description=&quot;Allele Frequency, for each ALT allele, in the same order as listed&quot;&gt;
##INFO=&lt;ID=AN,Number=1,Type=Integer,Description=&quot;Total number of alleles in called genotypes&quot;&gt;
#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  FP0008-C   FP0009-C
Pf_M76611       19      .       A       T       33196.8 MissingVQSLOD;Mitochondrion     AC=9;AF=0.0003019;AN=29816;ANN=T|intragenic_variant|MODIFIER|PF3D7_MIT00100|PF3D7_MIT00100|gene_variant|PF3D7_MIT00100|||n.19A&gt;T||||||;AS_InbreedingCoeff=0.0437;AS_QD=20.96;BaseQRankSum=0.141;CDS;DP=559138;ExcessHet=0.0026;FS=2.892;InbreedingCoeff=0.52;MLEAC=8;MLEAF=0.0002683;MQ=60;MQRankSum=0;QD=20.96;ReadPosRankSum=0.856;RegionType=Mitochondrion;SOR=0.789;set=genotypegvcfs  GT:AD:DP:GQ:PL  0/0:52,0:52:99:0,120,1800  0/0:54,0:54:99:0,120,1800
</code></pre>
<p>A small subset of samples I am working on can be found here:
<a href="https://www.mediafire.com/file/6dkct6guk5zx83m/sample.vcf/file" rel="nofollow noreferrer">https://www.mediafire.com/file/6dkct6guk5zx83m/sample.vcf/file</a></p>
<p>I am trying to calculate allele frequency (AF) for each variant in the dataset. The population is assumed to be all the samples in the file.</p>
<p>In other words, even though the VCF file contains AF statistics from population data, I would like to know how to recalculate AF only from the VCF file. I believe I need to know how AC and AN are calculated first.</p>
<p>AF is calculated this way: allele count (AC) / allele number (AN)</p>
<p>I tried to calculate them from the information from GT (counting 0/0, 1/1 and 0/1 in all the 200 samples), but then I learned that this calculation is not suitable for mithocondrial data.</p>
<p>All lines in the VCF file contain an AF statistic:</p>
<pre><code>$ grep -v 'AF=' sample.vcf | grep -v '^#'
&lt;No output&gt;
</code></pre>
<p>The information is provided under INFO column, but how were they calculated? If I do not have AC and AN to begin with, what should I do?
Is there a way to calculate AC and AN from the other information provided in the dataset? such as the data mentioned in the FORMAT column?
GT: genotype
AD: unfiltered allele depth
DP : read depth at this position for this sample (Integer)
GQ : conditional genotype quality, encoded as a phred quality
PL : the phred-scaled genotype likelihoods rounded to the closest integer</p>
<p>Thanks</p>

CC BY-SA 4.0
2023-01-27T16:58:25.557Z
Can you give us some more context? What is AC? Allele count in the population? From GnomAD, maybe? What tool produced this file?
terdon
CC BY-SA 4.0
2023-01-27T17:56:33.070Z
@terdon Thanks for the reply. I believe that is what AC means. I have not seen any other definition.
eh329
<vcf><allele-frequency>
How to Recalculate Allele Frequency from a VCF File?
# How to Recalculate Allele Frequency from a VCF File?

(Asked by: username_0 on 2023-01-26)

I have a VCF file with 200 samples (mitochondrial genome of Plasmodium falciparum). Here is a pic to take a look at:

And a few relevant lines from the actual file:
##INFO=
##INFO=
##INFO=
#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  FP0008-C   FP0009-C
Pf_M76611       19      .       A       T       33196.8 MissingVQSLOD;Mitochondrion     AC=9;AF=0.0003019;AN=29816;ANN=T|intragenic_variant|MODIFIER|PF3D7_MIT00100|PF3D7_MIT00100|gene_variant|PF3D7_MIT00100|||n.19A>T||||||;AS_InbreedingCoeff=0.0437;AS_QD=20.96;BaseQRankSum=0.141;CDS;DP=559138;ExcessHet=0.0026;FS=2.892;InbreedingCoeff=0.52;MLEAC=8;MLEAF=0.0002683;MQ=60;MQRankSum=0;QD=20.96;ReadPosRankSum=0.856;RegionType=Mitochondrion;SOR=0.789;set=genotypegvcfs  GT:AD:DP:GQ:PL  0/0:52,0:52:99:0,120,1800  0/0:54,0:54:99:0,120,1800

A small subset of samples I am working on can be found here:
https://www.mediafire.com/file/6dkct6guk5zx83m/sample.vcf/file
I am trying to calculate allele frequency (AF) for each variant in the dataset. The population is assumed to be all the samples in the file.
In other words, even though the VCF file contains AF statistics from population data, I would like to know how to recalculate AF only from the VCF file. I believe I need to know how AC and AN are calculated first.
AF is calculated this way: allele count (AC) / allele number (AN)
I tried to calculate them from the information from GT (counting 0/0, 1/1 and 0/1 in all the 200 samples), but then I learned that this calculation is not suitable for mithocondrial data.
All lines in the VCF file contain an AF statistic:
$ grep -v 'AF=' sample.vcf | grep -v '^#'


The information is provided under INFO column, but how were they calculated? If I do not have AC and AN to begin with, what should I do?
Is there a way to calculate AC and AN from the other information provided in the dataset? such as the data mentioned in the FORMAT column?
GT: genotype
AD: unfiltered allele depth
DP : read depth at this position for this sample (Integer)
GQ : conditional genotype quality, encoded as a phred quality
PL : the phred-scaled genotype likelihoods rounded to the closest integer
Thanks

**Tags:** vcf, allele-frequency
> **Comments:**
> **username_1**: Can you give us some more context? What is AC? Allele count in the population? From GnomAD, maybe? What tool produced this file?  
> **username_0**: @terdon Thanks for the reply. I believe that is what AC means. I have not seen any other definition.  

---

## 1 answers

### Answer by username_2 on 2023-01-30. Score: 2

It sounds like your main confusion is how to calculate AF. I'll answer that right away, but I'll also suggest a Python library that you can use to do this and much more.
You're correct, the AF is defined as the number of non-missing alternate alleles (AC) divided the number of non-missing alleles (AN). Calculating these are straightforward. Map each GT to its alternate alleles:
0/0 => 0
0/1 => 1
1/1 => 2
./0 => 0
./1 => 1
./. => 0

Now, compute the sum. That is the AC. For the AN, map each GT to its number of non-missing alleles:
0/0 => 2
0/1 => 2
1/1 => 2
./0 => 1
./1 => 1
./. => 0

Now, compute the sum. That is the AN. Now calculate AF as
AC / AN

All of this is kind of annoying to do though! I suggest you use Hail (disclaimer: I am a Hail developer).
import hail as hl

# Hail is not yet aware of Plasmodium falciparum reference genome
mt = hl.import_vcf('sample.vcf', reference_genome=None)

# let's peek at our data
mt.show(n_rows=10, n_cols=5)

# Let's get some information about the GT specifically:
mt.GT.summarize()

# We can use Hail's built in variant_qc method to calculate AF for us:
mt = hl.variant_qc(mt)

# Let's take a look at that
mt.variant_qc.show()

# We can overwrite the old AF field and export back to a VCF file if necessary
mt = mt.annotate_rows(
    # AF[1] is the frequency of the first alternate allele
    # Your dataset appears to have some multi-allelics. The other
    # alternate allele frequencies are available in AF[2], ...
    info = mt.info.annotate(
        AF = mt.variant_qc.AF[1]
    )
)
hl.export_vcf(mt, 'sample-updated.vcf.bgz')


This comment doesn't make sense to me:

I learned that this calculation is not suitable for mithocondrial data.

The definition of AF, AC, and AN doesn't change depending on the type of config (autosome, sex chromosome, mitochondrial). Perhaps the issue is ploidy encoding? A GT in a VCF can be:

1-ploidy, e.g. 0 (reference), 1 (first alternate), 2 (second alternate)
2-ploidy, e.g. 0/0 (homozygous reference), 0/1 (heterozygous), 1/1 (homozygous first alternate)
etc.

Sometimes people incorrectly encode a 1-ploidy config (for example the Y chromosome) as: 0/0 and 1/1 instead of 0 and 1. That's not your fault though and it doesn't mean the definition of AF, AC, and AN change. It means someone gave you malformed data.

> **Comments:**
> **username_0**: Thanks for the reply. You totally understood my question. Reagrding the part you said that does not make sense, what i meant was that cells can have multiple mts. From the info i gathered so far, due to this reason the simple equation of AF = AC / AN would not be appropriate for mt genome meaning using GT alone will not return the right result. The calculations used in softwares such as vcftools are only for non mt genome. I was trying to do it in Python anyway. I am glad you introduced this library. Does it work on Windows too?  
> **username_2**: We don't regularly test against windows, but, in principle, it should work. Let me know if you run into bugs and we'll fix them!  

---

