HaplotypeCaller doesn't output the VCF
2
0
Entering edit mode
3 months ago
Flexogore • 0

Hello everyone!

Still struggling to implement variant calling with my BAM file with the following command:

java -Xmx10g -jar gatk-package-4.1.9.0-local.jar HaplotypeCaller -R illumina.fasta -I final.bam -O illumina.vcf

BAM has been previously indexed and identified:

@SQ SN:GL000222.1   LN:186861   AS:human_g1k_v37_modified.ndx
@SQ SN:GL000200.1   LN:187035   AS:human_g1k_v37_modified.ndx
@SQ SN:GL000193.1   LN:189789   AS:human_g1k_v37_modified.ndx
@SQ SN:GL000194.1   LN:191469   AS:human_g1k_v37_modified.ndx
@SQ SN:GL000225.1   LN:211173   AS:human_g1k_v37_modified.ndx
@SQ SN:GL000192.1   LN:547496   AS:human_g1k_v37_modified.ndx
@RG ID:Sample1  LB:lib1M    PL:illumina SM:1M   PU:unit1
@PG ID:samtools PN:samtools VN:1.10 CL:samtools view -h illumina.initial.bam
@PG ID:samtools.1   PN:samtools PP:samtools VN:1.10 CL:samtools view -bo your_fixed.bam -
@PG ID:samtools.2   PN:samtools PP:samtools.1   VN:1.10 CL:samtools view -H final.bam

The command above seems to be right, however instead of VCF I keep getting such an output listed below (a fragment):

A00128:303:HW5CTDSXX:1:1272:4743:8234/2, A00128:303:HW5CTDSXX:2:2130:26304:21746/2, 
A00128:303:HW5CTDSXX:3:2569:18394:18223/2, A00128:303:HW5CTDSXX:3:1532:10890:10081/1, 
A00128:303:HW5CTDSXX:2:2139:4869:15499/2, A00128:303:HW5CTDSXX:4:2641:19226:19445/2, 
A00128:303:HW5CTDSXX:3:1631:3088:28510/1, A00128:303:HW5CTDSXX:2:2628:13801:23296/2, 
A00128:303:HW5CTDSXX:3:1111:23357:28134/2, A00128:303:HW5CTDSXX:2:1546:11930:7905/2, 
A00128:303:HW5CTDSXX:2:2545:32479:36636/2, A00128:303:HW5CTDSXX:3:2368:16749:30874/2, 
A00128:303:HW5CTDSXX:1:1578:31991:1720/1, A00128:303:HW5CTDSXX:1:1609:25446:6136/1, 
A00128:303:HW5CTDSXX:2:1434:12029:24173/1, A00128:303:HW5CTDSXX:3:2109:28782:36464/2, 
A00128:303:HW5CTDSXX:4:2649:12671:18239/2, A00128:303:HW5CTDSXX:4:2106:27344:22482/1, 
A00128:303:HW5CTDSXX:3:1121:18945:16924/2, A00128:303:HW5CTDSXX:2:2348:4255:23234/2, 
A00128:303:HW5CTDSXX:1:2461:25066:3912/2, A00128:303:HW5CTDSXX:4:2455:13105:34867/2, 
A00128:303:HW5CTDSXX:3:2368:32289:1391/1, A00128:303:HW5CTDSXX:1:1255:14606:10942/2, 
A00128:303:HW5CTDSXX:4:2358:14895:3771/1, A00128:303:HW5CTDSXX:2:1232:21377:21856/1, 
A00128:303:HW5CTDSXX:3:2202:3875:20541/1, A00128:303:HW5CTDSXX:1:1173:16568:12493/2, 
A00128:303:HW5CTDSXX:3:2347:13349:1939/1, A00128:303:HW5CTDSXX:2:1633:20835:5979/2, 
A00128:303:HW5CTDSXX:1:1210:29496:2440/2, A00128:303:HW5CTDSXX:2:1678:6497:10144/2, 
A00128:303:HW5CTDSXX:2:2169:21377:6136/1, A00128:303:HW5CTDSXX:4:2367:31467:29152/1, 
A00128:303:HW5CTDSXX:1:1478:11460:4022/1, A00128:303:HW5CTDSXX:4:2242:6994:16736/2, 
A00128:303:HW5CTDSXX:1:2220:9426:34601/1, A00128:303:HW5CTDSXX:1:2209:12292:6277/1, 
A00128:303:HW5CTDSXX:4:2377:28248:3599/2, A00128:303:HW5CTDSXX:2:2667:18285:36605/1, 
A00128:303:HW5CTDSXX:3:1543:15591:28244/1, A00128:303:HW5CTDSXX:2:1367:31783:14481/1, 
A00128:303:HW5CTDSXX:2:1367:31783:14481/2, A00128:303:HW5CTDSXX:4:2320:27019:28995/1, 
A00128:303:HW5CTDSXX:4:2320:27019:28995/2]
reads contigs = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, X, Y, MT, GL000207.1, GL000226.1, GL000229.1, 
GL000231.1, GL000210.1, GL000239.1, GL000235.1, GL000201.1, GL000247.1, GL000245.1, GL000197.1, GL000203.1, GL000246.1, 
GL000249.1, GL000196.1, GL000248.1, GL000244.1, GL000238.1, GL000202.1, GL000234.1, GL000232.1, GL000206.1, GL000240.1, 
GL000236.1, GL000241.1, GL000243.1, GL000242.1, GL000230.1, GL000237.1, GL000233.1, GL000204.1, GL000198.1, GL000208.1, 
GL000191.1, GL000227.1, GL000228.1, GL000214.1, GL000221.1, GL000209.1, GL000218.1, GL000220.1, GL000213.1, GL000211.1, 
GL000199.1, GL000217.1, GL000216.1, GL000215.1, GL000205.1, GL000219.1, GL000224.1, GL000223.1, GL000195.1, GL000212.1, 
GL000222.1, GL000200.1, GL000193.1, GL000194.1, GL000225.1, GL000192.1]

***********************************************************************
Set the system property GATK_STACKTRACE_ON_USER_EXCEPTION (--java-options '- 
DGATK_STACKTRACE_ON_USER_EXCEPTION=true') to print the stack trace.

I don't get any error messages so I am wondering what's the possible reason. Do you have any ideas?

BAM VCF GATK HaplotypeCaller • 396 views
ADD COMMENT
1
Entering edit mode
3 months ago

The reference used for mapping the reads is not the same than the one you're using for haplotype caller.

e.g: Incompatible contigs No overlapping contigs found.; GATK BaseRecalibrator error: input files reference and features have incompatible contigs ; error using bqsr and incompatible contigs ; etc... etc...

ADD COMMENT
0
Entering edit mode

I've used FASTA extracted from BAM via samtools but seems it didn't work. Are there any other options to get the reference without having it itself?

ADD REPLY
0
Entering edit mode
3 months ago
guyho • 0

I do the following steps for variant calling (Tool used):

Starting with fastq files

Align to genome with (BWA-mem)

Sort (samtools sort)

Add read groups (AddOrReplaceReadGroups)

mark duplicates (MarkDuplicatesSpark )

base recalibration (BaseRecalibrator)

base recalibration (ApplyBQSR )

Variant calling (HaplotypeCaller)

Integrate vcfs (GenotypeGVCFs )

ADD COMMENT
0
Entering edit mode

this doesn't answer the question.

ADD REPLY
0
Entering edit mode

HaplotypeCaller assumes sorted bam with read-group tag. If these requirements are not fulfilled it might cause the strange behavior. The question did not specify how the bam was created.

ADD REPLY

Login before adding your answer.

Traffic: 2241 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