Tagged with #samtools
1 documentation article | 0 announcements | 4 forum discussions



Created 2013-07-02 00:16:14 | Updated 2015-04-23 21:52:27 | Tags: install rscript igv picard gsalib samtools r ggplot2 rstudio
Comments (43)

Objective

Install all software packages required to follow the GATK Best Practices.

Prerequisites

To follow these instructions, you will need to have a basic understanding of the meaning of the following words and command-line operations. If you are unfamiliar with any of the following, you should consult a more experienced colleague or your systems administrator if you have one. There are also many good online tutorials you can use to learn the necessary notions.

  • Basic Unix environment commands
  • Binary / Executable
  • Compiling a binary
  • Adding a binary to your path
  • Command-line shell, terminal or console
  • Software library

You will also need to have access to an ANSI compliant C++ compiler and the tools needed for normal compilations (make, shell, the standard library, tar, gunzip). These tools are usually pre-installed on Linux/Unix systems. On MacOS X, you may need to install the MacOS Xcode tools. See https://developer.apple.com/xcode/ for relevant information and software downloads. The XCode tools are free but an AppleID may be required to download them.

Starting with version 2.6, the GATK requires Java Runtime Environment version 1.7. All Linux/Unix and MacOS X systems should have a JRE pre-installed, but the version may vary. To test your Java version, run the following command in the shell:

java -version 

This should return a message along the lines of ”java version 1.7.0_25” as well as some details on the Runtime Environment (JRE) and Virtual Machine (VM). If you have a version other than 1.7.x, be aware that you may run into trouble with some of the more advanced features of the Picard and GATK tools. The simplest solution is to install an additional JRE and specify which you want to use at the command-line. To find out how to do so, you should seek help from your systems administrator.

Software packages

  1. BWA
  2. SAMtools
  3. Picard
  4. Genome Analysis Toolkit (GATK)
  5. IGV
  6. RStudio IDE and R libraries ggplot2 and gsalib

Note that the version numbers of packages you download may be different than shown in the instructions below. If so, please adapt the number accordingly in the commands.


1. BWA

Read the overview of the BWA software on the BWA project homepage, then download the latest version of the software package.

  • Installation

Unpack the tar file using:

tar xvzf bwa-0.7.12.tar.bz2 

This will produce a directory called bwa-0.7.12 containing the files necessary to compile the BWA binary. Move to this directory and compile using:

cd bwa-0.7.12
make

The compiled binary is called bwa. You should find it within the same folder (bwa-0.7.12 in this example). You may also find other compiled binaries; at time of writing, a second binary called bwamem-lite is also included. You can disregard this file for now. Finally, just add the BWA binary to your path to make it available on the command line. This completes the installation process.

  • Testing

Open a shell and run:

bwa 

This should print out some version and author information as well as a list of commands. As the Usage line states, to use BWA you will always build your command lines like this:

bwa <command> [options] 

This means you first make the call to the binary (bwa), then you specify which command (method) you wish to use (e.g. index) then any options (i.e. arguments such as input files or parameters) used by the program to perform that command.


2. SAMtools

Read the overview of the SAMtools software on the SAMtools project homepage, then download the latest version of the software package.

  • Installation

Unpack the tar file using:

tar xvzf samtools-0.1.2.tar.bz2 

This will produce a directory called samtools-0.1.2 containing the files necessary to compile the SAMtools binary. Move to this directory and compile using:

cd samtools-0.1.2 
make 

The compiled binary is called samtools. You should find it within the same folder (samtools-0.1.2 in this example). Finally, add the SAMtools binary to your path to make it available on the command line. This completes the installation process.

  • Testing

Open a shell and run:

samtools 

This should print out some version information as well as a list of commands. As the Usage line states, to use SAMtools you will always build your command lines like this:

samtools <command> [options] 

This means you first make the call to the binary (samtools), then you specify which command (method) you wish to use (e.g. index) then any options (i.e. arguments such as input files or parameters) used by the program to perform that command. This is a similar convention as used by BWA.


3. Picard

Read the overview of the Picard software on the Picard project homepage, then download the latest version of the software package.

  • Installation

Unpack the zip file using:

