Name: Gap2Seq
Owner: Broad Institute
Description: Customization of Gap2Seq, a gap filling and insertion genotyping tool, for use with broadinstitute/viral-ngs; NOT the official Gap2Seq repo. .
Forked from: rikuu/Gap2Seq
Created: 2017-09-23 17:47:15.0
Updated: 2017-09-23 22:26:40.0
Pushed: 2017-09-25 02:38:59.0
Size: 2853
Language: C++
GitHub Committers
| User | Most Recent Commit | # Commits |
|---|
Other Committers
| User | Most Recent Commit | # Commits |
|---|
Gap2Seq is a program for filling gaps in scaffolds produced by genome assembly tools using short read data such as reads produced by Illumina sequencing.
Since version 3.0, it can also genotype insertions from variant calls produced by variant calling tools.
L. Salmela, K. Sahlin, V. Mäkinen, and A.I. Tomescu: Gap filling as exact path length problem. In Proc. RECOMB 2015, LNBI 9029, Springer 2015, pp. 281-292.
L. Salmela, A.I. Tomescu: Safely filling gaps with partial solutions common to all solutions. In Proc. WABI 2016, LNBI 9838, Springer 2016, xiv, short abstract.
R. Walve, L. Salmela, V. Mäkinen: Variant Genotyping with Gap Filling. (Submitted)
Gap2Seq has been tested on systems running Linux on a X86_64 architecture. Gap2Seq uses GATB library (http://gatb-core.gforge.inria.fr/index.html) for the de Bruijn graph implementation and htslib (http://www.htslib.org) for reading alignments for read filtering. The libraries are included in the Gap2Seq package.
Compiling Gap2Seq requires gcc version 4.5 or newer and cmake.
Unpack the Gap2Seq package. For compiling Gap2Seq run
r build; cd build; cmake ..; make
The main script Gap2Seq and the binaries Gap2Seq-core, GapMerger, GapCutter and ReadFilter can then be found in the build directory.
Seq [parameters]
ired parameters:
--filled <FASTA file> output file for filled scaffolds
--reads <FASTA/Q files> short reads, several files can be specified as a list separated by ','
--scaffolds <FASTA/Q file>
--gaps <FASTA/Q file>
--bed <BED file>
onal parameters:
--libraries <lib conf>
--threads <int> number of threads to use [default 1]
int> kmer length for DBG [default 31]
x-mem <float> maximum memory usage of DP table computation in gigabytes (excluding DBG) [default 20]
z <int> number of nucleotides to ignore on gap fringes [default 10]
st-error <int> maximum error in gap estimates [default 500]
lid <int> threshold for solid k-mers for building the DBG [default 2]
-upper If specified, all filled bases are in upper case.
que If specified, only gaps with a unique path of best length are filled.
t-only If specified, only paths that have optimal length are considered.
--help show this help message and exit
This example shows how to run Gap2Seq on the GAGE S. aureus data.
Download the GAGE data sets from
http://gage.cbcb.umd.edu/data/Staphylococcus_aureus/Data.original.tgz http://gage.cbcb.umd.edu/data/Staphylococcus_aureus/Assembly.tgz
Unpack the data files.
Run Gap2Seq (here we run it for the SGA scaffolds)
Seq --scaffolds Assembly/SGA/genome.scf.fasta --filled Assembly/SGA/genome.scf.fill.fasta --reads Data/original/frag_1.fastq,Data/original/frag_2.fastq,Data/original/shortjump_1.fastq,Data/original/shortjump_2.fastq
The filled scaffolds are then in the file Assembly/SGA/genome.scf.fill.fasta.
Align, sort, and index the read libraries to the scaffolds with e.g. BWA MEM.
index Assembly/SGA/genome.scf.fasta
mem Assembly/SGA/genome.scf.fasta Data/original/frag_1.fastq Data/original/frag_2.fastq | samtools sort -O bam - > Data/original/frag.bam
ools index Data/original/frag.bam
mem Assembly/SGA/genome.scf.fasta Data/original/shortjump_1.fastq Data/original/shortjump_2.fastq | samtools sort -O bam - > Data/original/shortjump.bam
ools index Data/original/shortjump.bam
Create a read library configuration file, a tab-separated list with a single read library per line: bam, read length, mean insert size, std dev, threshold
-e "Data/original/frag.bam\t101\t180\t18\t40" > libraries.txt
-e "Data/original/shortjump.bam\t37\t3500\t350\t40" >> libraries.txt
Run Gap2Seq.
Seq --scaffolds Assembly/SGA/genome.scf.fasta --filled Assembly/SGA/genome.scf.fill.fasta --libraries libraries.txt
When classifying filled bases into safe and unsafe bases, all paths within the allowed interval are now considered. In the previous version, only paths with optimal length were considered. The old behaviour can still be invoked using the parameter -best-only.
Gap2Seq.sh is replaced with a Python script, which accepts gaps/scaffolds in FASTA/FASTQ and VCF formats and reads in FASTA/FASTQ and SAM/BAM formats.
Optional per-gap read filtering when run with new script.
Gap2Seq binary is renamed Gap2Seq-core and can still be used instead of the new script.
Flanks of length between k and k+fuz are now used rather than the hard limit of k+fuz.
Switched to GATB-core 1.2.2.
Parallelization is now on gap level when run with the Gap2Seq.sh script.
Safe bases inserted into gaps are outputted in upper case and unsafe bases are outputted in lower case.
Optimized version of the algorithm.
Output now indicates on which gaps search was aborted because of the memory limit.
Reorganized parallel execution.
Proper synchronization for access to the memuse hash table.
Switched to GATB 1.0.5.
The maximum memory limitation option is now total for all threads. This is then divided evenly to all threads.
Memory usage tracking now includes all major data structures excluding the DBG.