Entering edit mode
3.6 years ago
wang.yiguan
▴
10
Hi,
I am trying to get the genotype quality(GQ) for all calls using BCFtools, including the homozygous REF(0/0) call.
Here is my code:
bcftools mpileup -Ou --thread 48 -a DP,AD,ADF,ADR \
--max-depth 100 \
--min-MQ 20 \
--regions 4 \
--fasta-ref my.fa \
sample_1.bam | \
bcftools call -m -Ob -f GQ --gvcf 100 --threads 24 -o sample_1_gvcf.bcf
bcftools views sample_1_gvcf.bcf > sample_1_gvcf.gvcf
But I find the output GVCF file only give GQ for variants(the third one in the following example).
For the homozygous REF(the first and second in the following example), there isn't GQ annotation.
4 1776 . G . 285 . DP=45;MQSB=1;MQ0F=0;AN=2;DP4=14,31,0,0;MQ=27;MinDP=45 GT:DP:ADF:ADR:AD 0/0:45:14:31:45
4 1777 . C . 285 . DP=45;MQSB=1;MQ0F=0;AN=2;DP4=14,31,0,0;MQ=27;MinDP=45 GT:DP:ADF:ADR:AD 0/0:45:14:31:45
4 1778 . A T 126 . DP=45;VDB=0.147712;SGB=-0.690438;RPB=0.978049;MQB=1;MQSB=1;BQB=0.990184;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=9,19,5,12;MQ=27 GT:PL:DP:ADF:ADR:AD:GQ 0/1:159,0,222:45:9,5:19,12:28,17:127
I am wondering how could I get the GQ for the non-variant sites.
Thanks,