tar xjf picard-tools-1.130.zip 

This will produce a directory called picard-tools-1.130 containing the Picard jar files. Picard tools are distributed as a pre-compiled Java executable (jar file) so there is no need to compile them.

Note that it is not possible to add jar files to your path to make the tools available on the command line; you have to specify the full path to the jar file in your java command, which would look like this:

java -jar ~/my_tools/jars/picard.jar <Toolname> [options]

This syntax will be explained in a little more detail further below.

However, you can set up a shortcut called an "environment variable" in your shell profile configuration to make this easier. The idea is that you create a variable that tells your system where to find a given jar, like this:

PICARD = "~/my_tools/jars/picard.jar"

So then when you want to run a Picard tool, you just need to call the jar by its shortcut, like this:

java -jar $PICARD <Toolname> [options]

The exact way to set this up depends on what shell you're using and how your environment is configured. We like this overview and tutorial which explains how it all works; but if you are new to the command line environment and you find this too much too deal with, we recommend asking for help from your institution's IT support group.

This completes the installation process.

  • Testing

Open a shell and run:

java -jar picard.jar -h 

This should print out some version and usage information about the AddOrReplaceReadGroups.jar tool. At this point you will have noticed an important difference between BWA and Picard tools. To use BWA, we called on the BWA program and specified which of its internal tools we wanted to apply. To use Picard, we called on Java itself as the main program, then specified which jar file to use, knowing that one jar file = one tool. This applies to all Picard tools; to use them you will always build your command lines like this:

java -jar picard.jar <ToolName> [options] 

This means you first make the call to Java itself as the main program, then specify the picard.jar file, then specify which tool you want, and finally you pass whatever other arguments (input files, parameters etc.) are needed for the analysis.

Note that the command-line syntax of Picard tools has recently changed from java -jar <ToolName>.jar to java -jar picard.jar <ToolName>. We are using the newer syntax in this document, but some of our other documents may not have been updated yet. If you encounter any documents using the old syntax, let us know and we'll update them accordingly. If you are already using an older version of Picard, either adapt the commands or better, upgrade your version!

Next we will see that GATK tools are called in essentially the same way, although the way the options are specified is a little different. The reasons for how tools in a given software package are organized and invoked are largely due to the preferences of the software developers. They generally do not reflect strict technical requirements, although they can have an effect on speed and efficiency.


4. Genome Analysis Toolkit (GATK)

Hopefully if you're reading this, you're already acquainted with the purpose of the GATK, so go ahead and download the latest version of the software package.

In order to access the downloads, you need to register for a free account on the GATK support forum. You will also need to read and accept the license agreement before downloading the GATK software package. Note that if you intend to use the GATK for commercial purposes, you will need to purchase a license. See the licensing page for an overview of the commercial licensing conditions.

  • Installation

Unpack the tar file using:

tar xjf GenomeAnalysisTK-3.3-0.tar.bz2 

This will produce a directory called GenomeAnalysisTK-3.3-0 containing the GATK jar file, which is called GenomeAnalysisTK.jar, as well as a directory of example files called resources. GATK tools are distributed as a single pre-compiled Java executable so there is no need to compile them. Just like we discussed for Picard, it's not possible to add the GATK to your path, but you can set up a shortcut to the jar file using environment variables as described above.

This completes the installation process.

  • Testing

Open a shell and run:

java -jar GenomeAnalysisTK.jar -h 

This should print out some version and usage information, as well as a list of the tools included in the GATK. As the Usage line states, to use GATK you will always build your command lines like this:

java -jar GenomeAnalysisTK.jar -T <ToolName> [arguments] 

This means that just like for Picard, you first make the call to Java itself as the main program, then specify the GenomeAnalysisTK.jar file, then specify which tool you want, and finally you pass whatever other arguments (input files, parameters etc.) are needed for the analysis.


5. IGV

The Integrated Genomics Viewer is a genome browser that allows you to view BAM, VCF and other genomic file information in context. It has a graphical user interface that is very easy to use, and can be downloaded for free (though registration is required) from this website. We encourage you to read through IGV's very helpful user guide, which includes many detailed tutorials that will help you use the program most effectively.


