[raonyguimaraes]

I rerunned all my analysis several times since the alignment using BWA -I option till UnifiedGenotyper and I'm still getting this output from UnifiedGenotyper:

GenomeAnalysisTK.jar -T UnifiedGenotyper -l INFO -I output/exome.real.dedup.recal.bam -R ../.
./input/b37/human_g1k_v37.fasta -B:intervals,BED ../../input/bed/exome_plus10.merged.bed -B:dbsnp,VCF ../../input/dbsnp/dbsnp-134.vcf -glm BOTH -stand_call_
conf 50.0 -stand_emit_conf 20.0 -dcov 300 -A AlleleBalance -A DepthOfCoverage -A FisherStrand -o output/exome.raw.vcf -log logs/gatk/UnifiedGenotyper.log

Visited bases 3102559836
Callable bases 2864301370
Confidently called bases 1142400
% callable bases of all loci 92.321
% confidently called bases of all loci 0.037
% confidently called bases of callable loci 0.040
Actual calls made 350263

Since it's the same number of variants I'm getting from DNANexus I'm starting to believe this is the right number before filtering the variants ...

I decided to download the BED File from UCSC Table Browser as suggested on the manual using "Exon plus 10bp", and tried to use this file with UnifiedGenotyper since I was using a bedfile from "SeqCap EZ Human Exome Library v2.0".

Running this part again I was receiving a message saying that there were overlaps on the intervals, so decided to use bedtools to merge the intervals with the command:

mergeBed -i exome_plus10.bed > exome_plus10.merged.bed

Does anyone else had to do it ?

After doing this and trying again i'm receiving the message:

exome_plus10.merged.bed and reference have incompatible contigs: No overlapping contigs found.
exome_plus10.merged.bed contigs = [chr1, chr10, chr11, chr12, chr13,
chr14, chr15, chr16, chr17, chr17_ctg5_hap1, chr17_gl000205_random, chr18, chr19, chr19_gl000209_random, chr1_gl000191_random, chr2, chr20, chr21, chr22, ch
r3, chr4, chr4_ctg9_hap1, chr4_gl000193_random, chr4_gl000194_random, chr5, chr6, chr6_apd_hap1, chr6_cox_hap2, chr6_dbb_hap3, chr6_mann_hap4, chr6_mcf_hap5
, chr6_qbl_hap6, chr6_ssto_hap7, chr7, chr7_gl000195_random, chr8, chr9, chrUn_gl000211, chrUn_gl000212, chrUn_gl000218, chrUn_gl000219, chrUn_gl000220, chr
Un_gl000222, chrUn_gl000223, chrUn_gl000228, chrX, chrY]
##### ERROR reference 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]

Sortbed didn't help me so I wrote an script in python to parse this bedfile and put everything as it should be chr1/1, chr2/2 and so on ... Does anyone else had to do the same?

I checked the quality of my reads with FASTQC and they looked ok, so I didn't do any clean on my reads before using BWA->GATK.

What I could use to clean my illumina reads ? NGS Backbone, SeqClean, CleanSeq, Prinseq, FastQX ? Does anyone improved the number of calls by doing it ?

For the Variant Quality Score Recalibration they suggest that "in order to achieve the best exome results one needs to use an exome SNP callset with at least 30 samples."

Does anyone tried to merge other exomes and got better results from it ?

Still looking for my 20k variants

[emilyjia2000]

Thanks both of you. Now it works.

[raonyguimaraes]

Please use everything from resource bundle ... http://www.broadinstitute.org/gsa/wi...esource_bundle

[Heisman]

Emily,

Make sure the various chromosome headers are in the same order in all files you are using (and make sure they match). They should also be in the same order as the way your files are sorted in whatever way you sort, I think.

[emilyjia2000]

Thanks raonyguimaraes, you are right. After I put dictionary file and reference in the same directory, it works in that way, but another error message came out:
##### ERROR ------------------------------------------------------------------------------------------
##### ERROR A USER ERROR has occurred (version 1.2-26-g43b0c98):
##### ERROR The invalid arguments or inputs must be corrected before the GATK can proceed
##### ERROR Please do not post this error to the GATK forum
##### ERROR
##### ERROR See the documentation (rerun with -h) for this tool to view allowable command-line arguments.
##### ERROR Visit our wiki for extensive documentation http://www.broadinstitute.org/gsa/wiki
##### ERROR Visit our forum to view answers to commonly asked questions http://getsatisfaction.com/gsa
##### ERROR
##### ERROR MESSAGE: Input files reads and reference have incompatible contigs: Order of contigs differences, which is unsafe.
##### ERROR reads contigs = [chrM, chr1, chr2, chr3, chr4, chr5, chr6, chr7, chr8, chr9, chr10, chr11, chr12, chr13, chr14, chr15, chr16, chr17, chr18, chr19, chr20, chr21, chr22, chrX, chrY]
##### ERROR reference contigs = [chr1, chr2, chr3, chr4, chr5, chr6, chr7, chrX, chr8, chr9, chr10, chr11, chr12, chr13, chr14, chr15, chr16, chr17, chr18, chr20, chrY, chr19, chr22, chr21, chr6_ssto_hap7, chr6_mcf_hap5, chr6_cox_hap2, chr6_mann_hap4, chr6_apd_hap1, chr6_qbl_hap6, chr6_dbb_hap3, chr17_ctg5_hap1, chr4_ctg9_hap1, chr1_gl000192_random, chrUn_gl000225, chr4_gl000194_random, chr4_gl000193_random, chr9_gl000200_random, chrUn_gl000222, chrUn_gl000212, chr7_gl000195_random, chrUn_gl000223, chrUn_gl000224, chrUn_gl000219, chr17_gl000205_random, chrUn_gl000215, chrUn_gl000216, chrUn_gl000217, chr9_gl000199_random, chrUn_gl000211, chrUn_gl000213, chrUn_gl000220, chrUn_gl000218, chr19_gl000209_random, chrUn_gl000221, chrUn_gl000214, chrUn_gl000228, chrUn_gl000227, chr1_gl000191_random, chr19_gl000208_random, chr9_gl000198_random, chr17_gl000204_random, chrUn_gl000233, chrUn_gl000237, chrUn_gl000230, chrUn_gl000242, chrUn_gl000243, chrUn_gl000241, chrUn_gl000236, chrUn_gl000240, chr17_gl000206_random, chrUn_gl000232, chrUn_gl000234, chr11_gl000202_random, chrUn_gl000238, chrUn_gl000244, chrUn_gl000248, chr8_gl000196_random, chrUn_gl000249, chrUn_gl000246, chr17_gl000203_random, chr8_gl000197_random, chrUn_gl000245, chrUn_gl000247, chr9_gl000201_random, chrUn_gl000235, chrUn_gl000239, chr21_gl000210_random, chrUn_gl000231, chrUn_gl000229, chrM, chrUn_gl000226, chr18_gl000207_random]
##### ERROR ------------------------------------------------------------------------------------------

