VQSR first step do not generating plots
0
0
Entering edit mode
21 months ago
Yoosef ▴ 60

Hello. After running my codes for VQSR's first step, I receive all files (including R scripte for plots) except the pdf file that contains plots. I do not have any errors from the vqsr. My code:

java -jar gatk-package-4.3.0.0-local.jar VariantRecalibrator \
-O /home/yousef/Desktop/VQSR_hg38/haplotype_hg38_SNP_Recal.vcf \
--resource:hapmap,known=false,training=true,truth=true,prior=15.0 '/media/yousef/EEFC0BDBFC0B9CC9/Sequencing/Next_generation_sequencing/Index_files/VQSR_data/Coverted_to_VCF/Coverted_broad_hg38_v0_hapmap_3.3.hg38.vcf' \
--resource:omni,known=false,training=true,truth=false,prior=12.0 '/media/yousef/EEFC0BDBFC0B9CC9/Sequencing/Next_generation_sequencing/Index_files/VQSR_data/Coverted_to_VCF/resources_broad_hg38_v0_1000G_omni2.5.hg38.vcf' \
--resource:1000G,known=false,training=true,truth=false,prior=10.0 '/media/yousef/EEFC0BDBFC0B9CC9/Sequencing/Next_generation_sequencing/Index_files/VQSR_data/Coverted_to_VCF/resources_broad_hg38_v0_1000G_phase1.snps.high_confidence.hg38.vcf' \
--resource:dbsnp,known=true,training=false,truth=false,prior=2.0 '/media/yousef/EEFC0BDBFC0B9CC9/Sequencing/Next_generation_sequencing/Index_files/VQSR_data/resources_broad_hg38_v0_Homo_sapiens_assembly38.dbsnp138.vcf' \
--tranches-file  /home/yousef/Desktop/VQSR_hg38/Tranches.txt \
-an QD -an SOR -an FS -an SOR -an ReadPosRankSum \
-V /home/yousef/Desktop/VQSR_hg38/haplotype_hg38.vcf \
--max-gaussians 4 \
-R /media/yousef/EEFC0BDBFC0B9CC9/Sequencing/Next_generation_sequencing/Index_files/hg38_FASTA/hg38.fa \
--rscript-file /home/yousef/Desktop/VQSR_hg38/Rscript.R \
--mode SNP \
-tranche 90.0 -tranche 93.0 -tranche 95.0 -tranche 97.0 -tranche 98.0 -tranche 99.9 -tranche 100.0

The results from software

