Can't properly use GATK Mutect2 to call somatic variants from my tumor/normal samples pair (BAM to VCF)

Question

[Samuel]

Hi,

I am trying to peform a variant calling out of my BAM file, to get a VCF using Mutect2 from GATK.

I have two files : my tumor sample BAM, and my paired control BAM (from blood). I tried to use GATK best practices documentation and this command :

gatk Mutect2 -R ${ref} -I results/${bam} --panel-of-normals vcf/1000g_pon.hg38.vcf.gz --germline-resource vcf/af-only-gnomad.hg38.vcf.gz --f1r2-tar-gz vcf_res/${f1r2} -L vcf/${intervals} -O vcf_res/${vcf_gz}
  gatk Mutect2 \
     -R /home/sbesseau/scratch/LOTUS/genome.fa \
     -I child_KTN102-2.bam \
     -I child_KTN102-b.bam \
     -normal " " \
     --germline-resource af-only-gnomad.hg38.vcf.gz \
     --panel-of-normals 1000g_pon.hg38.vcf.gz \
     -O somatic.vcf.gz

Using this command, I wanted to eliminate potential germline kwown variants (gnomad), and panel of normal variants (from population). These data come from : https://www.bcgsc.ca/downloads/morinlab/reference/ Is this idea correct? I want to only obtain somatic variants, with their VAF in population if possible.

After obtaining the final vcf, i get this, which seems too be bad (0 read in vcf) :

INFO Mutect2 - 24133700 read(s) filtered by: MappingQualityReadFilter
0 read(s) filtered by: MappingQualityAvailableReadFilter
0 read(s) filtered by: MappingQualityNotZeroReadFilter
0 read(s) filtered by: MappedReadFilter
0 read(s) filtered by: NotSecondaryAlignmentReadFilter
0 read(s) filtered by: NotDuplicateReadFilter
0 read(s) filtered by: PassesVendorQualityCheckReadFilter
0 read(s) filtered by: NonChimericOriginalAlignmentReadFilter
0 read(s) filtered by: NonZeroReferenceLengthAlignmentReadFilter
0 read(s) filtered by: ReadLengthReadFilter
0 read(s) filtered by: GoodCigarReadFilter
222092798 read(s) filtered by: WellformedReadFilter
246226498 total reads filtered
17:24:35.476 INFO ProgressMeter - KI270757.1:71101 131.9 10332600 78337.2
17:24:35.476 INFO ProgressMeter - Traversal complete. Processed 10332600 total regions in 131.9 minutes.

Thank you very much, I"m having a lot of trouble trying to use Mutect2.

[GenoMax]

Please include the FilterMutectCalls command as well.

[Samuel]

gatk FilterMutectCalls -V test.vcf -O filtered.vcf -R /home/sbesseau/scratch/LOTUS/genome.fa

but I can't work because the first command creates a vcf with 0 variant : how to filter out germline mutations (and PON) and add POPAF?

[GenoMax]

Looks like you are not providing name of the "normal" sample in your command line. Use the sample name that is in the read group for the blood file.

[Samuel]

I don't understand where to find this name sorry, I'm a bit new to genomics. Can you precise?

[GenoMax]

GATK requires read groups to be present in your BAM files. You must have created these when you ran the alignments. What "sample name" did you use there? It should be the name you specified for SM : https://gatk.broadinstitute.org/hc/en-us/articles/360035890671-Read-groups

[Samuel]

samtools addreplacerg -r '@RG\tID:samplename\tSM:samplename' child_KTN102-2.bam -o child3_KTN102-2.bam

I used this for the test, so the SM tag should be samplename, but when I use :

  gatk Mutect2 \
     -R /home/sbesseau/scratch/LOTUS/genome.fa \
     -I child2_KTN102-2.bam \
     -I child_KTN102-b.bam \
     -normal samplename \
     --germline-resource af-only-gnomad.hg38.vcf.gz \
     --panel-of-normals 1000g_pon.hg38.vcf.gz \
     -O test2.vcf

I get a java error :