6. RStudio IDE and R libraries ggplot2 and gsalib

Download the latest version of RStudio IDE. The webpage should automatically detect what platform you are running on and recommend the version most suitable for your system.

  • Installation

Follow the installation instructions provided. Binaries are provided for all major platforms; typically they just need to be placed in your Applications (or Programs) directory. Open RStudio and type the following command in the console window:

install.packages("ggplot2") 

This will download and install the ggplot2 library as well as any other library packages that ggplot2 depends on for its operation. Note that some users have reported having to install two additional package themselves, called reshape and gplots, which you can do as follows:

install.packages("reshape")
install.packages("gplots")

Finally, do the same thing to install the gsalib library:

install.packages("gsalib")

This will download and install the gsalib library.

Important note

If you are using a recent version of ggplot2 and a version of GATK older than 3.2, you may encounter an error when trying to generate the BQSR or VQSR recalibration plots. This is because until recently our scripts were still using an older version of certain ggplot2 functions. This has been fixed in GATK 3.2, so you should either upgrade your version of GATK (recommended) or downgrade your version of ggplot2. If you experience further issues generating the BQSR recalibration plots, please see this tutorial.

No posts found with the requested search criteria.

Created 2015-08-14 19:03:41 | Updated | Tags: bam samtools
Comments (5)

Hi,

I’ve been trying to track down an issue that cropped up when we were validating our pipeline on a newer system. We have a test that produces different output each time it’s run (it seems to cycle randomly between five different outputs), but only on ubuntu 14. The same test produces the same result every time on ubuntu 12, and when run on a Mac OS X desktop.

I was able to create an isolated test using only SamReader and SamLocusIterator that simply iterates over the BAM and writes out every locus to a text file. This file exhibits the same behavior, e.g. it cycles between five different outputs.

The test code is here:

http://paste.openstack.org/show/412884/

Note that this will generate a lot of output, so it’s best to run with a small BAM file. My test file of 170M produces an 800M output file.

I tested using JDK 8u51 on all machines, and am using the latest picard tools (1.138) although it originally showed up using a much earlier version of picard tools (1.93) and JDK 7.

I’m not a bioinformatician so I don’t really know what the expected behavior is, although I couldn’t find the source of this variation by looking at the SAM tools source. If this is expected behavior that is OK, but I’d like confirmation of it before considering that the validation is complete. Better yet would be a way to disable or control this behavior to allow for consistent test results.

Thanks,

Phil Shapiro


Created 2015-04-14 07:56:45 | Updated 2015-04-14 07:57:44 | Tags: samtools combinegvcfs bcftools
Comments (2)

Hi GATK team, I've tried to use the gvcf produced by samtools/bcftools. Here is a snapshot of my Makefile for the calling part:

define gvcf_samtools_haloplex


$(1): $(2) $(3)
        mkdir -p $$(dir $$@)   &&       $${samtools.exe}  mpileup  --uncompressed --BCF --output-tags DP,DV,DP4,SP  --redo-BAQ --adjust-MQ 50 --min-ireads 3 --gap-frac 0.00
2  --min-MQ 30 --fasta-ref $$(REF)              --max-depth 10000               --max-idepth 10000              $$(if $$(filter %.bed,$$^),--positions $$(filter %.bed,$$^),
)               --count-orphans                 $$<  |  $${bcftools.exe} call  --gvcf 20  --output-type v --variants-only --multiallelic-caller --format-fields GQ,GP   -o $
$(addsuffix .tmp.vcf,$$@) -  &&         mv $$(addsuffix .tmp.vcf,$$@)  $$@


endef

when the gvcf files are used by CombineGVCFs, I get the following exception:

java -jarGenomeAnalysisTK.jar -T CombineGVCFs             -R /commun/data/pubdb/broadinstitute.org/bundle/1.5/b37/human_g1k_v37.fasta --variant OUTDIR/Projects/bigone/VCF/__DELETE__bigone2015.bigone.samtools.gvcf.list              -S SILENT -o OUTDIR/Projects/bigone/VCF/__DELETE__bigone2015.bigone.samtools.combined.g.vcf.tmp.vcf
dynamic mode

