Attention:
The frontline support team will be slow on the forum because we are occupied with the GATK Workshop on March 21st and 22nd 2019. We will be back and more available to answer questions on the forum on March 25th 2019.

Haplotype caller-generated VCF cannot be validated

fojackson8fojackson8 Cambridge, UKMember

Hi, I have VCFs generated by Haplotype Caller and Unified Genotyper on the NA12878 genome. Running Variant Recalibrator brings up errors in the VCF, so I'm attempting to understand them using GATK Validate Variants.

Should the Validate Variants tool run through a Haplotype caller-generated VCF no problem? It's fine with the Unified Genotyper VCF, but the Haplotype Caller VCF returns a number of errors. Specifically:

ERROR MESSAGE: File /home/felix/Genomics/VCFs/HCoutput1.vcf fails strict validation: one or more of the ALT allele(s) for the record at position 1:825826 are not observed at all in the sample genotypes,

Then if I rerun with --validationTypeToExclude ALLELES, I get:

ERROR MESSAGE: File /home/felix/Genomics/VCFs/HCoutput1.vcf fails strict validation: the Allele Number (AN) tag is incorrect for the record at position 1:825826, 2 vs. 0

I'm not sure how this has come about, because the Unified Genotyper-generated VCF on the exact same input BAM has no errors.

Any pointers would be greatly appreciated.

