[raonyguimaraes]

Following ulz_peter's original doc, I have some problem when doing the SNP-calling.

java -Xmx4g -jar /path/GenomeAnalysisTK-1.1-35-ge253f6f/GenomeAnalysisTK.jar \
-glm BOTH \
-R hg18.fa \
-T UnifiedGenotyper \
-I myinput.marked.realigned.fixed.recal.bam \
-D dbsnp132_hg18.txt \
-o myoutput.snps.vcf \
-metrics snps.metrics \
-stand_call_conf 50.0 \
-stand_emit_conf 10.0 \
-dcov 1000 \
-A DepthOfCoverage \
-A AlleleBalance \
-L hg18_exonIntervals.bed

This "-L" option does not work.
I got the hg18_exonIntervals.bed from UCSC as ulz_peter's original doc shows.
I run the SNP-calling without the "-L" line.
Then the variant quality score recalibration step does not work, generating an empty output.tranches file.

Can somebody help me out? Thanks a lot.

Did you try with -B:targetIntervals,BED hg18_exonIntervals.bed?

By the way I couldn't figure out how to use this on version 1.2

[ulz_peter]

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...

[liu_xt005]

Problem with SNP-calling

Following ulz_peter's original doc, I have some problem when doing the SNP-calling.

java -Xmx4g -jar /path/GenomeAnalysisTK-1.1-35-ge253f6f/GenomeAnalysisTK.jar \
-glm BOTH \
-R hg18.fa \
-T UnifiedGenotyper \
-I myinput.marked.realigned.fixed.recal.bam \
-D dbsnp132_hg18.txt \
-o myoutput.snps.vcf \
-metrics snps.metrics \
-stand_call_conf 50.0 \
-stand_emit_conf 10.0 \
-dcov 1000 \
-A DepthOfCoverage \
-A AlleleBalance \
-L hg18_exonIntervals.bed

This "-L" option does not work.
I got the hg18_exonIntervals.bed from UCSC as ulz_peter's original doc shows.
I run the SNP-calling without the "-L" line.
Then the variant quality score recalibration step does not work, generating an empty output.tranches file.

Can somebody help me out? Thanks a lot.

[ulz_peter]

thanks a lot ulz_peter!
Could you please clarify which steps of your pipeline are specifically for diploid genomes in order I can customize for my purposes?

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...

[Heisman]

I am not affiliated with VAAST in any way, but I have used it extensively and absolutely love it. I don't want to derail this thread in any way but I can certainly answer questions about it.

[raonyguimaraes]

I think we should all give a try to VAAST as well

http://www.yandell-lab.org/software/vaast.html

A probabilistic disease-gene finder for personal genomes.
Yandell M, Huff C, Hu H, Singleton M, Moore B, Xing J, Jorde LB, Reese MG.
Source

Department of Human Genetics, Eccles Institute of Human Genetics, University of Utah and School of Medicine, Salt Lake City, UT 84112, USA. [email protected]

VAAST (the Variant Annotation, Analysis & Search Tool) is a probabilistic search tool for identifying damaged genes and their disease-causing variants in personal genome sequences. VAAST builds on existing amino acid substitution (AAS) and aggregative approaches to variant prioritization, combining elements of both into a single unified likelihood framework that allows users to identify damaged genes and deleterious variants with greater accuracy, and in an easy-to-use fashion. VAAST can score both coding and noncoding variants, evaluating the cumulative impact of both types of variants simultaneously. VAAST can identify rare variants causing rare genetic diseases, and it can also use both rare and common variants to identify genes responsible for common diseases. VAAST thus has a much greater scope of use than any existing methodology. Here we demonstrate its ability to identify damaged genes using small cohorts (n = 3) of unrelated individuals, wherein no two share the same deleterious variants, and for common, multigenic diseases using as few as 150 cases.

[Michael.James.Clark]

Thanks you ulz_peter and raonyguimaraes. I'm starting doing NGS quite soon, and being a newbie, this pipeline also seems very similar to what my instructor recommended me to do (BWA, SamTools/Varscan, and Annovar). She also said something about using Sift and Polyphen to predict the effect of the mutation on the gene functionality (continuous score that is benign below a certain threshold, and destructive above it). Anyone knows what I'm talking about?

Annovar can annotate with SIFT and Polyphen now.

[Orr Shomroni]

Thank you guys

Thank you ulz_peter and raonyguimaraes. I'm starting doing NGS quite soon, and being a newbie, this pipeline is quite helpful. Also it seems very similar to what my instructor recommended I should do (BWA, SamTools/Varscan, and Annovar). She also said something about using Sift and Polyphen to predict the effect of the mutation on the gene functionality (continuous score that is benign below a certain threshold, and destructive above it). Anyone familiar with those techniques?

[Michael.James.Clark]

Very cool! I was planning to put together a little Google Site going through how I analyze exome-seq that's very similar to this. Now I'm not sure I should bother!

[mirabilia]

I haven't tried it with prokaryiotic samples, but it should work actually (bwa, picard and samtools definitely work with prokaryotic data, not too sure about the GATK though...)

You need to adjust it though, for example index your own reference sequences and analysis depends on what sequence variation you'd expect (this pipeline works for diploid genomes only, though you might use some parts of it for different purposes)

Hope that helps.

thanks a lot ulz_peter!
Could you please clarify which steps of your pipeline are specifically for diploid genomes in order I can customize for my purposes?

[hanifk]

thanks very much

[ulz_peter]

I haven't tried it with prokaryiotic samples, but it should work actually (bwa, picard and samtools definitely work with prokaryotic data, not too sure about the GATK though...)

You need to adjust it though, for example index your own reference sequences and analysis depends on what sequence variation you'd expect (this pipeline works for diploid genomes only, though you might use some parts of it for different purposes)

Hope that helps.

[mirabilia]

hi folks,
thanks for sharing your expertise...it's a great help for a quite newbie like me.
I'm wondering if this analysis pipeline is suitable also for prokaryotic case or needs some adjustments. In case, could you suggest me some references?

thx!

[ulz_peter]

On the Unifier Genotyper I'm using the following parameters:

# # #Standard Raw VCF
java -Xmx15g -jar $GATK_DIR/GenomeAnalysisTK.jar -T UnifiedGenotyper \
-l INFO \
-I $OUT_DIR/exome.real.dedup.recal.bam \
-R $REFERENCE \
-B:intervals,BED $EXON_CAPTURE_FILE \
-B:dbsnp,VCF $DBSNP \
-glm BOTH \
-stand_call_conf 50.0 \
-stand_emit_conf 20.0 \
-dcov 300 \
-A AlleleBalance \
-A DepthOfCoverage \
-A FisherStrand \
-o $OUT_DIR/exome.raw.vcf \
-log $LOG_DIR/UnifiedGenotyper.log \
-nt 4

The company where this where done guarantees 30X of coverage ... (http://www.otogenetics.com/human_exome_page.htm)

I know this number should reduce after Variant Recalibrator ... I just want to know how many variants people are getting on this step.

By filtering out mutations you mean using the BED File to call only at the target regions ? If so, yes !

I use stand_emit_conf 10.0 and we usually get ~60k SNPs

[ulz_peter]

I just adapted the manual to fit it in the Wiki How to section:

Any changes, recommendations, complaints, etc. welcome:

http://seqanswers.com/wiki/How-to/exome_analysis