INFO  09:39:55,978 HelpFormatter - --------------------------------------------------------------------------------
INFO  09:39:55,981 HelpFormatter - The Genome Analysis Toolkit (GATK) v3.3-0-g37228af, Compiled 2014/10/24 01:07:22
INFO  09:39:55,981 HelpFormatter - Copyright (c) 2010 The Broad Institute
INFO  09:39:55,982 HelpFormatter - For support and documentation go to http://www.broadinstitute.org/gatk
INFO  09:39:55,989 HelpFormatter - Program Args: -T CombineGVCFs -R /commun/data/pubdb/broadinstitute.org/bundle/1.5/b37/human_g1k_v37.fasta --variant OUTDIR/Projects/bigone/VCF/__DELETE__bigone2015.bigone.samtools.gvcf.list -S SILENT -o OUTDIR/Projects/bigone/VCF/__DELETE__bigone2015.bigone.samtools.combined.g.vcf.tmp.vcf
INFO  09:39:55,994 HelpFormatter - Executing as lindenb@node005 on Linux 2.6.32-431.17.1.el6.x86_64 amd64; Java HotSpot(TM) 64-Bit Server VM 1.7.0_60-b19.
INFO  09:39:55,995 HelpFormatter - Date/Time: 2015/04/14 09:39:55
INFO  09:39:55,995 HelpFormatter - --------------------------------------------------------------------------------
INFO  09:39:55,996 HelpFormatter - --------------------------------------------------------------------------------
INFO  09:39:58,104 GenomeAnalysisEngine - Strictness is SILENT
INFO  09:39:58,269 GenomeAnalysisEngine - Downsampling Settings: Method: BY_SAMPLE, Target Coverage: 1000
INFO  09:40:18,824 GenomeAnalysisEngine - Preparing for traversal
INFO  09:40:18,859 GenomeAnalysisEngine - Done preparing for traversal
INFO  09:40:18,860 ProgressMeter - [INITIALIZATION COMPLETE; STARTING PROCESSING]
INFO  09:40:18,861 ProgressMeter -                 | processed |    time |    per 1M |           |   total | remaining
INFO  09:40:18,862 ProgressMeter -        Location |     sites | elapsed |     sites | completed | runtime |   runtime
##### ERROR ------------------------------------------------------------------------------------------
##### ERROR A USER ERROR has occurred (version 3.3-0-g37228af):
##### ERROR
##### ERROR This means that one or more arguments or inputs in your command are incorrect.
##### ERROR The error message below tells you what is the problem.
##### ERROR
##### ERROR If the problem is an invalid argument, please check the online documentation guide
##### ERROR (or rerun your command with --help) to view allowable command-line arguments for this tool.
##### ERROR
##### ERROR Visit our website and forum for extensive documentation and answers to
##### ERROR commonly asked questions http://www.broadinstitute.org/gatk
##### ERROR
##### ERROR Please do NOT post this error to the GATK forum unless you have really tried to fix it yourself.
##### ERROR
##### ERROR MESSAGE: Badly formed genome loc: Parameters to GenomeLocParser are incorrect:The stop position 32768065 is less than start 32768066 in contig 1
##### ERROR ------------------------------------------------------------------------------------------

here is a snapshot of the gvcf contaning 32768065:

 cat __DELETE__bigone2015.bigone.samtools.gvcf.list | xargs cut -f2,4,5,8 | grep -whF -A1 -B1 32768066