Answers

  • SheilaSheila Broad InstituteMember, Broadie, Moderator admin

    @fojackson8
    Hi,

    Can you please post the exact command and version you used to get the invalid VCF? Please also post the record at position 1:825826.

    Thanks,
    Sheila

  • fojackson8fojackson8 Cambridge, UKMember

    Sure. Command:

    java -jar GenomeAnalysisTK.jar -T ValidateVariants -V ~/Genomics/VCFs/HCoutput1.vcf -R ~/Genomics/_1_human_g1k_v37_decoy.fasta

    The version is GATK v3.7-0-gcfedb67.

    The record is:

    1 825826 rs139918331 C CA 625.74 . AC=1;AF=0.500;AN=2;BaseQRankSum=0.628;ClippingRankSum=0.000;DB;DP=41;ExcessHet=3.0103;FS=3.142;MLEAC=1;MLEAF=0.500;MQ=56.91;MQRankSum=5.478;QD=20.86;ReadPosRankSum=0.624;SOR=0.105 GT ./.

    Thanks,
    Felix

  • SheilaSheila Broad InstituteMember, Broadie, Moderator admin

    @fojackson8
    Hi Felix,

    Hmm. It looks like the issue is AN=2, but the genotype has 0 called alleles. I think AN should be 0. Can you confirm you used the latest version of HaplotypeCaller to create the VCF? Please also post the exact command you ran to generate the VCF. I may need you to submit a bug report.

    Thanks,
    Sheila

  • fojackson8fojackson8 Cambridge, UKMember

    Hi, yes the GATK version is 3.70, which I think is the latest version. I'm running haplotype caller on a genomics platform using a docker image, so the command is a bit different, but here it is anyway:

    java -Xmx16384M -jar /gatk/GenomeAnalysisTK.jar -T HaplotypeCaller -R /sbgenomics/Projects/4c32289d-6d41-4bd3-9d60-bbce2cf364f9/human_g1k_v37_decoy.fasta -I /sbgenomics/Projects/4c32289d-6d41-4bd3-9d60-bbce2cf364f9/_2_HG001-NA12878-50x.base_recalibrated.bam -o output.vcf --dbsnp /sbgenomics/Projects/4c32289d-6d41-4bd3-9d60-bbce2cf364f9/dbsnp_137.b37.vcf -nt 1 -nct 4

    Thanks,

    Felix

  • Geraldine_VdAuweraGeraldine_VdAuwera Cambridge, MAMember, Administrator, Broadie admin
    What do you mean by "the command is a bit different"? Is this what gets run or does the pipeline generate a command that's different from this? And does the pipeline you're using do any kind of processing after calling with HaplotypeCaller? If you post the header of the vcf (minus SQ lines) we can try to dig into that.
  • fojackson8fojackson8 Cambridge, UKMember
    edited March 2017

    Ah no, it's the command is as above - if the above command looks alright then that's good. No it doesn't do any further processing, I'm trying to run Variant Recalibrator on it, but it GATK tools keep throwing errors with it.

    Here's the header:

    fileformat=VCFv4.2

    FILTER=<ID=LowQual,Description="Low quality">

    FORMAT=<ID=AD,Number=R,Type=Integer,Description="Allelic depths for the ref and alt alleles in the order listed">

    FORMAT=<ID=DP,Number=1,Type=Integer,Description="Approximate read depth (reads with MQ=255 or with bad mates are filtered)">

    FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Genotype Quality">

    FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">

    FORMAT=<ID=PL,Number=G,Type=Integer,Description="Normalized, Phred-scaled likelihoods for genotypes as defined in the VCF specification">

    GATKCommandLine.HaplotypeCaller=<ID=HaplotypeCaller,Version=3.7-0-gcfedb67,Date="Thu Feb 16 11:53:30 UTC 2017",Epoch=1487246010641,CommandLineOptions="analysis_type=HaplotypeCaller input_file=[/sbgenomics/Projects/4c32289d-6d41-4bd3-9d60-bbce2cf364f9/_2_HG001-NA12878-50x.base_recalibrated.bam] showFullBamList=false read_buffer_size=null read_filter=[] disable_read_filter=[] intervals=null excludeIntervals=null interval_set_rule=UNION interval_merging=ALL interval_padding=0 reference_sequence=/sbgenomics/Projects/4c32289d-6d41-4bd3-9d60-bbce2cf364f9/human_g1k_v37_decoy.fasta nonDeterministicRandomSeed=false disableDithering=false maxRuntime=-1 maxRuntimeUnits=MINUTES downsampling_type=BY_SAMPLE downsample_to_fraction=null downsample_to_coverage=500 baq=OFF baqGapOpenPenalty=40.0 refactor_NDN_cigar_string=false fix_misencoded_quality_scores=false allow_potentially_misencoded_quality_scores=false useOriginalQualities=false defaultBaseQualities=-1 performanceLog=null BQSR=null quantize_quals=0 static_quantized_quals=null round_down_quantized=false disable_indel_quals=false emit_original_quals=false preserve_qscores_less_than=6 globalQScorePrior=-1.0 secondsBetweenProgressUpdates=10 validation_strictness=SILENT remove_program_records=false keep_program_records=false sample_rename_mapping_file=null unsafe=null disable_auto_index_creation_and_locking_when_reading_rods=false no_cmdline_in_header=false sites_only=false never_trim_vcf_format_field=false bcf=false bam_compression=null simplifyBAM=false disable_bam_indexing=false generate_md5=false num_threads=1 num_cpu_threads_per_data_thread=4 num_io_threads=0 monitorThreadEfficiency=false num_bam_file_handles=null read_group_black_list=null pedigree=[] pedigreeString=[] pedigreeValidationType=STRICT allow_intervals_with_unindexed_bam=false generateShadowBCF=false variant_index_type=DYNAMIC_SEEK variant_index_parameter=-1 reference_window_stop=0 phone_home= gatk_key=null tag=NA logging_level=INFO log_to_file=null help=false version=false out=/sbgenomics/Projects/4c32289d-6d41-4bd3-9d60-bbce2cf364f9/workspace/981b02b3-6e37-4b3f-ac0e-ef4783a89481/gatk-haplotype-caller/output.vcf likelihoodCalculationEngine=PairHMM heterogeneousKmerSizeResolution=COMBO_MIN dbsnp=(RodBinding name=dbsnp source=/sbgenomics/Projects/4c32289d-6d41-4bd3-9d60-bbce2cf364f9/dbsnp_137.b37.vcf) dontTrimActiveRegions=false maxDiscARExtension=25 maxGGAARExtension=300 paddingAroundIndels=150 paddingAroundSNPs=20 comp=[] annotation=[] excludeAnnotation=[] group=[StandardAnnotation, StandardHCAnnotation] debug=false useFilteredReadsForAnnotations=false emitRefConfidence=NONE bamOutput=null bamWriterType=CALLED_HAPLOTYPES emitDroppedReads=false disableOptimizations=false annotateNDA=false useNewAFCalculator=false heterozygosity=0.001 indel_heterozygosity=1.25E-4 heterozygosity_stdev=0.01 standard_min_confidence_threshold_for_calling=10.0 standard_min_confidence_threshold_for_emitting=30.0 max_alternate_alleles=6 max_genotype_count=1024 max_num_PL_values=100 input_prior=[] sample_ploidy=2 genotyping_mode=DISCOVERY alleles=(RodBinding name= source=UNBOUND) contamination_fraction_to_filter=0.0 contamination_fraction_per_sample_file=null p_nonref_model=null exactcallslog=null output_mode=EMIT_VARIANTS_ONLY allSitePLs=false gcpHMM=10 pair_hmm_implementation=VECTOR_LOGLESS_CACHING pair_hmm_sub_implementation=ENABLE_ALL always_load_vector_logless_PairHMM_lib=false phredScaledGlobalReadMismappingRate=45 noFpga=false sample_name=null kmerSize=[10, 25] dontIncreaseKmerSizesForCycles=false allowNonUniqueKmersInRef=false numPruningSamples=1 recoverDanglingHeads=false doNotRecoverDanglingBranches=false minDanglingBranchLength=4 consensus=false maxNumHaplotypesInPopulation=128 errorCorrectKmers=false minPruning=2 debugGraphTransformations=false allowCyclesInKmerGraphToGeneratePaths=false graphOutput=null kmerLengthForReadErrorCorrection=25 minObservationsForKmerToBeSolid=20 GVCFGQBands=[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 70, 80, 90, 99] indelSizeToEliminateInRefModel=10 min_base_quality_score=10 includeUmappedReads=false useAllelesTrigger=false doNotRunPhysicalPhasing=true keepRG=null justDetermineActiveRegions=false dontGenotype=false dontUseSoftClippedBases=false captureAssemblyFailureBAM=false errorCorrectReads=false pcr_indel_model=CONSERVATIVE maxReadsInRegionPerSample=10000 minReadsPerAlignmentStart=10 mergeVariantsViaLD=false activityProfileOut=null activeRegionOut=null activeRegionIn=null activeRegionExtension=null forceActive=false activeRegionMaxSize=null bandPassSigma=null maxReadsInMemoryPerSample=30000 maxTotalReadsInMemory=10000000 maxProbPropagationDistance=50 activeProbabilityThreshold=0.002 min_mapping_quality_score=20 filter_reads_with_N_cigar=false filter_mismatching_base_and_quals=false filter_bases_not_stored=false">

    INFO=<ID=AC,Number=A,Type=Integer,Description="Allele count in genotypes, for each ALT allele, in the same order as listed">

    INFO=<ID=AF,Number=A,Type=Float,Description="Allele Frequency, for each ALT allele, in the same order as listed">

    INFO=<ID=AN,Number=1,Type=Integer,Description="Total number of alleles in called genotypes">

    INFO=<ID=BaseQRankSum,Number=1,Type=Float,Description="Z-score from Wilcoxon rank sum test of Alt Vs. Ref base qualities">

    INFO=<ID=ClippingRankSum,Number=1,Type=Float,Description="Z-score From Wilcoxon rank sum test of Alt vs. Ref number of hard clipped bases">

    INFO=<ID=DB,Number=0,Type=Flag,Description="dbSNP Membership">

    INFO=<ID=DP,Number=1,Type=Integer,Description="Approximate read depth; some reads may have been filtered">

    INFO=<ID=DS,Number=0,Type=Flag,Description="Were any of the samples downsampled?">

    INFO=<ID=ExcessHet,Number=1,Type=Float,Description="Phred-scaled p-value for exact test of excess heterozygosity">

    INFO=<ID=FS,Number=1,Type=Float,Description="Phred-scaled p-value using Fisher's exact test to detect strand bias">

    INFO=<ID=HaplotypeScore,Number=1,Type=Float,Description="Consistency of the site with at most two segregating haplotypes">

    INFO=<ID=InbreedingCoeff,Number=1,Type=Float,Description="Inbreeding coefficient as estimated from the genotype likelihoods per-sample when compared against the Hardy-Weinberg expectation">

    INFO=<ID=MLEAC,Number=A,Type=Integer,Description="Maximum likelihood expectation (MLE) for the allele counts (not necessarily the same as the AC), for each ALT allele, in the same order as listed">

    INFO=<ID=MLEAF,Number=A,Type=Float,Description="Maximum likelihood expectation (MLE) for the allele frequency (not necessarily the same as the AF), for each ALT allele, in the same order as listed">

    INFO=<ID=MQ,Number=1,Type=Float,Description="RMS Mapping Quality">

    INFO=<ID=MQRankSum,Number=1,Type=Float,Description="Z-score From Wilcoxon rank sum test of Alt vs. Ref read mapping qualities">

    INFO=<ID=QD,Number=1,Type=Float,Description="Variant Confidence/Quality by Depth">

    INFO=<ID=ReadPosRankSum,Number=1,Type=Float,Description="Z-score from Wilcoxon rank sum test of Alt vs. Ref read position bias">

    INFO=<ID=SOR,Number=1,Type=Float,Description="Symmetric Odds Ratio of 2x2 contingency table to detect strand bias">

    contig=<ID=1,length=249250621,assembly=b37>

    contig=<ID=2,length=243199373,assembly=b37>

    contig=<ID=3,length=198022430,assembly=b37>

    contig=<ID=4,length=191154276,assembly=b37>

    contig=<ID=5,length=180915260,assembly=b37>

    contig=<ID=6,length=171115067,assembly=b37>

    contig=<ID=7,length=159138663,assembly=b37>

    contig=<ID=8,length=146364022,assembly=b37>

    contig=<ID=9,length=141213431,assembly=b37>

    contig=<ID=10,length=135534747,assembly=b37>

    contig=<ID=11,length=135006516,assembly=b37>

    contig=<ID=12,length=133851895,assembly=b37>

    contig=<ID=13,length=115169878,assembly=b37>

    contig=<ID=14,length=107349540,assembly=b37>

    contig=<ID=15,length=102531392,assembly=b37>

    contig=<ID=16,length=90354753,assembly=b37>

    contig=<ID=17,length=81195210,assembly=b37>

    contig=<ID=18,length=78077248,assembly=b37>

    contig=<ID=19,length=59128983,assembly=b37>

    contig=<ID=20,length=63025520,assembly=b37>

    contig=<ID=21,length=48129895,assembly=b37>

    contig=<ID=22,length=51304566,assembly=b37>

    contig=<ID=X,length=155270560,assembly=b37>

    contig=<ID=Y,length=59373566,assembly=b37>

    contig=<ID=MT,length=16569,assembly=b37>

    contig=<ID=GL000207.1,length=4262,assembly=b37>

    contig=<ID=GL000226.1,length=15008,assembly=b37>

    contig=<ID=GL000229.1,length=19913,assembly=b37>

    contig=<ID=GL000231.1,length=27386,assembly=b37>

    contig=<ID=GL000210.1,length=27682,assembly=b37>

    contig=<ID=GL000239.1,length=33824,assembly=b37>

    contig=<ID=GL000235.1,length=34474,assembly=b37>

    contig=<ID=GL000201.1,length=36148,assembly=b37>

    contig=<ID=GL000247.1,length=36422,assembly=b37>

    contig=<ID=GL000245.1,length=36651,assembly=b37>

    contig=<ID=GL000197.1,length=37175,assembly=b37>

    contig=<ID=GL000203.1,length=37498,assembly=b37>

    contig=<ID=GL000246.1,length=38154,assembly=b37>

    contig=<ID=GL000249.1,length=38502,assembly=b37>

    contig=<ID=GL000196.1,length=38914,assembly=b37>

    contig=<ID=GL000248.1,length=39786,assembly=b37>

    contig=<ID=GL000244.1,length=39929,assembly=b37>

    contig=<ID=GL000238.1,length=39939,assembly=b37>

    contig=<ID=GL000202.1,length=40103,assembly=b37>

    contig=<ID=GL000234.1,length=40531,assembly=b37>

    contig=<ID=GL000232.1,length=40652,assembly=b37>

    contig=<ID=GL000206.1,length=41001,assembly=b37>

    contig=<ID=GL000240.1,length=41933,assembly=b37>

    contig=<ID=GL000236.1,length=41934,assembly=b37>

    contig=<ID=GL000241.1,length=42152,assembly=b37>

    contig=<ID=GL000243.1,length=43341,assembly=b37>

    contig=<ID=GL000242.1,length=43523,assembly=b37>

    contig=<ID=GL000230.1,length=43691,assembly=b37>

    contig=<ID=GL000237.1,length=45867,assembly=b37>

    contig=<ID=GL000233.1,length=45941,assembly=b37>

    contig=<ID=GL000204.1,length=81310,assembly=b37>

    contig=<ID=GL000198.1,length=90085,assembly=b37>

    contig=<ID=GL000208.1,length=92689,assembly=b37>

    contig=<ID=GL000191.1,length=106433,assembly=b37>

    contig=<ID=GL000227.1,length=128374,assembly=b37>

    contig=<ID=GL000228.1,length=129120,assembly=b37>

    contig=<ID=GL000214.1,length=137718,assembly=b37>

    contig=<ID=GL000221.1,length=155397,assembly=b37>

    contig=<ID=GL000209.1,length=159169,assembly=b37>

    contig=<ID=GL000218.1,length=161147,assembly=b37>

    contig=<ID=GL000220.1,length=161802,assembly=b37>

    contig=<ID=GL000213.1,length=164239,assembly=b37>

    contig=<ID=GL000211.1,length=166566,assembly=b37>

    contig=<ID=GL000199.1,length=169874,assembly=b37>

    contig=<ID=GL000217.1,length=172149,assembly=b37>

    contig=<ID=GL000216.1,length=172294,assembly=b37>

    contig=<ID=GL000215.1,length=172545,assembly=b37>

    contig=<ID=GL000205.1,length=174588,assembly=b37>

    contig=<ID=GL000219.1,length=179198,assembly=b37>

    contig=<ID=GL000224.1,length=179693,assembly=b37>

    contig=<ID=GL000223.1,length=180455,assembly=b37>

    contig=<ID=GL000195.1,length=182896,assembly=b37>

    contig=<ID=GL000212.1,length=186858,assembly=b37>

    contig=<ID=GL000222.1,length=186861,assembly=b37>

    contig=<ID=GL000200.1,length=187035,assembly=b37>

    contig=<ID=GL000193.1,length=189789,assembly=b37>

    contig=<ID=GL000194.1,length=191469,assembly=b37>

    contig=<ID=GL000225.1,length=211173,assembly=b37>

    contig=<ID=GL000192.1,length=547496,assembly=b37>

    contig=<ID=NC_007605,length=171823,assembly=b37>

    contig=<ID=hs37d5,length=35477943,assembly=b37>

    reference=file:///sbgenomics/Projects/4c32289d-6d41-4bd3-9d60-bbce2cf364f9/human_g1k_v37_decoy.fasta

    CHROM POS ID REF ALT QUAL FILTER INFO FORMAT HG001

    Sorry for the poor formating, and thanks for your help so far, it's much appreciated.

  • SheilaSheila Broad InstituteMember, Broadie, Moderator admin

    @fojackson8
    Hi,

    If you can confirm this happens without -nct 4, then I will need you to submit a bug report. Instructions are here. Please just include ~1000 bases before and after the site that gives the incorrect AN with HaplotypeCaller.

    Thanks
    Sheila

Sign In or Register to comment.