Test-drive the GATK tools and Best Practices pipelines on Terra


Check out this blog post to learn how you can get started with GATK and try out the pipelines in preconfigured workspaces (with a user-friendly interface!) without having to install anything.

AD field doesn't update after merging GVCF files with GenotypeGVCFs

leobokileoboki Atlanta, Georgia, USAMember

Hello,

I'm using GATK version 3.1-1 and following Best Practices to call SNPs and INDELs for my exome sequencing data. I noticed that the AD field didn't update properly after I merged per-sample GVCF files into the per-cohort VCF file.

As an example, one SNP on mitochondrion at location 185 has following entries in per-sample GVCF files:
chrM 185 . G T, 1254.77 . DP=23;MLEAC=2,0;MLEAF=1.00,0.00;MQ=52.36;MQ0=0 GT:AD:DP:GQ:PL:SB 1/1:0,23,0:23:93:1288,93,0,1288,93,1288:0,0,0,0
chrM 185 . G T, 410.77 . DP=7;MLEAC=2,0;MLEAF=1.00,0.00;MQ=34.70;MQ0=0 GT:AD:DP:GQ:PL:SB 1/1:0,7,0:7:33:444,33,0,444,33,444:0,0,0,0
chrM 185 . G A, 259.84 . DP=5;MLEAC=2,0;MLEAF=1.00,0.00;MQ=37.32;MQ0=0 GT:AD:DP:GQ:PL:SB 1/1:0,5,0:5:21:293,21,0,293,21,293:0,0,0,0
chrM 185 . G A, 205.90 . DP=3;MLEAC=2,0;MLEAF=1.00,0.00;MQ=60.00;MQ0=0 GT:AD:DP:GQ:PL:SB 1/1:0,3,0:3:18:239,18,0,239,18,239:0,0,0,0
chrM 185 . G T, 1316.77 . DP=21;MLEAC=2,0;MLEAF=1.00,0.00;MQ=35.14;MQ0=0 GT:AD:DP:GQ:PL:SB 1/1:0,21,0:21:99:1350,99,0,1350,99,1350:0,0,0,0
chrM 185 . G T, 171.03 . DP=3;MLEAC=2,0;MLEAF=1.00,0.00;MQ=60.00;MQ0=0 GT:AD:DP:GQ:PL:SB 1/1:0,3,0:3:15:204,15,0,204,15,204:0,0,0,0
All other samples didn't have evidence for this SNP.

After merging with GenotypeGVCFs, the same SNP has the entry in the per-cohort VCF file:
chrM 185 . G T,A 3561.78 . AC=8,4;AF=0.098,0.049;AN=82;DP=565;FS=0.000;InbreedingCoeff=0.8981;MLEAC=8,3;MLEAF=0.098,0.037;MQ=52.36;MQ0=0;QD=32.42 GT:AD:DP:GQ:PL 0/0:.:13:30:0,30,450,30,450,450 0/0:.:16:39:0,39,585,39,585,585 1/1:0,23,0:23:93:1288,93,0,1288,93,1288 0/0:.:8:21:0,21,315,21,315,315 0/0:.:14:27:0,27,405,27,405,405 0/0:.:11:24:0,24,360,24,360,360 0/0:.:2:6:0,6,75,6,75,75 0/0:.:9:18:0,18,270,18,270,270 0/0:.:8:21:0,21,269,21,269,269 0/0:.:18:36:0,36,540,36,540,540 0/0:.:10:21:0,21,315,21,315,315 0/0:.:23:60:0,60,900,60,900,900 0/0:.:45:99:0,102,1530,102,1530,1530 0/0:.:50:99:0,120,1800,120,1800,1800 0/0:.:8:21:0,21,266,21,266,266 0/0:.:16:36:0,36,540,36,540,540 1/1:0,7,0:7:33:444,33,0,444,33,444 0/0:.:5:6:0,6,90,6,90,90 0/0:.:4:9:0,9,135,9,135,135
0/0:.:7:21:0,21,249,21,249,249 0/0:.:12:27:0,27,405,27,405,405 0/0:.:23:60:0,60,900,60,900,900 2/2:0,5,0:5:21:293,293,293,21,21,0 0/0:.:8:21:0,21,278,21,278,278 0/0:.:11:24:0,24,360,24,360,360 0/0:.:8:24:0,24,292,24,292,292 2/2:0,3,0:3:18:239,239,239,18,18,0 0/0:.:12:27:0,27,405,27,405,405 0/0:.:12:21:0,21,315,21,315,315 1/1:0,21,0:21:99:1350,99,0,1350,99,1350 0/0:.:13:27:0,27,405,27,405,405 0/0:.:3:9:0,9,112,9,112,112 1/1:0,3,0:3:15:204,15,0,204,15,204 0/0:.:7:21:0,21,209,21,209,209 0/0:.:13:27:0,27,405,27,405,405 0/0:.:5:6:0,6,90,6,90,90 0/0:.:7:21:0,21,260,21,260,260 0/0:.:10:24:0,24,360,24,360,360 0/0:.:40:93:0,93,1395,93,1395,1395 0/0:.:38:84:0,84,1260,84,1260,1260 0/0:.:14:33:0,33,495,33,495,495

Both alternative alleles were included in the ALT field and the genotype field (GT) was updated correspondingly. For genotype 2/2, I think the correct AD field should be 0,0,5 and 0,0,3 in my aforementioned example. However, these fields remained the same as in the original per-sample GVCF files (i.e., 0,5,0 and 0,3,0). I'm wondering whether my understanding about the AD field is correct or not. And also, if it is an erroneous reporting, whether it'll influence the downstream VQSR and other filtering/annotations or not.

Thanks,
Leo

Best Answer

Answers

Sign In or Register to comment.