32757721        G       .       END=32757891
32768025        A       .       END=32768065
32768066        T       .       END=32768065
32768066        TC      T       INDEL;IDV=36;IMF=0.349515;DP=100;VDB=5.37817e-15;SGB=-0.693054;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=5,0,28,0;MQ=49
32768067        C       T       DP=61;VDB=1.14947e-07;SGB=-0.688148;RPB=0.0443614;MQB=0.986249;BQB=0.0581566;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=10,0,15,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=107;VDB=1.27981e-36;SGB=-0.693147;RPB=0.145801;MQB=0.659702;BQB=0.0567542;MQ0F=0;AC=2;AN=2;DP4=11,0,51,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=82;VDB=6.85164e-21;SGB=-0.693146;RPB=0.825581;MQB=0.697674;BQB=0.883721;MQ0F=0;AC=2;AN=2;DP4=2,0,43,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=84;VDB=7.44212e-14;SGB=-0.693145;RPB=0.030182;MQB=0.575018;BQB=0.0531329;MQ0F=0;AC=2;AN=2;DP4=6,0,41,0;MQ=48
--
32768025        A       .       END=32768064
32768065        T       C       DP=122;VDB=4.86965e-13;SGB=-0.692067;RPB=0.409269;MQB=0.707596;BQB=0.00055632;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=30,0,20,0;MQ=48
32768065        T       .       END=32768066
32768067        C       T       DP=106;VDB=1.4115e-10;SGB=-0.693054;RPB=0.0329847;MQB=0.883119;BQB=0.647624;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=19,0,28,0;MQ=49
--
32757679        C       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=157;VDB=2.99447e-25;SGB=-0.693147;RPB=0.446601;MQB=0.783611;BQB=0.323437;MQ0F=0;AC=2;AN=2;DP4=9,0,77,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=91;VDB=1.91056e-12;SGB=-0.693147;RPB=0.921275;MQB=0.857999;BQB=0.380388;MQ0F=0;AC=2;AN=2;DP4=6,0,49,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=35;VDB=1.38853e-08;SGB=-0.69168;RPB=0.524801;MQB=0.979651;BQB=0.524801;MQ0F=0;AC=2;AN=2;DP4=4,0,19,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=76;VDB=2.26381e-08;SGB=-0.693136;RPB=0.0423819;MQB=0.639856;BQB=0.0341585;MQ0F=0;AC=2;AN=2;DP4=7,0,35,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=63;VDB=2.20811e-09;SGB=-0.692976;RPB=0.170874;MQB=0.722875;BQB=0.0252673;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=5,0,26,0;MQ=49
--
32757684        C       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=40;VDB=6.04185e-10;SGB=-0.689466;RPB=0.25;MQB=0.9375;BQB=0.25;MQ0F=0;AC=2;AN=2;DP4=2,0,16,0;MQ=49
--
32757679        C       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=207;VDB=2.1795e-29;SGB=-0.693147;RPB=0.00231647;MQB=0.538424;BQB=0.00483428;MQ0F=0;AC=2;AN=2;DP4=13,0,99,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=43;VDB=1.23981e-07;SGB=-0.683931;MQ0F=0;AC=2;AN=2;DP4=0,0,13,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=28;VDB=3.51673e-06;SGB=-0.688148;RPB=0.333333;MQB=0.8;BQB=0.466667;MQ0F=0;AC=2;AN=2;DP4=2,0,15,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=38;VDB=3.77934e-07;SGB=-0.691153;RPB=0.173145;MQB=0.877714;BQB=0.32969;MQ0F=0;AC=2;AN=2;DP4=4,0,18,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=38;VDB=4.5048e-09;SGB=-0.692562;RPB=0.729381;MQB=0.881726;BQB=0.0429872;MQ0F=0;AC=2;AN=2;DP4=3,0,22,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768065
32768066        TC      T       INDEL;IDV=30;IMF=0.309278;DP=92;VDB=4.68847e-13;SGB=-0.693097;MQ0F=0;AC=2;AN=2;DP4=6,0,30,0;MQ=48
32768067        C       T       DP=61;VDB=0.00740154;SGB=-0.616816;RPB=0.231146;MQB=0.924584;BQB=0.924584;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=3,0,6,0;MQ=49
--
32757892        C       T       DP=2;VDB=0.02;SGB=-0.453602;MQ0F=0;AC=2;AN=2;DP4=0,0,0,2;MQ=50
32768025        A       .       END=32768066
32768067        C       T       DP=43;VDB=2.52435e-11;SGB=-0.692831;RPB=0.913904;MQB=0.953497;BQB=0.604728;MQ0F=0;AC=2;AN=2;DP4=3,0,24,0;MQ=48
--
32768025        A       .       END=32768050
32768066        TC      T       INDEL;IDV=4;IMF=0.153846;DP=26;VDB=1.68071e-06;SGB=-0.556411;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=7,0,4,0;MQ=48
32768067        C       T       DP=22;VDB=0.154358;SGB=-0.511536;MQ0F=0;AC=2;AN=2;DP4=0,0,3,0;MQ=48
--
32757688        T       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=106;VDB=5.33136e-17;SGB=-0.693147;MQ0F=0;AC=2;AN=2;DP4=0,0,50,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=52;VDB=2.64783e-06;SGB=-0.691153;RPB=0.78077;MQB=0.922371;BQB=0.0374746;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=3,0,18,0;MQ=48
--
32768025        A       .       END=32768050
32768051        T       .       END=32768066
32768067        C       T       DP=27;VDB=1.18758e-07;SGB=-0.690438;RPB=1;MQB=1;BQB=1;MQ0F=0;AC=2;AN=2;DP4=1,0,17,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768065
32768067        C       T       DP=35;VDB=0.00415585;SGB=-0.683931;RPB=0.0851108;MQB=0.832553;BQB=0.0731397;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=3,0,13,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=96;VDB=4.21707e-19;SGB=-0.693147;RPB=0.429291;MQB=0.986088;BQB=0.00988274;MQ0F=0;AC=2;AN=2;DP4=7,0,47,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768068        T       .       END=32768071
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=35;VDB=5.3147e-08;SGB=-0.69168;RPB=0.605263;MQB=0.736842;BQB=0.157895;MQ0F=0;AC=2;AN=2;DP4=2,0,19,0;MQ=49
--
32757645        A       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=111;VDB=1.28269e-05;SGB=-0.69168;RPB=0.0468862;MQB=0.524801;BQB=0.573544;MQ0F=0;AC=2;AN=2;DP4=4,0,19,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=29;VDB=0.000228416;SGB=-0.676189;RPB=0.646383;MQB=0.782349;BQB=0.092944;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=3,0,11,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=61;VDB=2.08113e-07;SGB=-0.693054;RPB=0.333887;MQB=0.808136;BQB=0.724199;MQ0F=0;AC=2;AN=2;DP4=5,0,28,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=26;VDB=1.62429e-07;SGB=-0.676189;RPB=1;MQB=1;BQB=1;MQ0F=0;AC=2;AN=2;DP4=1,0,11,0;MQ=50
--
32768025        A       .       END=32768050
32768066        TC      T       INDEL;IDV=4;IMF=0.173913;DP=22;VDB=0.0270423;SGB=-0.556411;MQ0F=0;AC=2;AN=2;DP4=1,0,4,0;MQ=48
32768068        T       .       END=32768070
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=57;VDB=9.28657e-08;SGB=-0.691153;RPB=0.0307499;MQB=0.969233;BQB=0.0900137;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=8,0,18,0;MQ=49
--
32757721        G       .       END=32757891
32768025        A       .       END=32768065
32768066        TC      T       INDEL;IDV=10;IMF=0.204082;DP=49;VDB=1.28199e-08;SGB=-0.651104;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=14,0,8,0;MQ=47
32768067        C       T       DP=39;VDB=0.0505719;SGB=-0.636426;RPB=1;MQB=1;BQB=1;MQ0F=0;AC=2;AN=2;DP4=1,0,7,0;MQ=47
--
32757645        A       .       END=32757891
32768025        A       .       END=32768065
32768066        TC      T       INDEL;IDV=23;IMF=0.164286;DP=136;VDB=5.53751e-13;SGB=-0.692067;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=13,0,20,0;MQ=48
32768067        C       T       DP=113;VDB=0.00011222;SGB=-0.686358;RPB=0.0294684;MQB=0.816272;BQB=0.38558;MQ0F=0;AC=2;AN=2;DP4=4,0,14,0;MQ=49
--
32768025        A       .       END=32768062
32768064        T       .       END=32768066
32768067        C       T       DP=31;VDB=2.8024e-07;SGB=-0.680642;RPB=0.861658;MQB=0.767432;BQB=0.913864;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=4,0,12,0;MQ=49
--
32757646        G       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=186;VDB=2.94273e-44;SGB=-0.693147;RPB=0.000575688;MQB=0.642812;BQB=0.0982708;MQ0F=0;AC=2;AN=2;DP4=14,0,106,0;MQ=47
--
32757722        G       .       END=32757891
32768026        A       .       END=32768065
32768066        T       .       END=32768065
32768066        TC      T       INDEL;IDV=42;IMF=0.488372;DP=83;VDB=5.40575e-21;SGB=-0.693145;MQ0F=0;AC=2;AN=2;DP4=7,0,40,0;MQ=49
32768067        C       T       DP=41;VDB=0.000382795;SGB=-0.651104;RPB=0.0172892;MQB=0.562066;BQB=0.959099;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=11,0,8,0;MQ=48
--
32768026        A       .       END=32768064
32768065        T       C       DP=46;VDB=4.85251e-07;SGB=-0.686358;RPB=0.670716;MQB=0.958545;BQB=0.000568099;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=16,0,14,0;MQ=49
32768066        T       .       END=32768066
32768068        T       .       END=32768073
--
32757722        G       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=34;VDB=1.38418e-05;SGB=-0.688148;RPB=0.7;MQB=0.933333;BQB=0.933333;MQ0F=0;AC=2;AN=2;DP4=2,0,15,0;MQ=48
--
32757722        G       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=75;VDB=2.16734e-11;SGB=-0.693136;RPB=0.00921987;MQB=0.811936;BQB=0.82152;MQ0F=0;AC=2;AN=2;DP4=6,0,35,0;MQ=48
--
32757646        G       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=46;VDB=1.10418e-07;SGB=-0.693079;RPB=0.181374;MQB=0.964868;BQB=0.981262;MQ0F=0;AC=2;AN=2;DP4=5,0,29,0;MQ=47
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=43;VDB=1.27519e-12;SGB=-0.692976;RPB=1;MQB=1;BQB=1;MQ0F=0;AC=2;AN=2;DP4=1,0,26,0;MQ=49
--
32757680        C       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=46;VDB=2.72057e-11;SGB=-0.692717;RPB=0.903753;MQB=0.354771;BQB=0.971625;MQ0F=0;AC=2;AN=2;DP4=5,0,23,0;MQ=49
--
32768026        A       .       END=32768060
32768066        T       .       END=32768066
32768067        C       T       DP=30;VDB=2.04238e-07;SGB=-0.690438;RPB=0.987474;MQB=0.892753;BQB=0.283511;MQ0F=0;AC=2;AN=2;DP4=3,0,17,0;MQ=48
--
32757646        G       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=62;VDB=2.44263e-09;SGB=-0.693054;MQ0F=0;AC=2;AN=2;DP4=0,0,28,0;MQ=48
--
32757646        G       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=134;VDB=7.06062e-26;SGB=-0.693147;RPB=0.000880621;MQB=0.588454;BQB=0.670928;MQ0F=0;AC=2;AN=2;DP4=10,0,75,0;MQ=48
--
32757680        C       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=52;VDB=9.35919e-11;SGB=-0.692831;RPB=0.114218;MQB=0.845852;BQB=0.651439;MQ0F=0;AC=2;AN=2;DP4=3,0,24,0;MQ=48
--
32768025        A       .       END=32768050
32768051        T       .       END=32768066
32768067        C       T       DP=29;VDB=8.52499e-06;SGB=-0.691153;RPB=1;MQB=1;BQB=1;MQ0F=0;AC=2;AN=2;DP4=1,0,18,0;MQ=48
--
32757646        G       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=220;VDB=0;SGB=-0.693147;RPB=0.167013;MQB=0.427009;BQB=0.162397;MQ0F=0;AC=2;AN=2;DP4=19,0,125,0;MQ=48
--
32757646        G       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=50;VDB=1.56275e-05;SGB=-0.690438;RPB=0.170712;MQB=0.860378;BQB=0.801124;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=5,0,17,0;MQ=48
--
32768025        A       .       END=32768050
32768051        T       .       END=32768066
32768067        C       T       DP=35;VDB=1.06013e-06;SGB=-0.686358;MQ0F=0;AC=2;AN=2;DP4=0,0,14,0;MQ=49
--
32757888        G       .       END=32757891
32768065        T       C       DP=18;VDB=0.0221621;SGB=-0.511536;RPB=1;MQB=1;BQB=1;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=1,0,3,0;MQ=50
32768066        TC      T       INDEL;IDV=5;IMF=0.277778;DP=18;VDB=0.00021572;SGB=-0.616816;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=2,0,6,0;MQ=49
32768188        T       .       END=32768247
--
32757680        C       .       END=32757891
32768026        A       .       END=32768066
32768067        C       T       DP=93;VDB=1.22543e-26;SGB=-0.693147;RPB=0.405658;MQB=0.999512;BQB=0.69141;MQ0F=0;AC=2;AN=2;DP4=4,0,53,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=51;VDB=1.59314e-10;SGB=-0.692067;RPB=1;MQB=1;BQB=1;MQ0F=0;AC=2;AN=2;DP4=1,0,20,0;MQ=48
--
32757888        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=51;VDB=0.00037052;SGB=-0.683931;RPB=0.48045;MQB=0.832553;BQB=0.641825;MQ0F=0;AC=2;AN=2;DP4=3,0,13,0;MQ=48
--
32757679        C       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=72;VDB=0.000573512;SGB=-0.686358;RPB=1;MQB=1;BQB=1;MQ0F=0;AC=2;AN=2;DP4=1,0,14,0;MQ=48
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=48;VDB=0.0303686;SGB=-0.616816;RPB=0;MQB=0.666667;BQB=0.333333;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=2,0,6,0;MQ=48
--
32757645        A       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=63;VDB=9.36785e-08;SGB=-0.692352;RPB=0.0501991;MQB=0.673273;BQB=0.606114;MQ0F=0;AC=2;AN=2;DP4=4,0,21,0;MQ=48
--
32768025        A       .       END=32768050
32768051        T       .       END=32768065
32768066        T       .       END=32768065
32768066        TC      T       INDEL;IDV=9;IMF=0.166667;DP=51;VDB=0.000464403;SGB=-0.636426;MQ0F=0;AC=2;AN=2;DP4=1,0,7,0;MQ=48
32768067        C       T       DP=40;VDB=3.93392e-06;SGB=-0.689466;RPB=0.602752;MQB=0.602752;BQB=0.855345;MQ0F=0;AC=2;AN=2;DP4=3,0,16,0;MQ=47
--
32757721        G       .       END=32757891
32768025        A       .       END=32768066
32768067        C       T       DP=56;VDB=3.41847e-09;SGB=-0.693097;RPB=0.995312;MQB=0.968732;BQB=0.183337;MQ0F=0;ICB=1;HOB=0.5;AC=1;AN=2;DP4=7,0,30,0;MQ=48

