Notice:
If you happen to see a question you know the answer to, please do chime in and help your fellow community members. We appreciate your help!

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.

mutect2: quality column in the vcf file is empty

Hi,

I ran mutect2 (version 4.1.2.0) and in the output vcf file the quality column was empty. Can you please help me with this.

Thanks

Tagged:

Answers

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    HI @efratklig

    Can you please post the exact command you used and some records of the output vcf.

  • efratkligefratklig IsraelMember

    Hi @bhanuGandham

    The command is:

    …/GenomeAnalysisTK-4.1.2.0/gatk --java-options "-Xmx8g -XX:ParallelGCThreads=4" Mutect2 -R …/Homo_sapiens_assembly38.fasta -I …/Sample_1_pos.rg.sorted.markDups.recal.bam -I …/Sample_1_germline.rg.sorted.markDups.recal.bam -normal Sample_1_germline -L …/truseq-exome-targeted-regions-manifest-v1-2_converted_to_hg38_sorted_byChrAsDict_1-22_M_alt.bed -O …/Sample_1_somatic.chr_1-22.vcf --germline-resource …/af-only-gnomad.hg38.vcf --panel-of-normals …/1000g_pon.hg38.vcf --interval-padding 100 && …/GenomeAnalysisTK-4.1.2.0/gatk --java-options "-Xmx8g -XX:ParallelGCThreads=4" Mutect2 -R …/Homo_sapiens_assembly38.fasta -I …/Sample_1_pos.rg.sorted.markDups.recal.bam -I …/Sample_1_germline.rg.sorted.markDups.recal.bam -normal Sample_1_germline -L …/truseq-exome-targeted-regions-manifest-v1-2_converted_to_hg38_sorted_byChrAsDict_X-Y.bed --germline-resource …/af-only-gnomad.hg38.vcf --panel-of-normals …/1000g_pon.hg38.vcf -O …/Sample_1_somatic.chr_XY.vcf --interval-padding 100 && grep '^#' …/Sample_1_somatic.chr_1-22.vcf > …/Sample_1_somatic.vcf && zgrep -v '^#' …/Sample_1_somatic.chr_1-22.vcf >> …/Sample_1_somatic.vcf && grep -v '^#' …/Sample_1_somatic.chr_XY.vcf >> …/Sample_1_somatic.vcf

    After it I ran the FilterMutectCalls:

    …/GenomeAnalysisTK-4.1.2.0/gatk --java-options "-Xmx8g -XX:ParallelGCThreads=4" FilterMutectCalls -R …/Homo_sapiens_assembly38.fasta -V …/Sample_1_somatic.vcf -O …/Sample_1_filtered.vcf

    Some records of the output vcf:

    CHROM POS ID REF ALT QUAL FILTER INFO FORMAT Sample_1_pos Sample_1_germline

    chr1 1204500 . T G . base_qual;normal_artifact;strand_bias;weak_evidence CONTQ=93;DP=137;ECNT=1;GERMQ=93;MBQ=20,17;MFRL=161,171;MMQ=60,60;MPOS=19;NALOD=-1.233e-01;NLOD=8.98;POPAF=2.68;SEQQ=1;STRANDQ=1;TLOD=4.06 GT:AD:AF:DP:F1R2:F2R1:SB 0/1:58,11:0.122:69:19,1:29,1:28,30,0,11 0/0:55,8:0.075:63:22,0:26,0:27,28,0,8

    Thanks

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    HI @efratklig

    For questions like this, please refer to vcf the spec first. You will find most of your information there with regards to vcf file formats and meaning of various values in the vcf fields. In VCF, "." at QUAL means a missing value – i.e. the QUAL is unknown.

  • efratkligefratklig IsraelMember

    Hi @bhanuGandham ,

    Thanks for your answer, I saw it, but my question is why the QUAL is unknown for all of the variants in the vcf file, specially for variants with a filter "PASS"

    Thanks

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin
    edited September 20

    Hi @efratklig

    Can you please post a few records of variants that have QUAL "." and Filter "PASS".

    Post edited by bhanuGandham on
  • efratkligefratklig IsraelMember

    Hi @bhanuGandham ,

    Examples of variants that have QUAL "." and Filter "PASS":

    CHROM POS ID REF ALT QUAL FILTER INFO FORMAT Sample_1_pos Sample_1_germline

    chr1 31658423 . C T . PASS CONTQ=93;DP=143;ECNT=1;GERMQ=93;MBQ=36,34;MFRL=168,171;MMQ=60,60;MPOS=29;NALOD=1.75;NLOD=16.47;POPAF=6.00;SEQQ=93;STRANDQ=93;TLOD=52.01 GT:AD:AF:DP:F1R2:F2R1:SB 0/1:61,22:0.270:83:34,9:25,13:38,23,16,6 0/0:55,0:0.017:55:37,0:17,0:38,17,0,0
    chr1 39770451 . C G . PASS CONTQ=93;DP=281;ECNT=1;GERMQ=93;MBQ=26,37;MFRL=160,163;MMQ=60,60;MPOS=28;NALOD=2.04;NLOD=32.29;POPAF=6.00;SEQQ=93;STRANDQ=93;TLOD=131.64 GT:AD:AF:DP:F1R2:F2R1:SB 0/1:115,53:0.316:168:67,31:45,18:62,53,37,16 0/0:108,0:9.081e-03:108:58,0:46,0:60,48,0,0
    chr7 107211262 . A G . PASS CONTQ=93;DP=366;ECNT=1;GERMQ=93;MBQ=36,36;MFRL=164,163;MMQ=60,60;MPOS=24;NALOD=2.13;NLOD=39.54;POPAF=6.00;SEQQ=93;STRANDQ=93;TLOD=332.55 GT:AD:AF:DP:F1R2:F2R1:SB 0/1:105,118:0.529:223:50,52:52,62:25,80,27,91 0/0:132,0:7.360e-03:132:62,0:68,0:38,94,0,0

    Thanks

  • davidbendavidben BostonMember, Broadie, Dev ✭✭✭

    @efratklig Mutect2 does not emit a QUAL. The _____Q annotations are essentially qual scores for individual filters. For example SEQQ=93; STRANDQ=30 means we think the probability of it being a strand artifact is 1 in 1000 but it is almost certainly not a sequencing error.

  • efratkligefratklig IsraelMember

    Hi @davidben,

    Thanks for your answer. So as I understand the Qual column is empty in Mutect2 and instead we must look on the values of SEQQ (quality of the sequence?) and STRANDQ(quality of the strand?). In my examples when the values of SEQQ and STRANDQ were equal to 93 the filter was "PASS" and when they were equal to 1 the filter was not "PASS" , and only when the filter is pass the varient is probably not a sequencing error, does it right?

    Thanks

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    Hi @efratklig

    SEQQ is Phred-scaled quality that alt alleles are not sequencing errors
    STRANDQ is Phred-scaled quality of strand bias artifact

    As David mentioned Mutect2 does not ever emit a QUAL field. FilterMutectCalls passes variants by estimating relative likelihoods of being an alt allele in a real subclone vs one of several kinds of artifact. Some of those qualities are emitted in the INFO field (they are the properties that end in 'Q'). For more information on documentation of Mutect2 and filtering is maintained here: https://github.com/broadinstitute/gatk/blob/master/docs/mutect/mutect.pdf

  • efratkligefratklig IsraelMember

    Hi @bhanuGandham,

    Thanks a lot for your answer. Another thing in the same topic, is it recommended to use the optional input from CalculateContamination in the filtering, and is there a site for downloading the file: common-biallelic.vcf?

    Thanks

Sign In or Register to comment.