Bug Bulletin: The GenomeLocPArser error in SplitNCigarReads has been fixed; if you encounter it, use the latest nightly build.

Gaussian Mixture model plot interpret

mahyarheymahyarhey BostonPosts: 37Member

I got this plot after VariantRecalibration for 42 samples in a VCF file. As it can bee seen in the plot there is no "known" variants detected. What is the problem? Which walker do you recommend to solve this issue? thanks

Best Answer

Answers

  • Geraldine_VdAuweraGeraldine_VdAuwera Posts: 6,286Administrator, GATK Developer admin

    What was your command line?

    Geraldine Van der Auwera, PhD

  • mahyarheymahyarhey BostonPosts: 37Member

    I used the following commands:

    bsub -q short -W 12:0 -R "rusage[mem=32000]" -N -o /hms/scratch1/mahyar/error.log java -jar GenomeAnalysisTK-2.8-1-g932cd3a/GenomeAnalysisTK.jar \ -T VariantRecalibrator \ -R /hms/scratch1/mahyar/ucsc.hg19.fasta \ --input /hms/scratch1/mahyar/Overal-42post-RGSM-allsites.vcf \ --resource:dbsnp,VCF,known=false,training=true,truth=true,prior=6.0 /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/dbsnp_137.hg19.vcf \ -an QD -an HaplotypeScore -an MQRankSum -an ReadPosRankSum -an FS -an MQ \ --mode SNP \ -rf BadCigar \ --recal_file /hms/scratch1/mahyar/All42_post_VQRS.recal \ --tranches_file /hms/scratch1/mahyar/All42_post_VQRS.tranches \ --rscript_file /hms/scratch1/mahyar/All42_post_VQRS_plots.R

  • Geraldine_VdAuweraGeraldine_VdAuwera Posts: 6,286Administrator, GATK Developer admin

    The reason you have no known variants in the plots is because you're not providing any set of knowns (you have known=false for the one resource you provide).

    The bigger problem here is that you're not following our Best Practices for variant recalibration. This command will give you very poor results. Please read the documentation on the Best Practices to learn how you should do this.

    Geraldine Van der Auwera, PhD

  • mahyarheymahyarhey BostonPosts: 37Member

    Do I need use all 4 resources (e.g. hapmap, omni, 1000G, dbsnp) for the VariantRecalibrator or only one resource is enough? I used "dbsnp" resource, because I used it for calling variants via UnifiedGenotyper!

  • mahyarheymahyarhey BostonPosts: 37Member

    I run again VariantRecalibrator only for one sample and got the following error: What is the issue?

    INFO 14:43:50,332 HelpFormatter - -------------------------------------------------------------------------------- INFO 14:43:50,336 HelpFormatter - The Genome Analysis Toolkit (GATK) v2.8-1-g932cd3a, Compiled 2013/12/06 16:47:15 INFO 14:43:50,336 HelpFormatter - Copyright (c) 2010 The Broad Institute INFO 14:43:50,336 HelpFormatter - For support and documentation go to http://www.broadinstitute.org/gatk INFO 14:43:50,342 HelpFormatter - Program Args: -T VariantRecalibrator -R /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/ucsc.hg19.fasta --input /groups/body/Mahyar/Danny/TEST-post/12V_post/12V_post_afterQC.vcf -resource:hapmap,known=false,training=true,truth=true,prior=15.0 /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/hapmap_3.3.hg19.vcf -resource:omni,known=false,training=true,truth=true,prior=12.0 /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/1000G_omni2.5.hg19.vcf -resource:dbsnp,known=true,training=false,truth=false,prior=2.0 /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/dbsnp_137.hg19.vcf -an DP -an QD -an HaplotypeScore -an MQRankSum -an ReadPosRankSum -an FS -an MQ -mode SNP -tranche 100.0 -tranche 99.9 -tranche 99.0 -tranche 90.0 -rf BadCigar --recal_file /hms/scratch1/mahyar/Danny/data/TEST/12V_post_VQRS.recal --tranches_file /hms/scratch1/mahyar/Danny/data/TEST/12V_post_VQRS.tranches --rscript_file /hms/scratch1/mahyar/Danny/data/TEST/12V_post_VQRS_plots.R INFO 14:43:50,342 HelpFormatter - Date/Time: 2014/01/24 14:43:50 INFO 14:43:50,342 HelpFormatter - -------------------------------------------------------------------------------- INFO 14:43:50,342 HelpFormatter - -------------------------------------------------------------------------------- INFO 14:43:50,372 ArgumentTypeDescriptor - Dynamically determined type of /groups/body/Mahyar/Danny/TEST-post/12V_post/12V_post_afterQC.vcf to be VCF INFO 14:43:50,405 ArgumentTypeDescriptor - Dynamically determined type of /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/hapmap_3.3.hg19.vcf to be VCF INFO 14:43:50,420 ArgumentTypeDescriptor - Dynamically determined type of /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/1000G_omni2.5.hg19.vcf to be VCF INFO 14:43:50,425 ArgumentTypeDescriptor - Dynamically determined type of /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/dbsnp_137.hg19.vcf to be VCF INFO 14:43:51,790 GenomeAnalysisEngine - Strictness is SILENT INFO 14:43:51,910 GenomeAnalysisEngine - Downsampling Settings: Method: BY_SAMPLE, Target Coverage: 1000 INFO 14:43:51,946 RMDTrackBuilder - Loading Tribble index from disk for file /groups/body/Mahyar/Danny/TEST-post/12V_post/12V_post_afterQC.vcf INFO 14:43:52,028 RMDTrackBuilder - Loading Tribble index from disk for file /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/hapmap_3.3.hg19.vcf INFO 14:43:52,115 RMDTrackBuilder - Loading Tribble index from disk for file /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/1000G_omni2.5.hg19.vcf WARN 14:43:52,166 RMDTrackBuilder - Index file /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/1000G_omni2.5.hg19.vcf.idx is out of date (index older than input file), deleting and updating the index file INFO 14:43:52,191 RMDTrackBuilder - Creating Tribble index in memory for file /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/1000G_omni2.5.hg19.vcf INFO 14:44:02,429 RMDTrackBuilder - Writing Tribble index to disk for file /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/1000G_omni2.5.hg19.vcf.idx INFO 14:44:03,872 RMDTrackBuilder - Loading Tribble index from disk for file /groups/body/JDM_RNA_Seq-2012/GATK/bundle-2.3/ucsc.hg19/dbsnp_137.hg19.vcf INFO 14:44:04,093 GenomeAnalysisEngine - Preparing for traversal INFO 14:44:04,107 GenomeAnalysisEngine - Done preparing for traversal INFO 14:44:04,107 ProgressMeter - [INITIALIZATION COMPLETE; STARTING PROCESSING] INFO 14:44:04,108 ProgressMeter - Location processed.sites runtime per.1M.sites completed total.runtime remaining INFO 14:44:04,133 TrainingSet - Found hapmap track: Known = false Training = true Truth = true Prior = Q15.0 INFO 14:44:04,134 TrainingSet - Found omni track: Known = false Training = true Truth = true Prior = Q12.0 INFO 14:44:04,134 TrainingSet - Found dbsnp track: Known = true Training = false Truth = false Prior = Q2.0 INFO 14:44:34,112 ProgressMeter - chr1:189685554 3.15e+06 30.0 s 9.0 s 6.0% 8.3 m 7.8 m INFO 14:45:04,114 ProgressMeter - chr2:134069987 6.95e+06 60.0 s 8.0 s 12.2% 8.2 m 7.2 m INFO 14:45:34,115 ProgressMeter - chr3:86393896 1.07e+07 90.0 s 8.0 s 18.5% 8.1 m 6.6 m INFO 14:46:04,117 ProgressMeter - chr4:82027426 1.45e+07 120.0 s 8.0 s 24.6% 8.1 m 6.1 m INFO 14:46:34,119 ProgressMeter - chr5:83285029 1.82e+07 2.5 m 8.0 s 30.8% 8.1 m 5.6 m INFO 14:47:04,120 ProgressMeter - chr6:98389806 2.21e+07 3.0 m 8.0 s 37.0% 8.1 m 5.1 m INFO 14:47:34,124 ProgressMeter - chr7:131997281 2.61e+07 3.5 m 8.0 s 43.5% 8.0 m 4.5 m INFO 14:48:04,125 ProgressMeter - chr9:17446123 3.01e+07 4.0 m 7.0 s 49.6% 8.1 m 4.1 m INFO 14:48:34,126 ProgressMeter - chr10:101411137 3.41e+07 4.5 m 7.0 s 56.8% 7.9 m 3.4 m INFO 14:49:04,129 ProgressMeter - chr12:26405210 3.80e+07 5.0 m 7.0 s 63.0% 7.9 m 2.9 m INFO 14:49:34,130 ProgressMeter - chr13:112799222 4.19e+07 5.5 m 7.0 s 70.1% 7.9 m 2.4 m INFO 14:50:04,131 ProgressMeter - chr16:13840622 4.57e+07 6.0 m 7.0 s 77.3% 7.8 m 105.0 s INFO 14:50:34,133 ProgressMeter - chr18:48814933 4.96e+07 6.5 m 7.0 s 83.8% 7.8 m 75.0 s INFO 14:51:04,134 ProgressMeter - chr22:17655260 5.36e+07 7.0 m 7.0 s 90.8% 7.7 m 42.0 s INFO 14:51:25,743 VariantDataManager - DP: mean = 27.12 standard deviation = 47.63 INFO 14:51:25,769 VariantDataManager - QD: mean = 25.49 standard deviation = 10.51 INFO 14:51:25,789 VariantDataManager - HaplotypeScore: mean = 44.88 standard deviation = 243.46 INFO 14:51:25,803 VariantDataManager - MQRankSum: mean = 0.01 standard deviation = 0.90 INFO 14:51:25,838 VariantDataManager - ReadPosRankSum: mean = 0.01 standard deviation = 1.21 INFO 14:51:25,856 VariantDataManager - FS: mean = 1.27 standard deviation = 3.41 INFO 14:51:25,871 VariantDataManager - MQ: mean = 49.98 standard deviation = 0.53 INFO 14:51:26,035 VariantDataManager - Annotations are now ordered by their information content: [MQ, HaplotypeScore, DP, QD, FS, ReadPosRankSum, MQRankSum] INFO 14:51:26,053 VariantDataManager - Training with 30280 variants after standard deviation thresholding. INFO 14:51:26,058 GaussianMixtureModel - Initializing model with 100 k-means iterations... INFO 14:51:30,143 VariantRecalibratorEngine - Finished iteration 0. INFO 14:51:31,724 VariantRecalibratorEngine - Finished iteration 5. Current change in mixture coefficients = 0.60233 INFO 14:51:33,133 VariantRecalibratorEngine - Finished iteration 10. Current change in mixture coefficients = 0.10546 INFO 14:51:34,135 ProgressMeter - chrY:59358202 5.66e+07 7.5 m 7.0 s 98.7% 7.6 m 6.0 s INFO 14:51:34,540 VariantRecalibratorEngine - Finished iteration 15. Current change in mixture coefficients = 0.07049 INFO 14:51:35,947 VariantRecalibratorEngine - Finished iteration 20. Current change in mixture coefficients = 0.14396 INFO 14:51:37,378 VariantRecalibratorEngine - Finished iteration 25. Current change in mixture coefficients = 0.04686 INFO 14:51:38,819 VariantRecalibratorEngine - Finished iteration 30. Current change in mixture coefficients = 0.01551 INFO 14:51:40,256 VariantRecalibratorEngine - Finished iteration 35. Current change in mixture coefficients = 0.00538 INFO 14:51:41,694 VariantRecalibratorEngine - Finished iteration 40. Current change in mixture coefficients = 0.00206 INFO 14:51:41,982 VariantRecalibratorEngine - Convergence after 41 iterations! INFO 14:51:42,259 VariantDataManager - Training with worst 0 scoring variants --> variants with LOD <= -5.0000. INFO 14:51:44,126 GATKRunReport - Uploaded run statistics report to AWS S3

    ERROR ------------------------------------------------------------------------------------------
    ERROR stack trace

    java.lang.NullPointerException at org.broadinstitute.sting.gatk.walkers.variantrecalibration.VariantRecalibratorEngine.generateModel(VariantRecalibratorEngine.java:83) at org.broadinstitute.sting.gatk.walkers.variantrecalibration.VariantRecalibrator.onTraversalDone(VariantRecalibrator.java:359) at org.broadinstitute.sting.gatk.walkers.variantrecalibration.VariantRecalibrator.onTraversalDone(VariantRecalibrator.java:139) at org.broadinstitute.sting.gatk.executive.Accumulator$StandardAccumulator.finishTraversal(Accumulator.java:129) at org.broadinstitute.sting.gatk.executive.LinearMicroScheduler.execute(LinearMicroScheduler.java:116) at org.broadinstitute.sting.gatk.GenomeAnalysisEngine.execute(GenomeAnalysisEngine.java:313) at org.broadinstitute.sting.gatk.CommandLineExecutable.execute(CommandLineExecutable.java:113) at org.broadinstitute.sting.commandline.CommandLineProgram.start(CommandLineProgram.java:245) at org.broadinstitute.sting.commandline.CommandLineProgram.start(CommandLineProgram.java:152) at org.broadinstitute.sting.gatk.CommandLineGATK.main(CommandLineGATK.java:91)

    ERROR ------------------------------------------------------------------------------------------
    ERROR A GATK RUNTIME ERROR has occurred (version 2.8-1-g932cd3a):
    ERROR
    ERROR This might be a bug. Please check the documentation guide to see if this is a known problem.
    ERROR If not, please post the error message, with stack trace, to the GATK forum.
    ERROR Visit our website and forum for extensive documentation and answers to
    ERROR commonly asked questions http://www.broadinstitute.org/gatk
    ERROR
    ERROR MESSAGE: Code exception (see stack trace for error itself)
    ERROR ------------------------------------------------------------------------------------------
Sign In or Register to comment.