I calculated genotype likelihoods from a BAM file with the command: "samtools mpileup -g -f filename" and got the BCF file.
Then I used "bcftools call -c -v filename" to get the VCF.
In the VCF file, I see the following two:
1 589051 . G A 7.79993 . DP=1;SGB=-0.379885;MQ0F=0;AF1=1;AC1=2;DP4=0,0,1,0;MQ=37;FQ=-29.991 GT:PL 1/1:37,3,0
1 700436 . G A 4.77219 . DP=1;SGB=-0.379885;MQ0F=0;AF1=1;AC1=2;DP4=0,0,0,1;MQ=37;FQ=-29.9923 GT:PL 0/1:33,3,0
In the top SNP I get that the most likely genotype is 1/1 (given by the lowest Phred score of 0). Therefore the called genotype is 1/1.
However in the below one even though the loweest Phred score of 0 belongs to the 1/1 genotype, the called genotype is given by 0/1.
What am I missing?
I am new so this is probably obvious to many.
Thanks.
both variants doesn't mean much things: they both have only one read(DP=1) overlapping the variant....
I know. I am just trying to understand how the genotype likelihood part works.
Do you have a comment about that?
Thanks.