is it a problem with GATK or samtools ?

thanks !

Pierre


Created 2013-08-13 19:56:44 | Updated 2013-08-14 14:33:51 | Tags: unifiedgenotyper haplotypecaller samtools
Comments (2)

Hello All,

I'm a new user of haplotypecaller. I am currently comparing it to unifiedgenotyper and samtools. When I do a count of SNPs for all DPs, it looks like the is a multiple of 3 issues. There seems to be a lot more SNPs with a DP of 3, 6, 9, ... compared to 1, 2, 4, 5, 7, 8 ....

DP #ofSNP
1 43
2 41
3 105
4 24
5 18
6 232
7 38
8 66
9 3749

I was wondering if this was a normal behavior and what caused it.

Yannick


Created 2013-08-01 03:12:03 | Updated | Tags: filter gatk samtools
Comments (1)

Hi all,

I am comparing the variant calls from samtools and GATK. For samtools, I have been using quality score cut-offs 100 for SNP and 1000 for INDEL (quite stringent) and as a result, many variants are excluded after filtering. In case of GATK, I have been using our default setting (99% sensitivity for SNPs and 95% sensitivity for INDEL) and included only the variants with FILTER field "PASS". I was wondering if there is any more stringent filters that I can apply and that could be equivalent to samtools QS thresholds since it does not look like this is a fair comparison. Any of your suggestions will be appreciated.