Notice:
If you happen to see a question you know the answer to, please do chime in and help your fellow community members. We encourage our fourm members to be more involved, jump in and help out your fellow researchers with their questions. GATK forum is a community forum and helping each other with using GATK tools and research is the cornerstone of our success as a genomics research community.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.

PON: Mutect2 include MNPs and crash GenomicsDBImport

manolismanolis ✭✭Member ✭✭
edited April 20 in Ask the GATK team

GATK v4.1.1.0, linux server, bash

Hi,

using the following command I create the input files for the pon:

${gatk}  Mutect2 \
-R ${hg38} \
-I ${bqsr_bam} \
-O ${sample_pon} \
-L ${intervals}

When I run GenomicsDBImport I have an error which reports the presence of MNPs in the ${sample_pon} files

Then I exclude the MNPs with:

${gatk} SelectVariants \
-R ${hg38} \
-V ${sample_pon} \
-O ${clean_sample_pon} \
-xl-select-type MNP

If I go to count the variants before and after SelectVariants I have the following counts

sample_pon_1
62223
clean_sample_pon_1
61395

sample_pon_2
66974
clean_sample_pon_2
66013

How can I deactivate the call of MNPs in Mutect2 in order to use GenomicsDBImport?
Is it the default behavior of Mutect2/pon-pipeline or I'm doing something wrong?

Many thanks

Post edited by manolis on

Best Answer

Answers

  • SkyWarriorSkyWarrior ✭✭✭ TurkeyMember ✭✭✭

    Neither GenomicsDBImport does support mutect2 nor does mutect2 have a gvcf mode. How can you import them in the first place?

  • manolismanolis ✭✭ Member ✭✭

    Hi @SkyWarrior !

    In the last version, v4.1.1.0, according to the pon pipeline you can use M2+GDBI+CreateSomaticPanelOfNormals in case you have to process a large number of samples, this is what I understand, hoping correctly...

           gatk GenomicsDBImport --genomicsdb-workspace-path pon_db -R ${ref_fasta} -V ${sep=' -V ' input_vcfs} -L ${intervals}
    
            gatk --java-options "-Xmx${command_mem}g"  CreateSomaticPanelOfNormals -R ${ref_fasta} --germline-resource ${gnomad} \
                -V gendb://pon_db -O ${output_vcf_name}.vcf ${create_pon_extra_args}
    

    I found my error, in HC the default value of "--max-mnp-distance" is 0 but in the M2 tool is 1.

    I ran M2 with the value 0 and then I filtered the MNPs with SelecteVariants... now, before and after SelectVariants the number of variants is the same... just for a fast check.

    @bhanuGandham, please, is correct the solution of my post?

    Best

  • SkyWarriorSkyWarrior ✭✭✭ TurkeyMember ✭✭✭
    edited April 21

    Wow this is uncalled for. Changing a default parameter is really dangerous. I guess this must be reverted in the next version.

    Also another user is reporting a similar problem with haplotypecaller and combinegvcfs. Is it possible the MNPs are somehow emitted due to some parameter change and nobody is really aware of. I have to check my own GVCFs. I have not done joint genotyping for a long time.

  • manolismanolis ✭✭ Member ✭✭
    edited April 21

    Personally I ran the CNN pipeline and now I want to run the GDBI-VQSR pipeline because I want to compare the outputs.

    Moreover, if I'm going to create a pon file with Mutect2 "--max-mnp-distance 0", then when I'm going to use this pon with a tumor sample (somatic SNPs/Indels discovery) again in the Mutect2 I have to set up 0, correct?

  • SkyWarriorSkyWarrior ✭✭✭ TurkeyMember ✭✭✭

    That sounds the most logical way.

  • bhanuGandhambhanuGandham admin Cambridge MAMember, Administrator, Broadie, Moderator admin

    Hi @SkyWarrior and @manolis,

    I agree that is weird. Let me check with the dev team and get back to you shortly.

  • manolismanolis ✭✭ Member ✭✭
  • syer89syer89 Member
    Hi @bhanuGandham , I have same problem with v4.1.2.0 GenomicsDBImport for creating PON and I wanted to confirm this bug is still active in the new version as well, and we continue using --max-mnp-distance 0 with Mutect2 ?
  • AdelaideRAdelaideR admin Unconfirmed, Member, Broadie, Moderator admin

    It appears that this issue has not been closed. You can follow the progress on github here

  • syer89syer89 Member
    edited May 20
    Hi @AdelaideR @davidben . That seems like a different issue (unless I got it wrong). When creating PON it's advised to use Mutect2 "--max-mnp-distance 0" as GenomicsDBImport cannot handle otherwise. When compared the results between default (1) and 0 and results look like expected and not like in the issue above.
    When now analyzing a tumor sample with Mutect2 with this PON downstream, is it tested that it works with default --max-mnp-distance or should we use the same "--max-mnp-distance 0" for sample processing as well ? Please advise.
  • bhanuGandhambhanuGandham admin Cambridge MAMember, Administrator, Broadie, Moderator admin

    Hi @syer89

    Hi @bhanuGandham , I have same problem with v4.1.2.0 GenomicsDBImport for creating PON and I wanted to confirm this bug is still active in the new version as well, and we continue using --max-mnp-distance 0 with Mutect2 ?

    That is correct.

    When now analyzing a tumor sample with Mutect2 with this PON downstream, is it tested that it works with default --max-mnp-distance or should we use the same "--max-mnp-distance 0" for sample processing as well ?

    Downstream it works with the default --max-mnp-distance and does not need "--max-mnp-distance 0" for sample processing.

    I hope this helps clarify the question.

Sign In or Register to comment.