FilterMutectCalls --unique-alt-read-count argument
New GATK user here. Using mutect2 to detect somatic mutations in a CRISPR/Cas9 edited population of plant clones. Sequencing is targeted illumina short read PE using bait probes for pull-down. I know this is a very unique use of this program but it seemed ideally suited to the task.
Running GATK 4.1.4.1 on terra and GATK 4.1.9.0 locally. Running mutect2 in normal/tumor mode. Mutect2 runs to completion no problem.
My issue concerns FilterMutectCalls. I would like to use the --unique-alt-read-count argument to increase the stringency and eliminate sites that are scored PASS but have only a few unique reads supporting the alternative allele. When I run FilterMutectCalls using --unique-alt-read-count with values of 1,10,50,100,1000, there is no change in the sites that reach a PASS call. Manual inspection of PASS sites in a browser shows some sites should be eliminated based on this argument. When I use the --min-allele-fraction argument I am successfully able to reduce PASS calls based on the allele frequency. Running mutect2 locally in GATK 4.1.9.0 and FilterMutectCalls locally in the same version produces the same result. When running FilterMutectCalls with --unique-alt-read-count I receive no errors and all the expected output files.
How do I get the --unique-alt-read-count parameter to function in FilterMutectCalls?
Thanks very much for your time!
-
Hi Greg Goralogia,
Thanks for writing in about this, we will try to figure out how to get it to work. Could you share some of the examples you found with the sites that should have been eliminated?
Thank you,
Genevieve
-
Hi Genevieve,
Thanks very much for helping me with this. Below is the rough command I am using in FilterMutectCalls. Just as additional information am using a list of intervals corresponding to targeted sequencing regions during the mutect2 run.
java -jar (GATKlocal).jar FilterMutectCalls \-V sampleX(tumor-normal_pair).vcf \-R reference.fasta \-O sampleXfilteredUniquereadcount100.vcf \--unique-alt-read-count 100
final output lines look like this
11:20:34.058 INFO FilterMutectCalls - Initializing engine
11:20:34.446 INFO FeatureManager - Using codec VCFCodec to read file sampleX.vcf
11:20:34.543 INFO FilterMutectCalls - Done initializing engine
11:20:34.670 INFO ProgressMeter - Starting traversal
11:20:34.672 INFO ProgressMeter - Current Locus Elapsed Minutes Variants Processed Variants/Minute
11:20:34.673 INFO FilterMutectCalls - Starting pass 0 through the variants
11:20:35.838 INFO FilterMutectCalls - Finished pass 0 through the variants
11:20:35.924 INFO FilterMutectCalls - Starting pass 1 through the variants
11:20:36.318 INFO FilterMutectCalls - Finished pass 1 through the variants
11:20:36.335 INFO FilterMutectCalls - Starting pass 2 through the variants
11:20:36.660 INFO FilterMutectCalls - Finished pass 2 through the variants
11:20:36.664 INFO FilterMutectCalls - Starting pass 3 through the variants
11:20:37.316 INFO FilterMutectCalls - Finished pass 3 through the variants
11:20:37.329 INFO FilterMutectCalls - No variants filtered by: AllowAllVariantsVariantFilter
11:20:37.329 INFO FilterMutectCalls - 0 read(s) filtered by: AllowAllReadsReadFilter
11:20:37.330 INFO ProgressMeter - Chr18:1087311 0.0 8200 185171.2
11:20:37.330 INFO ProgressMeter - Traversal complete. Processed 8200 total variants in 0.0 minutes.
11:20:37.385 INFO FilterMutectCalls - Shutting down engine
[February 17, 2021 11:20:37 AM PST] org.broadinstitute.hellbender.tools.walkers.mutect.filtering.FilterMutectCalls done. Elapsed time: 0.06 minutes.
Below I'll put the filtered .vcf files, the PASS entries in both number and site location are the same between values of 10 and 100, suggesting this argument isn't working at all for me. I'll also put two screenshots of IGV browser at two sites (with the normal above and tumor sample below), one which got a PASS despite only a single read supporting an insertion. The second has the trappings of a somatic mutation in the tumor sample (C to A blue/green at the left of second screenshot) but should be eliminated by a --unique-alt-read-count value of 10 or 100. As I said in the original post using the argument --min-allele-fraction (0.1/0.2/0.3) I get a reduction of PASS sites with increasing stringency but not with this argument.
I hope I'm not making a silly mistake with the commands or accepted values etc. I also hope the images, entries are informative; let me know if I can send anything else helpful.
Thanks again very much for your help!
-
Hi Greg Goralogia,
Could you check to see if your VCF files have the annotation UniqueAltReadCount? I cannot see from what you have shared, and checking this value would be informative if these positions should be getting a PASS or not.
Thank you,
Genevieve
-
Hi Genevieve,
Thanks again for looking at that. Looking at the .vcf file the unique alt read count line has a value of 1 when it should have been set to 10 or 100 by the argument. I'll try to play around with the command format to see if I can get FilterMutectCalls to apply the value change.
If you have any advice if I am inputting it wrong, that would be great. (final line is output.vcf \--unique-alt-read-count 100)
Thanks again for your time,
Greg
-
Hi Greg,
Is the UNIQ_ALT_READ_COUNT in the INFO column of your variants? That would be where each position is annotated with how many alt reads are unique. The line you showed is from the header, so it does not show information about your specific variants.
Thank you,
Genevieve
-
Hi Genevieve,
Sorry about that; Hope the attached image is where that would be. I don't see any UNIQ ARC in the info column. See the attached image.
Thanks again,
Greg
-
Sorry that might have ben too low of resolution.
-
Hi Greg,
Thanks for sending that. It seems like the reason it did not filter is that you do not have values for this annotation. I am going to check in with my team about the expected behavior of the filter and get you more information.
Best,
Genevieve
-
Hi Greg Goralogia,
This is an advanced filter in FilterMutectCalls used for variant calling with UMIs, we don't recommend it for new users. If you do want to get it to work, you would call the UNIQ_ALT_READ_COUNT annotation with the -A argument in Mutect2.
I spoke with our Mutect2 developers regarding what you should try instead, and they recommended sticking to the default settings with FilterMutectCalls because that is where you are going to get the best results. There are many different parameters in FilterMutectCalls, but most should really only be used in advanced cases because they could cause issues rather than help your results.
One option you could try is the option --f-score-beta, which adjusts the recall vs precision weight.
Hope this helps!
Best,
Genevieve
-
Hi Genevieve,
Thanks very much for looking into that, I really appreciate it. I'll give both the strategies you listed above a try, but as you say likely fall back on the default settings as the developers recommend.
Thanks again for your time!
Greg
Please sign in to leave a comment.
10 comments