Exception in thread "main" java.lang.IncompatibleClassChangeError: Inconsistent constant pool data in classfile for class org/broadinstitute/hellbender/transformers/ReadTransformer. Method 'org.broadinstitute.hellbender.utils.read.GATKRead lambda$identity$d67512bf$1(org.broadinstitute.hellbender.utils.read.GATKRead)' at index 65 is CONSTANT_MethodRef and should be CONSTANT_InterfaceMethodRef
        at org.broadinstitute.hellbender.transformers.ReadTransformer.identity(ReadTransformer.java:30)
        at org.broadinstitute.hellbender.engine.GATKTool.makePreReadFilterTransformer(GATKTool.java:345)
        at org.broadinstitute.hellbender.engine.AssemblyRegionWalker.traverse(AssemblyRegionWalker.java:276)
        at org.broadinstitute.hellbender.engine.GATKTool.doWork(GATKTool.java:1039)
        at org.broadinstitute.hellbender.cmdline.CommandLineProgram.runTool(CommandLineProgram.java:139)
        at org.broadinstitute.hellbender.cmdline.CommandLineProgram.instanceMainPostParseArgs(CommandLineProgram.java:191)
        at org.broadinstitute.hellbender.cmdline.CommandLineProgram.instanceMain(CommandLineProgram.java:210)
        at org.broadinstitute.hellbender.Main.runCommandLineProgram(Main.java:162)
        at org.broadinstitute.hellbender.Main.mainEntry(Main.java:205)
        at org.broadinstitute.hellbender.Main.main(Main.java:291)

[GenoMax]

Do both BAM files have proper read groups?

[Samuel]

Thank you for anwering! Oh I didn't think about it, so now that I used :

samtools addreplacerg -r '@RG\tID:tumorsample\tSM:tumorsample' child_KTN102-2.bam -o child2_KTN102-2.bam
samtools addreplacerg -r '@RG\tID:bloodsample\tSM:bloodsample' child_KTN102-b.bam -o blood_KTN102-b.bam

I confirm it :

samtools view -H child2_KTN102-2.bam | grep "^@RG" | grep -oP "(?<=SM:)[^\t$]+" | uniq

gives me "samplename" and

samtools view -H blood_KTN102-b.bam | grep "^@RG" | grep -oP "(?<=SM:)[^\t$]+" | uniq

gives me "bloodsample"

And then the command again :

  gatk Mutect2 \
     -R /home/sbesseau/scratch/LOTUS/genome.fa \
     -I child2_KTN102-2.bam \
     -I blood_KTN102-b.bam \
     -normal bloodsample \
     --germline-resource af-only-gnomad.hg38.vcf.gz \
     --panel-of-normals 1000g_pon.hg38.vcf.gz \
     -O minus_blood_PON_germ.vcf

but I still get :

Exception in thread "main" java.lang.IncompatibleClassChangeError: Inconsistent constant pool data in classfile for class org/broadinstitute/hellbender/transformers/ReadTransformer. Method 'org.broadinstitute.hellbender.utils.read.GATKRead lambda$identity$d67512bf$1(org.broadinstitute.hellbender.utils.read.GATKRead)' at index 65 is CONSTANT_MethodRef and should be CONSTANT_InterfaceMethodRef at org.broadinstitute.hellbender.transformers.ReadTransformer.identity(ReadTransformer.java:30) at org.broadinstitute.hellbender.engine.GATKTool.makePreReadFilterTransformer(GATKTool.java:345) at org.broadinstitute.hellbender.engine.AssemblyRegionWalker.traverse(AssemblyRegionWalker.java:276) at org.broadinstitute.hellbender.engine.GATKTool.doWork(GATKTool.java:1039) at org.broadinstitute.hellbender.cmdline.CommandLineProgram.runTool(CommandLineProgram.java:139) at org.broadinstitute.hellbender.cmdline.CommandLineProgram.instanceMainPostParseArgs(CommandLineProgram.java:191) at org.broadinstitute.hellbender.cmdline.CommandLineProgram.instanceMain(CommandLineProgram.java:210) at org.broadinstitute.hellbender.Main.runCommandLineProgram(Main.java:162) at org.broadinstitute.hellbender.Main.mainEntry(Main.java:205) at org.broadinstitute.hellbender.Main.main(Main.java:291)

Answer 1

[GenoMax]

Looking at the error and searching for it brings up

https://gatk.broadinstitute.org/hc/en-us/community/posts/9946093932059-Java-exception-index-65-is-CONSTANT-MethodRef-and-should-be-CONSTANT-InterfaceMethodRef

and

Unable to generate VCF

You will want to make sure the correct version of Java (Java 8 ) is being used.

[Samuel]

After struggling to create envs with various java versions, I figured it out! Thank you very much, you helped me to solve my problem!