This is a test version of Biostars. For the public version, visit https://www.biostars.org.
bcftools AA annotation failing

I am attempting to annotate a phased vcf file with ancestral alleles (AA), but it is not working.

My original vcf file (less the full header) format is:

#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO    FORMAT  Sample1  Sample2    ...
chr1    2602    .       C       A       .       .       AC=0;AF=0.0155039;CM=0;AN=22    GT      0|0   0|0      0|0     0|0     0|0     0|0     0|0     0|0     0|0     0|0     0|0
chr1    2614    .       G       A       .       .       AC=0;AF=0.0193798;CM=1.2e-05;AN=22      GT    0|0      0|0     0|0     0|0     0|0     0|0     0|0     0|0     0|0     0|0     0|0

My annotation file with the AAs I would like to bring over into the original vcf is:

#CHROM  POS     ID      REF     ALT     QUAL    FILTER  INFO
chr1    2602    .       .       .       .       .       AA=C
chr1    2614    .       .       .       .       .       AA=G
chr1    2649    .       .       .       .       .       AA=C
chr1    2830    .       .       .       .       .       AA=G

Here is the bcftools command:

bcftools annotate -a ancestral_alleles.vcf.gz -c CHROM,POS,INFO/AA -h hdr.txt -o annotated_ancestral.vcf original.vcf.gz

The output vcf file header includes:

##fileformat=VCFv4.2
##FILTER=<ID=PASS,Description="All filters passed">
##fileDate=07/08/2024 - 01:30:48
##source=shapeit4.1.3
##INFO=<ID=AF,Number=A,Type=Float,Description="Allele Frequency">
##INFO=<ID=AC,Number=1,Type=Integer,Description="Allele count">
##INFO=<ID=CM,Number=A,Type=Float,Description="Interpolated cM position">
##FORMAT=<ID=GT,Number=1,Type=String,Description="Phased genotypes">
##INFO=<ID=AN,Number=1,Type=Integer,Description="Total number of alleles in called genotypes">
##INFO=<ID=AA,Number=1,Type=String,Description="Ancestral Allele">
##bcftools_annotateVersion=1.19+htslib-1.19.1

While the file looks like:

chr1    2602    .       C       A       .       .       AC=0;AF=0.0155039;CM=0;AN=22    GT      0|0     0|0     0|0     0|0     0|0     0|0     0|0       0|0     0|0     0|0     0|0
chr1    2614    .       G       A       .       .       AC=0;AF=0.0193798;CM=1.2e-05;AN=22      GT      0|0     0|0     0|0     0|0     0|0     0|0       0|0     0|0     0|0     0|0     0|0

But I would like it to look something like:

chr1    2602    .       C       A       .       .       AC=0;AF=0.0155039;CM=0;AN=22;AA=C    GT      0|0     0|0     0|0     0|0     0|0     0|0     0|0       0|0     0|0     0|0     0|0
chr1    2614    .       G       A       .       .       AC=0;AF=0.0193798;CM=1.2e-05;AN=22;AA=G     GT      0|0     0|0     0|0     0|0     0|0     0|0       0|0     0|0     0|0     0|0     0|0
bcftools vcf ancestral-allele
• 1,407 views • link
written by Whirlingdaf • Whirlingdaf

I suppose it's not working because when you use a vcf as an annotation source, the REF (and ALT ?) allele is required ?

just try again but add the REF/ALT for the two variants above, just to check if that works.

• 0 views • link
written by Pierre Lindenbaum Pierre Lindenbaum

Ah! Maybe you are on to something. Do you expect this would probably also be the case for a tab separated (bgziped and indexed) file as well? E.g.:

chr1    2602    C
chr1    2614    G
chr1    2649    C
chr1    2830    G
• 0 views • link
written by Whirlingdaf • Whirlingdaf

0 answers

No answers yet.

Log in to answer this question.

More posts like this
  • vcf file chr notation
    written by Hriday • Hriday
    "I have a single VCF file named 'ALL.wgs.shapeit2_integrated_snvindels_v2a.GRCh38.27022019.sites.vcf.gz'. The issue at hand is that the file uses different chromosomal notation and lacks the 'chr' prefix. …
  • bcftools view to failed reader data
    written by Jihyun • Jihyun
    It's so confused for data editing with vcftools - bcftools fusion preprocessing. I want to view all results in bcftools, vcftool and zcat when I …
  • Problem to convert genotypes to plink format from bcftools
    written by Alejandro • Alejandro
    I had a vcf file with imputed genotypes from beagle like this: ##fileformat=VCFv4.2 ##filedate=20210504 ##source="beagle.18May20.d20.jar" ##INFO=<ID=AF,Number=A,Type=Float,Description="Estimated ALT Allele Frequencies"> ##INFO=<ID=DR2,Number=1,Type=Float,Description="Dosage R-Squared: estimated squared correlation between …
  • Split multiallelic SNPs to biallelic from vcf
    written by User000 User000
    Dear all, I have a particular vcf file like this, chrX 29 . G A,T . PASS AC=1,1;AN=3 GT:DP:HF:CILOW:CIUP:SDP 0/1/2:4839:0.003,0.001:0.002,0.0:0.005,0.003:14;0,4;2 I tried various tools to …
  • Annotating VCF with centimorgan information (cM)
    written by 4galaxy77 4galaxy77
    I often need to know the interpolated cumulative genetic distance in cM across the genome. I am aware of some software, such as shapeit4, which …
  • Filtering VCF file by INFO flag
    written by zhamouda • zhamouda
    Hi I'm new to bioinformatics and I get an empty file everytime I am trying to filter variants with HapMap3 "HMR". In my VCF file …
  • Get genotype GT from multi-sample VCF (per chromosome)
    written by zx8754 zx8754
    Example input: multi-sample VCF (adapted from [www.internationalgenome.org][1]): **Note:** my actual file is bgzipped and tabixed with ~2Mln variants (rows) and ~1000 samples (columns). ##fileformat=VCFv4.0 ##fileDate=20090805 …
  • About structure variants simulation
    written by Darren Lin • Darren Lin
    Guys, I am trying to use DWGSIM to generate some simulation paired-end reads with a VCF file of candidate mutations. I use the command like …
  • add coordinates to line of output extracted from vcf
    written by bioguy24 bioguy24
    I am trying to extrapolate/convert variants from a vcf file to genotype. I am close but can not seem to add the chr and start …
  • Grep Position List
    written by BAGeno BAGeno
    Hi, I have VCF file and list of position file. I want to grep those records from VCF which have position corresponding to position present …