19:28:49.256 INFO  NativeLibraryLoader - Loading libgkl_compression.so from jar:file:/media/yousef/EEFC0BDBFC0B9CC9/Sequencing/Next_generation_sequencing/Softwares/gatk-4.3.0.0/gatk-package-4.3.0.0-local.jar!/com/intel/gkl/native/libgkl_compression.so
19:28:49.569 INFO  VariantRecalibrator - ------------------------------------------------------------
19:28:49.569 INFO  VariantRecalibrator - The Genome Analysis Toolkit (GATK) v4.3.0.0
19:28:49.569 INFO  VariantRecalibrator - For support and documentation go to https://software.broadinstitute.org/gatk/
19:28:49.570 INFO  VariantRecalibrator - Executing as yousef@yousef-Katana-GF66-11UE on Linux v5.19.0-32-generic amd64
19:28:49.570 INFO  VariantRecalibrator - Java runtime: OpenJDK 64-Bit Server VM v11.0.17+8-post-Ubuntu-1ubuntu222.04
19:28:49.570 INFO  VariantRecalibrator - Start Date/Time: February 17, 2023 at 7:28:49 PM IRST
19:28:49.570 INFO  VariantRecalibrator - ------------------------------------------------------------
19:28:49.570 INFO  VariantRecalibrator - ------------------------------------------------------------
19:28:49.570 INFO  VariantRecalibrator - HTSJDK Version: 3.0.1
19:28:49.570 INFO  VariantRecalibrator - Picard Version: 2.27.5
19:28:49.570 INFO  VariantRecalibrator - Built for Spark Version: 2.4.5
19:28:49.570 INFO  VariantRecalibrator - HTSJDK Defaults.COMPRESSION_LEVEL : 2
19:28:49.570 INFO  VariantRecalibrator - HTSJDK Defaults.USE_ASYNC_IO_READ_FOR_SAMTOOLS : false
19:28:49.570 INFO  VariantRecalibrator - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_SAMTOOLS : true
19:28:49.570 INFO  VariantRecalibrator - HTSJDK Defaults.USE_ASYNC_IO_WRITE_FOR_TRIBBLE : false
19:28:49.570 INFO  VariantRecalibrator - Deflater: IntelDeflater
19:28:49.570 INFO  VariantRecalibrator - Inflater: IntelInflater
19:28:49.571 INFO  VariantRecalibrator - GCS max retries/reopens: 20
19:28:49.571 INFO  VariantRecalibrator - Requester pays: disabled
19:28:49.571 INFO  VariantRecalibrator - Initializing engine
19:28:49.674 INFO  FeatureManager - Using codec VCFCodec to read file file:///media/yousef/EEFC0BDBFC0B9CC9/Sequencing/Next_generation_sequencing/Index_files/VQSR_data/Coverted_to_VCF/Coverted_broad_hg38_v0_hapmap_3.3.hg38.vcf
19:28:49.733 INFO  FeatureManager - Using codec VCFCodec to read file file:///media/yousef/EEFC0BDBFC0B9CC9/Sequencing/Next_generation_sequencing/Index_files/VQSR_data/Coverted_to_VCF/resources_broad_hg38_v0_1000G_omni2.5.hg38.vcf
19:28:49.758 INFO  FeatureManager - Using codec VCFCodec to read file file:///home/yousef/Desktop/VQSR_hg38/haplotype_hg38.vcf
19:28:49.793 INFO  VariantRecalibrator - Done initializing engine
19:28:49.795 WARN  VariantDataManager - Ignoring duplicate annotations for recalibration %s.
19:28:49.796 INFO  TrainingSet - Found hapmap track:    Known = false   Training = true     Truth = true    Prior = Q15.0
19:28:49.796 INFO  TrainingSet - Found omni track:  Known = false   Training = true     Truth = false   Prior = Q12.0
19:28:49.819 INFO  ProgressMeter - Starting traversal
19:28:49.819 INFO  ProgressMeter -        Current Locus  Elapsed Minutes    Variants Processed  Variants/Minute
19:28:50.502 INFO  ProgressMeter - chr1_KI270711v1_random:20068              0.0                 26862        2359765.7
19:28:50.502 INFO  ProgressMeter - Traversal complete. Processed 26862 total variants in 0.0 minutes.
19:28:50.505 INFO  VariantDataManager - QD:      mean = 28.54    standard deviation = 6.38
19:28:50.509 INFO  VariantDataManager - SOR:     mean = 1.57     standard deviation = 1.33
19:28:50.513 INFO  VariantDataManager - FS:      mean = 0.25     standard deviation = 1.20
19:28:50.515 INFO  VariantDataManager - ReadPosRankSum:      mean = 0.38     standard deviation = 1.00
19:28:50.532 INFO  VariantDataManager - Annotation order is: [QD, SOR, ReadPosRankSum, FS]
19:28:50.534 INFO  VariantDataManager - Training with 10029 variants after standard deviation thresholding.
19:28:50.535 INFO  GaussianMixtureModel - Initializing model with 100 k-means iterations...
19:28:50.679 INFO  VariantRecalibratorEngine - Finished iteration 0.
19:28:50.759 INFO  VariantRecalibratorEngine - Finished iteration 5.    Current change in mixture coefficients = 0.13274
19:28:50.811 INFO  VariantRecalibratorEngine - Finished iteration 10.   Current change in mixture coefficients = 0.01978
19:28:50.844 INFO  VariantRecalibratorEngine - Finished iteration 15.   Current change in mixture coefficients = 0.00160
19:28:50.844 INFO  VariantRecalibratorEngine - Convergence after 15 iterations!
19:28:50.856 INFO  VariantRecalibratorEngine - Evaluating full set of 23560 variants...
19:28:50.983 INFO  VariantDataManager - Selected worst 231 scoring variants --> variants with LOD <= -5.0000.
19:28:50.983 INFO  GaussianMixtureModel - Initializing model with 100 k-means iterations...
19:28:50.984 INFO  VariantRecalibratorEngine - Finished iteration 0.
19:28:50.985 INFO  VariantRecalibratorEngine - Finished iteration 5.    Current change in mixture coefficients = 0.00158
19:28:50.985 INFO  VariantRecalibratorEngine - Convergence after 5 iterations!
19:28:50.986 INFO  VariantRecalibratorEngine - Evaluating full set of 23560 variants...
19:28:51.079 INFO  TrancheManager - Finding 7 tranches for 23560 variants
19:28:51.091 INFO  TrancheManager -   TruthSensitivityTranche threshold 90.00 => selection metric threshold 0.100
19:28:51.097 INFO  TrancheManager -   Found tranche for 90.000: 0.100 threshold starting with variant 8211; running score is 0.100 
19:28:51.097 INFO  TrancheManager -   TruthSensitivityTranche is TruthSensitivityTranche targetTruthSensitivity=90.00 minVQSLod=4.2772 known=(0 @ 0.0000) novel=(15349 @ 2.1687) truthSites(9229 accessible, 8306 called), name=anonymous]
19:28:51.097 INFO  TrancheManager -   TruthSensitivityTranche threshold 93.00 => selection metric threshold 0.070
19:28:51.100 INFO  TrancheManager -   Found tranche for 93.000: 0.070 threshold starting with variant 4029; running score is 0.070 
19:28:51.100 INFO  TrancheManager -   TruthSensitivityTranche is TruthSensitivityTranche targetTruthSensitivity=93.00 minVQSLod=3.3913 known=(0 @ 0.0000) novel=(19531 @ 2.1433) truthSites(9229 accessible, 8582 called), name=anonymous]
19:28:51.100 INFO  TrancheManager -   TruthSensitivityTranche threshold 95.00 => selection metric threshold 0.050
19:28:51.102 INFO  TrancheManager -   Found tranche for 95.000: 0.050 threshold starting with variant 1872; running score is 0.050 
19:28:51.102 INFO  TrancheManager -   TruthSensitivityTranche is TruthSensitivityTranche targetTruthSensitivity=95.00 minVQSLod=2.0420 known=(0 @ 0.0000) novel=(21688 @ 2.1259) truthSites(9229 accessible, 8767 called), name=anonymous]
19:28:51.102 INFO  TrancheManager -   TruthSensitivityTranche threshold 97.00 => selection metric threshold 0.030
19:28:51.104 INFO  TrancheManager -   Found tranche for 97.000: 0.030 threshold starting with variant 1439; running score is 0.030 
19:28:51.104 INFO  TrancheManager -   TruthSensitivityTranche is TruthSensitivityTranche targetTruthSensitivity=97.00 minVQSLod=1.4785 known=(0 @ 0.0000) novel=(22121 @ 2.1295) truthSites(9229 accessible, 8952 called), name=anonymous]
19:28:51.104 INFO  TrancheManager -   TruthSensitivityTranche threshold 98.00 => selection metric threshold 0.020
19:28:51.105 INFO  TrancheManager -   Found tranche for 98.000: 0.020 threshold starting with variant 1213; running score is 0.020 
19:28:51.105 INFO  TrancheManager -   TruthSensitivityTranche is TruthSensitivityTranche targetTruthSensitivity=98.00 minVQSLod=1.0573 known=(0 @ 0.0000) novel=(22347 @ 2.1340) truthSites(9229 accessible, 9044 called), name=anonymous]
19:28:51.105 INFO  TrancheManager -   TruthSensitivityTranche threshold 99.90 => selection metric threshold 0.001
19:28:51.108 INFO  TrancheManager -   Found tranche for 99.900: 0.001 threshold starting with variant 313; running score is 0.001 
19:28:51.108 INFO  TrancheManager -   TruthSensitivityTranche is TruthSensitivityTranche targetTruthSensitivity=99.90 minVQSLod=-1.7435 known=(0 @ 0.0000) novel=(23247 @ 2.1319) truthSites(9229 accessible, 9219 called), name=anonymous]
19:28:51.108 INFO  TrancheManager -   TruthSensitivityTranche threshold 100.00 => selection metric threshold 0.000
19:28:51.110 INFO  TrancheManager -   Found tranche for 100.000: 0.000 threshold starting with variant 0; running score is 0.000 
19:28:51.110 INFO  TrancheManager -   TruthSensitivityTranche is TruthSensitivityTranche targetTruthSensitivity=100.00 minVQSLod=-38.0002 known=(0 @ 0.0000) novel=(23560 @ 2.1196) truthSites(9229 accessible, 9229 called), name=anonymous]
19:28:51.111 INFO  VariantRecalibrator - Writing out recalibration table...
19:28:51.355 INFO  VariantRecalibrator - Writing out visualization Rscript file...
19:28:51.363 INFO  VariantRecalibrator - Building QD x SOR plot...
19:28:51.365 INFO  VariantRecalibratorEngine - Evaluating full set of 3600 variants...
19:28:51.424 INFO  VariantRecalibratorEngine - Evaluating full set of 3600 variants...
19:28:51.486 INFO  VariantRecalibrator - Building QD x ReadPosRankSum plot...
19:28:51.488 INFO  VariantRecalibratorEngine - Evaluating full set of 3600 variants...
19:28:51.529 INFO  VariantRecalibratorEngine - Evaluating full set of 3600 variants...
19:28:51.586 INFO  VariantRecalibrator - Building QD x FS plot...
19:28:51.587 INFO  VariantRecalibratorEngine - Evaluating full set of 3600 variants...
19:28:51.623 INFO  VariantRecalibratorEngine - Evaluating full set of 3600 variants...
19:28:51.672 INFO  VariantRecalibrator - Building QD x null plot...
19:28:51.756 INFO  VariantRecalibrator - Shutting down engine
[February 17, 2023 at 7:28:51 PM IRST] org.broadinstitute.hellbender.tools.walkers.vqsr.VariantRecalibrator done. Elapsed time: 0.04 minutes.
Runtime.totalMemory()=386924544
java.lang.ArrayIndexOutOfBoundsException: Index 4 out of bounds for length 4
    at org.broadinstitute.hellbender.tools.walkers.vqsr.VariantRecalibrator.createVisualizationScript(VariantRecalibrator.java:1067)
    at org.broadinstitute.hellbender.tools.walkers.vqsr.VariantRecalibrator.onTraversalSuccess(VariantRecalibrator.java:710)
    at org.broadinstitute.hellbender.engine.GATKTool.doWork(GATKTool.java:1099)
    at org.broadinstitute.hellbender.cmdline.CommandLineProgram.runTool(CommandLineProgram.java:140)
    at org.broadinstitute.hellbender.cmdline.CommandLineProgram.instanceMainPostParseArgs(CommandLineProgram.java:192)
    at org.broadinstitute.hellbender.cmdline.CommandLineProgram.instanceMain(CommandLineProgram.java:211)
    at org.broadinstitute.hellbender.Main.runCommandLineProgram(Main.java:160)
    at org.broadinstitute.hellbender.Main.mainEntry(Main.java:203)
    at org.broadinstitute.hellbender.Main.main(Main.java:289)
WGS WES VQSR • 1.5k views
ADD COMMENT
0
Entering edit mode

Can you try skipping the --max-gaussians parameter and see if that works?

ADD REPLY
0
Entering edit mode

I tested it , I did not make any changes in my result.

ADD REPLY
0
Entering edit mode

OK, maybe someone more familiar with GATK can help - I don't understand what's happening here.

ADD REPLY
0
Entering edit mode

I have asked it on GATK forum too. If I get any responses from them, I will share it with you too.

ADD REPLY

Login before adding your answer.

Traffic: 2500 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6