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)
Can you try skipping the
--max-gaussians
parameter and see if that works?I tested it , I did not make any changes in my result.
OK, maybe someone more familiar with GATK can help - I don't understand what's happening here.
I have asked it on GATK forum too. If I get any responses from them, I will share it with you too.