The chromosome in BAM file has been reordered, should I reorder reference file?
Thanks a lot

[raonyguimaraes]

If I'm not wrong ucsc.hg19.dict should be on the same directory of the reference file.

[emilyjia2000]

Thanks for all of your quick response. I used the command line:

java -Xmx4g -jar /path/to/GenomeAnalysisTK.jar -T RealignerTargetCreator -R hg19.fa -o output.interval -I /path/to/reorder_dedup.bam

I already copied the ucsc.hg19.dict in the same directory.

When I run this command, the error message:

##### ERROR ------------------------------------------------------------------------------------------
##### ERROR A USER ERROR has occurred (version 1.2-26-g43b0c98):
##### ERROR The invalid arguments or inputs must be corrected before the GATK can proceed
##### ERROR Please do not post this error to the GATK forum
##### ERROR
##### ERROR See the documentation (rerun with -h) for this tool to view allowable command-line arguments.
##### ERROR Visit our wiki for extensive documentation http://www.broadinstitute.org/gsa/wiki
##### ERROR Visit our forum to view answers to commonly asked questions http://getsatisfaction.com/gsa
##### ERROR
##### ERROR MESSAGE: Invalid command line: Failed to load reference dictionary
##### ERROR ------------------------------------------------------------------------------------------

[liu_xt005]

A quick comment: I would suggest people to do the SAM -> BAM conversion using Picard.
In fact, samtools sort generates 79 temporary files for one sample.

[Heisman]

When I tried to use GATK to do the local realignment according to ulz_peter's instruction, one error message occurred: Invalid command line: Failed to load reference dictionary. Could anybody let me know where to get this reference dictionary? how to use it in the command line?

Thanks in advance.

Look here: http://www.broadinstitute.org/gsa/wi...ference_genome

[raonyguimaraes]

Post your command line here... and the full error message

[emilyjia2000]

reference dictionary

When I tried to use GATK to do the local realignment according to ulz_peter's instruction, one error message occurred: Invalid command line: Failed to load reference dictionary. Could anybody let me know where to get this reference dictionary? how to use it in the command line?

Thanks in advance.

[ulz_peter]

You're absolutely right on the first one.
For the sampe -r option: I'd like to keep it in the bwa part as another step handling the BAM file takes a lot of time than just doing it once... but I could add that as an alternative

[MQ-BCBB]

Great resource!!!
Just a few comments:
In the picard/MarkDuplicates.jar, the option 'CREATE_INDEX=true' should be added.
With respect to adding read group information, instead of using the bwa sampe -r option, picard AddOrReplaceReadGroups.jar is an easier way to go as it tells you which options are required. Thanks for sharing!

[liu_xt005]

-L option problem solved

What is the error message when you specify the -L argument?

I actually stopped using the Variant quality Score recalibration as it often did not work out for me (I never work on more than 2 exomes at a time).

I out the version withouth the recalibration on the SEQanswers Wiki/How-To section. You may have a look there, as I will update that in the future and stop uploading newer versions of the PDF file...

Dear ulz_peter and raonyguimaraes,

Thanks VERY MUCH to both of you.
The problem seems to be solved by removing the random and hap intervals from the .bed file.

ulz_peter,
raonyguimaraes posted a similar pipeline by Gayle Philip.
SNPs and Indels are recalibrated/filtered separately and combined after.
I am trying that, and think that it is a good idea to exclude Indels from Gaussian models.

[mirabilia]

I didn't find it yet, but there was a statement on the GATK homepage that the options descirbed there (which are basically pretty mcuh the same as I use) only work for diploid genomes and expected shifts of allele frequency must be adressed. So the question is: what are you planning to do: find rare alleles within some strains, sequence a genetically homogeneous strain...

The first part of the analysis, until the local realignment around indels, is of course suitable for all genomes I have. Quality score recalibration and SNP calling need customization and I'm thinking on... Actually, I'd like to catch rare alleles in sequenced bacterial populations, like repeats, tandem repeats and structural variants, but I suspect it's a really difficult task to discriminate between real rare variations and noise introduced by such data.
Obviously any kind of suggestion, it's really appreciate!