changeset 0:25cd59a002d9 draft

planemo upload for repository https://github.com/galaxyproject/tools-iuc/tree/master/tools/genetrack commit e96df94dba60050fa28aaf55b5bb095717a5f260
author iuc
date Tue, 22 Dec 2015 17:03:27 -0500
parents
children df7ac50ade5d
files genetrack.py genetrack.xml genetrack_macros.xml genetrack_util.py static/images/genetrack.png test-data/genetrack_input2.gff test-data/genetrack_input3.scidx test-data/genetrack_input_unsorted4.gff test-data/genetrack_output2.gff test-data/genetrack_output3.gff test-data/genetrack_output4.gff tool_dependencies.xml
diffstat 12 files changed, 1142 insertions(+), 0 deletions(-) [+]
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/genetrack.py	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,65 @@
+"""
+genetrack.py
+
+Input: either scidx or gff format of reads
+Output: Called peaks in gff format
+"""
+import optparse
+import csv
+import os
+import genetrack_util
+
+CHUNK_SIZE = 10000000
+
+
+if __name__ == '__main__':
+    parser = optparse.OptionParser()
+    parser.add_option('-t', '--input_format', dest='input_format', type='string', help='Input format')
+    parser.add_option('-i', '--input', dest='inputs', type='string', action='append', nargs=2, help='Input datasets')
+    parser.add_option('-s', '--sigma', dest='sigma', type='int', default=5, help='Sigma.')
+    parser.add_option('-e', '--exclusion', dest='exclusion', type='int', default=20, help='Exclusion zone.')
+    parser.add_option('-u', '--up_width', dest='up_width', type='int', default=10, help='Upstream width of called peaks.')
+    parser.add_option('-d', '--down_width', dest='down_width', type='int', default=10, help='Downstream width of called peaks.')
+    parser.add_option('-f', '--filter', dest='filter', type='int', default=1, help='Absolute read filter.')
+    options, args = parser.parse_args()
+
+    os.mkdir('output')
+    for (dataset_path, hid) in options.inputs:
+        if options.input_format == 'gff':
+            # Make sure the reads for each chromosome are sorted by index.
+            input_path = genetrack_util.sort_chromosome_reads_by_index(dataset_path)
+        else:
+            # We're processing scidx data.
+            input_path = dataset_path
+        output_name = 's%se%su%sd%sF%s_on_data_%s' % (options.sigma,
+                                                      options.exclusion,
+                                                      options.up_width,
+                                                      options.down_width,
+                                                      options.filter,
+                                                      hid)
+        output_path = os.path.join('output', output_name)
+        reader = csv.reader(open(input_path, 'rU'), delimiter='\t')
+        writer = csv.writer(open(output_path, 'wt'), delimiter='\t')
+        width = options.sigma * 5
+        manager = genetrack_util.ChromosomeManager(reader)
+        while not manager.done:
+            cname = manager.chromosome_name()
+            # Should we process this chromosome?
+            data = manager.load_chromosome()
+            if not data:
+                continue
+            keys = genetrack_util.make_keys(data)
+            lo, hi = genetrack_util.get_range(data)
+            for chunk in genetrack_util.get_chunks(lo, hi, size=CHUNK_SIZE, overlap=width):
+                (slice_start, slice_end), process_bounds = chunk
+                window = genetrack_util.get_window(data, slice_start, slice_end, keys)
+                genetrack_util.process_chromosome(cname,
+                                                  window,
+                                                  writer,
+                                                  process_bounds,
+                                                  width,
+                                                  options.sigma,
+                                                  options.up_width,
+                                                  options.down_width,
+                                                  options.exclusion,
+                                                  options.filter)
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/genetrack.xml	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,173 @@
+<?xml version="1.0"?>
+<tool id="genetrack" name="GeneTrack" version="@WRAPPER_VERSION@.0">
+    <description>peak predictor</description>
+    <macros>
+        <import>genetrack_macros.xml</import>
+    </macros>
+    <expand macro="requirements" />
+    <command>
+        python $__tool_directory__/genetrack.py
+        --input_format $input_format_cond.input_format
+        #for $i in $input_format_cond.input:
+            --input "${i}" "${i.hid}"
+        #end for
+        --sigma $sigma
+        --exclusion $exclusion
+        --up_width $up_width
+        --down_width $down_width
+        --filter $filter
+    </command>
+    <inputs>
+        <conditional name="input_format_cond">
+            <param name="input_format" type="select" label="Format of files for conversion">
+                <option value="scidx" selected="True">ScIdx</option>
+                <option value="gff">Gff</option>
+            </param>
+            <when value="scidx">
+                <param name="input" type="data" format="scidx" multiple="True" label="Predict peaks on" />
+            </when>
+            <when value="gff">
+                <param name="input" type="data" format="gff" multiple="True" label="Predict peaks on" />
+            </when>
+        </conditional>
+        <param name="sigma" type="integer" value="5" min="1" label="Sigma to use when smoothing reads" help="Higher values increase computation but produce more smoothing." />
+        <param name="exclusion" type="integer" value="20" min="1" label="Peak exclusion zone" help="Exclusion zone around each peak that prevents others from being called." />
+        <param name="up_width" type="integer" value="10" min="0" label="Exclusion zone of upstream called peaks" />
+        <param name="down_width" type="integer" value="10" min="0" label="Exclusion zone of downstream called peaks" />
+        <param name="filter" type="integer" value="1" min="0" label="Absolute read filter" help="Removes peaks with lower peak height." />
+    </inputs>
+    <outputs>
+        <collection name="genetrack_output" type="list" label="Genetrack results on ${on_string}">
+            <discover_datasets pattern="(?P&lt;designation&gt;.*)" directory="output" ext="gff" visible="false" />
+        </collection>
+    </outputs>
+    <tests>
+        <test>
+            <param name="input" value="genetrack_input2.gff" ftype="gff" />
+            <param name="input_format" value="gff" />
+            <param name="sigma" value="5" />
+            <param name="exclusion" value="20" />
+            <param name="up_width" value="10" />
+            <param name="down_width" value="10" />
+            <param name="filter" value="3" />
+            <output_collection name="genetrack_output" type="list">
+                <element name="s5e20u10d10F3_on_data_1" file="genetrack_output2.gff" ftype="gff" />
+            </output_collection>
+        </test>
+        <test>
+            <param name="input" value="genetrack_input3.scidx" ftype="scidx" />
+            <param name="input_format" value="scidx" />
+            <param name="sigma" value="5" />
+            <param name="exclusion" value="20" />
+            <param name="up_width" value="10" />
+            <param name="down_width" value="10" />
+            <param name="filter" value="3" />
+            <output_collection name="genetrack_output" type="list">
+                <element name="s5e20u10d10F3_on_data_1" file="genetrack_output3.gff" ftype="gff" />
+            </output_collection>
+        </test>
+        <test>
+            <param name="input" value="genetrack_input_unsorted4.gff" ftype="gff" />
+            <param name="input_format" value="gff" />
+            <param name="sigma" value="5" />
+            <param name="exclusion" value="20" />
+            <param name="up_width" value="10" />
+            <param name="down_width" value="10" />
+            <param name="filter" value="3" />
+            <output_collection name="genetrack_output" type="list">
+                <element name="s5e20u10d10F3_on_data_1" file="genetrack_output4.gff" ftype="gff" />
+            </output_collection>
+        </test>
+    </tests>
+    <help>
+**What it does**
+
+GeneTrack separately identifies peaks on the forward "+” (W) and reverse “-” (C) strand.  The way that GeneTrack
+works is to replace each tag with a probabilistic distribution of occurrences for that tag at and around its mapped
+genomic coordinate.  The distance decay of the probabilistic distribution is set by adjusting the value of the
+tool's **Sigma to use when smoothing reads** parameter.  GeneTrack then sums the distribution over all mapped
+tags.  This results in a smooth continuous trace that can be globally broadened or tightened by adjusting the
+sigma value.  GeneTrack starts with the highest smoothed peak first, treating each strand separately if indicated
+by the data, then sets up an exclusion zone (centered over the peak) defined by the value of the **Peak exclusion
+zone** parameter (see figure).  The exclusion zone prevents any secondary peaks from being called on the same strand
+within that exclusion zone.  In rare cases, it may be desirable to set different exclusion zones upstream (more 5’)
+versus downstream (more 3’) of the peak.
+
+.. image:: $PATH_TO_IMAGES/genetrack.png
+
+GeneTrack continues through the data in order of peak height, until no other peaks are found, and in principle will
+call a peak at a single isolated tag, if no filter is set using the tool's **Absolute read filter** parameter.  A 
+filter value of 1 means that it will stop calling peaks when the tag count in the peak hits 1 (so single tag peaks
+will be excluded in this case).  GeneTrack outputs **chrom** (chromosome number), **strand** (+/W or -/C strand),
+**start** (lower coordinate of exclusion zone), **end** (higher coordinate of exclusion zone), and **value** (peak
+height).  Genetrack's GFF output reports the start (lower coordinate) and end (higher coordinate) of the exclusion
+zone.
+
+In principle, the width of the exclusion zone may be as large as the DNA region occupied by the native protein plus
+a steric exclusion zone between the protein and the exonuclease.  On the other hand the site might be considerably
+smaller if the protein is in a denatured state during exonuclease digestion (since it is pre-treated with SDS).
+
+In general, higher resolution data or smaller binding site size data should use smaller sigma values. Large binding
+site size data such as 147 bp nucleosomal DNA use a larger sigma value like 20 (-s 20).  For transcription factors
+mapped by ChIP-exo, sigma may initially be set at 5, and the exclusion zone set at 20 (-s 5 –e 20).  Sigma is typically
+varied between ~3 and ~20. Too high of a sigma value may merge two independent nearby binding events.  This may be
+desirable if closely bound factors are not distinguishable.  Too low of a sigma value will cause some tags that
+contribute to a binding event to be excluded, because they may not be located sufficiently close to the main peak.
+If alternative (mutually exclusive) binding is expected for two overlapping sites, and these sites are to be
+independently recorded, then an empirically determined smaller exclusion zone width is set.  Thus the value of sigma
+is set empirically for each mapped factor, depending upon the resolution and binding site size of the binding event.
+
+It might make sense to exclude peaks that have only a single tag, where -F 1 is used, or have their tags located on
+only a single coordinate (called Singletons, where stddev=0 in the output file).  However, low coverage datasets might
+be improved by including them, if additional analysis (e.g., motif discovery) validates them. In addition, idealized
+action of the exonuclease in ChIP-exo might place all tags for a peak on a single coordinate.
+
+-----
+
+**Options**
+
+ * **Sigma to use when smoothing reads** - Smooths clusters of tags via a Gaussian distribution.
+ * **Peak exclusion zone** - Exclusion zone around each peak, eliminating all other peaks on the same strand that are within a ± bp distance of the peak.
+ * **Exclusion zone of upstream called peaks** - Defines the exclusion zone centered over peaks upstream of a peak.
+ * **Exclusion zone of downstream called peaks** - Defines the exclusion zone centered over peaks downstream of a peak.
+ * **Filter** - Absolute read filter, restricts output to only peaks with larger peak height.
+ 
+-----
+
+**Output gff Columns**
+
+1. Chromosome
+2. Script
+3. Placeholder (no meaning)
+4. Start of peak exclusion zone (-e 20)
+5. End of peak exclusion zone
+6. Tag sum (not peak height or area under curve, which LionDB provides)
+7. Strand
+8. Placeholder (no meaning)
+9. Attributes (standard deviation of reads located within exclusion zone) = fuzziness of peak
+
+-----
+ 
+**Considerations**
+ 
+In principle, the width of the exclusion zone may be as large as the DNA region occupied by the native protein
+plus a steric exclusion zone between the protein and the exonuclease.  On the other hand the site might be considerably
+smaller if the protein is in a denatured state during exonuclease digestion (since it is pre-treated with SDS).
+ 
+In general, higher resolution data or smaller binding site size data should use smaller sigma values.  Large binding site
+size data such as 147 bp nucleosomal DNA use a larger sigma value like 20 (-s 20).  For transcription factors mapped by
+ChIP-exo, sigma may initially be set at 5, and the exclusion zone set at 20 (-s 5 –e 20).  Sigma is typically varied
+between ~3 and ~20.  Too high of a sigma value may merge two independent nearby binding events.  This may be desirable if
+closely bound factors are not distinguishable.  Too low of a sigma value will cause some tags that contribute to a binding
+event to be excluded, because they may not be located sufficiently close to the main peak.  If alternative (mutually
+exclusive) binding is expected for two overlapping sites, and these sites are to be independently recorded, then an
+empirically determined smaller exclusion zone width is set.  Thus, the value of sigma is set empirically for each mappedfactor depending upon the resolution and binding site size of the binding event.
+
+It might make sense to exclude peaks that have only a single tag, where -F 1 is used, or have their tags located on only
+a single coordinate (called Singletons, where stddev=0 in the output file).  However, low coverage datasets might be
+improved by including them, if additional analysis (e.g., motif discovery) validates them.  In addition, idealized action
+of the exonuclease in ChIP-exo might place all tags for a peak on a single coordinate.
+
+    </help>
+    <expand macro="citations" />
+</tool>
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/genetrack_macros.xml	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,22 @@
+<?xml version='1.0' encoding='UTF-8'?>
+<macros>
+    <token name="@WRAPPER_VERSION@">1.0</token>
+    <xml name="requirements">
+        <requirements>
+            <requirement type="package" version="2.3.0">anaconda</requirement>
+        </requirements>
+    </xml>
+    <xml name="stdio">
+        <stdio>
+            <exit_code range="1:"/>
+            <exit_code range=":-1"/>
+            <regex match="Error:"/>
+            <regex match="Exception:"/>
+        </stdio>
+    </xml>
+    <xml name="citations">
+        <citations>
+            <citation type="doi">10.1093/bioinformatics/btn119</citation>
+        </citations>
+    </xml>
+</macros>
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/genetrack_util.py	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,390 @@
+import bisect
+import math
+import numpy
+import re
+import subprocess
+import sys
+import tempfile
+
+GFF_EXT = 'gff'
+SCIDX_EXT = 'scidx'
+
+
+def noop(data):
+    return data
+
+
+def zeropad_to_numeric(data):
+    return re.sub(r'chr0(\d)', r'chr\1', data)
+
+
+def numeric_to_zeropad(data):
+    return re.sub(r'chr(\d([^\d]|$))', r'chr0\1', data)
+
+
+FORMATS = ['zeropad', 'numeric']
+IN_CONVERT = {'zeropad': zeropad_to_numeric, 'numeric': noop}
+OUT_CONVERT = {'zeropad': numeric_to_zeropad, 'numeric': noop}
+
+
+def conversion_functions(in_fmt, out_fmt):
+    """
+    Returns the proper list of functions to apply to perform a conversion
+    """
+    return [IN_CONVERT[in_fmt], OUT_CONVERT[out_fmt]]
+
+
+def convert_data(data, in_fmt, out_fmt):
+    for fn in conversion_functions(in_fmt, out_fmt):
+        data = fn(data)
+    return data
+
+
+class ChromosomeManager(object):
+    """
+    Manages a CSV reader of an index file to only load one chrom at a time
+    """
+
+    def __init__(self, reader):
+        self.done = False
+        self.reader = reader
+        self.processed_chromosomes = []
+        self.current_index = 0
+        self.next_valid()
+
+    def next(self):
+        self.line = self.reader.next()
+
+    def is_valid(self, line):
+        if len(line) not in [4, 5, 9]:
+            return False
+        try:
+            [int(i) for i in line[1:]]
+            self.format = SCIDX_EXT
+            return True
+        except ValueError:
+            try:
+                if len(line) < 6:
+                    return False
+                [int(line[4]), int(line[5])]
+                self.format = GFF_EXT
+                return True
+            except ValueError:
+                return False
+
+    def next_valid(self):
+        """
+        Advance to the next valid line in the reader
+        """
+        self.line = self.reader.next()
+        s = 0
+        while not self.is_valid(self.line):
+            self.line = self.reader.next()
+            s += 1
+        if s > 0:
+            # Skip initial line(s) of file
+            pass
+
+    def parse_line(self, line):
+        if self.format == SCIDX_EXT:
+            return [int(line[1]), int(line[2]), int(line[3])]
+        else:
+            return [int(line[3]), line[6], line[5]]
+
+    def chromosome_name(self):
+        """
+        Return the name of the chromosome about to be loaded
+        """
+        return self.line[0]
+
+    def load_chromosome(self, collect_data=True):
+        """
+        Load the current chromosome into an array and return it
+        """
+        cname = self.chromosome_name()
+        if cname in self.processed_chromosomes:
+            stop_err('File is not grouped by chromosome')
+        self.data = []
+        while self.line[0] == cname:
+            if collect_data:
+                read = self.parse_line(self.line)
+                if read[0] < self.current_index:
+                    msg = 'Reads in chromosome %s are not sorted by index. (At index %d)' % (cname, self.current_index)
+                    stop_err(msg)
+                self.current_index = read[0]
+                self.add_read(read)
+            try:
+                self.next()
+            except StopIteration:
+                self.done = True
+                break
+        self.processed_chromosomes.append(cname)
+        self.current_index = 0
+        data = self.data
+        # Don't retain reference anymore to save memory
+        del self.data
+        return data
+
+    def add_read(self, read):
+        if self.format == SCIDX_EXT:
+            self.data.append(read)
+        else:
+            index, strand, value = read
+            if value == '' or value == '.':
+                value = 1
+            else:
+                value = int(value)
+            if not self.data:
+                self.data.append([index, 0, 0])
+                current_read = self.data[-1]
+            if self.data[-1][0] == index:
+                current_read = self.data[-1]
+            elif self.data[-1][0] < index:
+                self.data.append([index, 0, 0])
+                current_read = self.data[-1]
+            else:
+                msg = 'Reads in chromosome %s are not sorted by index. (At index %d)' % (self.chromosome_name(), index)
+                stop_err(msg)
+            if strand == '+':
+                current_read[1] += value
+            elif strand == '-':
+                current_read[2] += value
+            else:
+                msg = 'Strand "%s" at chromosome "%s" index %d is not valid.' % (strand, self.chromosome_name(), index)
+                stop_err(msg)
+
+    def skip_chromosome(self):
+        """
+        Skip the current chromosome, discarding data
+        """
+        self.load_chromosome(collect_data=False)
+
+
+class Peak(object):
+    def __init__(self, index, pos_width, neg_width):
+        self.index = index
+        self.start = index - neg_width
+        self.end = index + pos_width
+        self.value = 0
+        self.deleted = False
+        self.safe = False
+
+    def __repr__(self):
+        return '[%d] %d' % (self.index, self.value)
+
+
+def gff_row(cname, start, end, score, source, type='.', strand='.', phase='.', attrs={}):
+    return (cname, source, type, start, end, score, strand, phase, gff_attrs(attrs))
+
+
+def gff_attrs(d):
+    if not d:
+        return '.'
+    return ';'.join('%s=%s' % item for item in d.items())
+
+
+def stop_err(msg):
+    sys.stderr.write(msg)
+    sys.exit(1)
+
+
+def is_int(i):
+    try:
+        int(i)
+        return True
+    except ValueError:
+        return False
+
+
+def make_keys(data):
+    return [read[0] for read in data]
+
+
+def make_peak_keys(peaks):
+    return [peak.index for peak in peaks]
+
+
+def get_window(data, start, end, keys):
+    """
+    Returns all reads from the data set with index between the two indexes
+    """
+    start_index = bisect.bisect_left(keys, start)
+    end_index = bisect.bisect_right(keys, end)
+    return data[start_index:end_index]
+
+
+def get_index(value, keys):
+    """
+    Returns the index of the value in the keys using bisect
+    """
+    return bisect.bisect_left(keys, value)
+
+
+def get_range(data):
+    lo = min([item[0] for item in data])
+    hi = max([item[0] for item in data])
+    return lo, hi
+
+
+def get_chunks(lo, hi, size, overlap=500):
+    """
+    Divides a range into chunks of maximum size size. Returns a list of
+    2-tuples (slice_range, process_range), each a 2-tuple (start, end).
+    process_range has zero overlap and should be given to process_chromosome
+    as-is, and slice_range is overlapped and should be used to slice the
+    data (using get_window) to be given to process_chromosome.
+    """
+    chunks = []
+    for start_index in range(lo, hi, size):
+        process_start = start_index
+        # Don't go over upper bound
+        process_end = min(start_index + size, hi)
+        # Don't go under lower bound
+        slice_start = max(process_start - overlap, lo)
+        # Don't go over upper bound
+        slice_end = min(process_end + overlap, hi)
+        chunks.append(((slice_start, slice_end), (process_start, process_end)))
+    return chunks
+
+
+def allocate_array(data, width):
+    """
+    Allocates a new array with the dimensions required to fit all reads in
+    the argument. The new array is totally empty. Returns the array and the
+    shift (number to add to a read index to get the position in the array it
+    should be at).
+    """
+    lo, hi = get_range(data)
+    rng = hi - lo
+    shift = width - lo
+    return numpy.zeros(rng+width*2, numpy.float), shift
+
+
+def normal_array(width, sigma, normalize=True):
+    """
+    Returns an array of the normal distribution of the specified width
+    """
+    sigma2 = float(sigma)**2
+
+    def normal_func(x):
+        return math.exp(-x * x / (2 * sigma2))
+
+    # width is the half of the distribution
+    values = map(normal_func, range(-width, width))
+    values = numpy.array(values, numpy.float)
+    # normalization
+    if normalize:
+        values = 1.0/math.sqrt(2 * numpy.pi * sigma2) * values
+    return values
+
+
+def call_peaks(array, shift, data, keys, direction, down_width, up_width, exclusion):
+    peaks = []
+
+    def find_peaks():
+        # Go through the array and call each peak
+        results = (array > numpy.roll(array, 1)) & (array > numpy.roll(array, -1))
+        indexes = numpy.where(results)
+        for index in indexes[0]:
+            pos = down_width or exclusion // 2
+            neg = up_width or exclusion // 2
+            # Reverse strand
+            if direction == 2:
+                # Swap positive and negative widths
+                pos, neg = neg, pos
+            peaks.append(Peak(int(index)-shift, pos, neg))
+    find_peaks()
+
+    def calculate_reads():
+        # Calculate the number of reads in each peak
+        for peak in peaks:
+            reads = get_window(data, peak.start, peak.end, keys)
+            peak.value = sum([read[direction] for read in reads])
+            # Flat list of indexes with frequency
+            indexes = [r for read in reads for r in [read[0]] * read[direction]]
+            peak.stddev = numpy.std(indexes)
+    calculate_reads()
+
+    def perform_exclusion():
+        # Process the exclusion zone
+        peak_keys = make_peak_keys(peaks)
+        peaks_by_value = peaks[:]
+        peaks_by_value.sort(key=lambda peak: -peak.value)
+        for peak in peaks_by_value:
+            peak.safe = True
+            window = get_window(peaks,
+                                peak.index-exclusion//2,
+                                peak.index+exclusion//2,
+                                peak_keys)
+            for excluded in window:
+                if excluded.safe:
+                    continue
+                i = get_index(excluded.index, peak_keys)
+                del peak_keys[i]
+                del peaks[i]
+    perform_exclusion()
+    return peaks
+
+
+def process_chromosome(cname, data, writer, process_bounds, width, sigma, down_width, up_width, exclusion, filter):
+    """
+    Process a chromosome. Takes the chromosome name, list of reads, a CSV
+    writer to write processes results to, the bounds (2-tuple) to write
+    results in, and options.
+    """
+    if not data:
+        return
+    keys = make_keys(data)
+    # Create the arrays that hold the sum of the normals
+    forward_array, forward_shift = allocate_array(data, width)
+    reverse_array, reverse_shift = allocate_array(data, width)
+    normal = normal_array(width, sigma)
+
+    def populate_array():
+        # Add each read's normal to the array
+        for read in data:
+            index, forward, reverse = read
+            # Add the normals to the appropriate regions
+            if forward:
+                forward_array[index+forward_shift-width:index+forward_shift+width] += normal * forward
+            if reverse:
+                reverse_array[index+reverse_shift-width:index+reverse_shift+width] += normal * reverse
+    populate_array()
+    forward_peaks = call_peaks(forward_array, forward_shift, data, keys, 1, down_width, up_width, exclusion)
+    reverse_peaks = call_peaks(reverse_array, reverse_shift, data, keys, 2, down_width, up_width, exclusion)
+    # Convert chromosome name in preparation for writing output
+    cname = convert_data(cname, 'zeropad', 'numeric')
+
+    def write(cname, strand, peak):
+        start = max(peak.start, 1)
+        end = peak.end
+        value = peak.value
+        stddev = peak.stddev
+        if value > filter:
+            # This version of genetrack outputs only gff files.
+            writer.writerow(gff_row(cname=cname,
+                                    source='genetrack',
+                                    start=start,
+                                    end=end,
+                                    score=value,
+                                    strand=strand,
+                                    attrs={'stddev': stddev}))
+
+    for peak in forward_peaks:
+        if process_bounds[0] < peak.index < process_bounds[1]:
+            write(cname, '+', peak)
+    for peak in reverse_peaks:
+        if process_bounds[0] < peak.index < process_bounds[1]:
+            write(cname, '-', peak)
+
+
+def sort_chromosome_reads_by_index(input_path):
+    """
+    Return a gff file with chromosome reads sorted by index.
+    """
+    # Will this sort produce different results across platforms?
+    output_path = tempfile.NamedTemporaryFile(delete=False).name
+    command = 'sort -k 1,1 -k 4,4n "%s" > "%s"' % (input_path, output_path)
+    p = subprocess.Popen(command, shell=True)
+    p.wait()
+    return output_path
Binary file static/images/genetrack.png has changed
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/genetrack_input2.gff	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,100 @@
+chr1	genetrack	.	17	37	918	+	.	stddev=5.96715849116
+chr1	genetrack	.	31	51	245	-	.	stddev=2.66582799529
+chr1	genetrack	.	40	60	2060	+	.	stddev=2.7859667372
+chr1	genetrack	.	62	82	1300	-	.	stddev=4.13061337623
+chr1	genetrack	.	73	93	397	+	.	stddev=0.0
+chr1	genetrack	.	89	109	521	+	.	stddev=0.747112137937
+chr1	genetrack	.	123	143	5129	+	.	stddev=3.01025384354
+chr1	genetrack	.	125	145	4659	-	.	stddev=3.8642622228
+chr1	genetrack	.	155	175	897	-	.	stddev=3.22709952671
+chr1	genetrack	.	171	191	956	-	.	stddev=4.95899971687
+chr1	genetrack	.	180	200	1527	+	.	stddev=4.62574275346
+chr1	genetrack	.	185	205	494	-	.	stddev=1.4255957
+chr1	genetrack	.	192	212	2538	+	.	stddev=5.04731591122
+chr1	genetrack	.	206	226	2087	-	.	stddev=3.6160253713
+chr1	genetrack	.	238	258	2496	+	.	stddev=2.11105291581
+chr1	genetrack	.	242	262	5047	-	.	stddev=3.62629343395
+chr1	genetrack	.	254	274	1525	+	.	stddev=4.46082441647
+chr1	genetrack	.	281	301	15	+	.	stddev=1.74610678049
+chr1	genetrack	.	302	322	626	-	.	stddev=0.0
+chr1	genetrack	.	308	328	1544	+	.	stddev=4.43066151722
+chr1	genetrack	.	334	354	533	+	.	stddev=1.34355443899
+chr1	genetrack	.	344	364	726	-	.	stddev=1.36767079956
+chr1	genetrack	.	347	367	286	+	.	stddev=0.0
+chr1	genetrack	.	358	378	792	-	.	stddev=1.47737416556
+chr1	genetrack	.	374	394	608	+	.	stddev=1.44652711793
+chr1	genetrack	.	389	409	126	-	.	stddev=0.471404520791
+chr1	genetrack	.	439	459	618	-	.	stddev=5.47536569145
+chr1	genetrack	.	441	461	1393	+	.	stddev=4.75587332865
+chr1	genetrack	.	461	481	754	-	.	stddev=3.28891288785
+chr1	genetrack	.	483	503	58	+	.	stddev=0.0
+chr1	genetrack	.	538	558	1015	-	.	stddev=0.0
+chr1	genetrack	.	728	748	39	-	.	stddev=0.0
+chr1	genetrack	.	757	777	23	+	.	stddev=0.0
+chr1	genetrack	.	799	819	607	+	.	stddev=0.0
+chr1	genetrack	.	844	864	665	+	.	stddev=0.0
+chr1	genetrack	.	877	897	468	+	.	stddev=0.0
+chr1	genetrack	.	903	923	107	-	.	stddev=0.0
+chr1	genetrack	.	944	964	2	-	.	stddev=0.0
+chr1	genetrack	.	1092	1112	740	+	.	stddev=0.0
+chr1	genetrack	.	1127	1147	940	-	.	stddev=3.96036497305
+chr1	genetrack	.	1183	1203	25	+	.	stddev=0.0
+chr1	genetrack	.	1291	1311	454	-	.	stddev=0.0
+chr1	genetrack	.	1329	1349	207	-	.	stddev=0.0
+chr1	genetrack	.	1484	1504	584	+	.	stddev=0.0
+chr1	genetrack	.	2075	2095	1181	+	.	stddev=0.0
+chr1	genetrack	.	2102	2122	481	+	.	stddev=0.0455486534308
+chr1	genetrack	.	2125	2145	199	-	.	stddev=0.0
+chr1	genetrack	.	2452	2472	1246	+	.	stddev=0.0
+chr1	genetrack	.	2602	2622	34	+	.	stddev=0.0
+chr1	genetrack	.	2833	2853	1062	+	.	stddev=1.01561431542
+chr1	genetrack	.	2838	2858	1144	-	.	stddev=1.09438744148
+chr1	genetrack	.	3011	3031	1212	-	.	stddev=0.0
+chr1	genetrack	.	3116	3136	555	-	.	stddev=0.0
+chr1	genetrack	.	3130	3150	17	+	.	stddev=0.0
+chr1	genetrack	.	3378	3398	525	-	.	stddev=0.0
+chr1	genetrack	.	3669	3689	845	+	.	stddev=0.0
+chr1	genetrack	.	3785	3805	23	-	.	stddev=0.0
+chr1	genetrack	.	3847	3867	316	-	.	stddev=0.0
+chr1	genetrack	.	3868	3888	491	+	.	stddev=0.0
+chr1	genetrack	.	4097	4117	536	-	.	stddev=0.0
+chr1	genetrack	.	4326	4346	482	+	.	stddev=0.0
+chr1	genetrack	.	4395	4415	3	+	.	stddev=0.0
+chr1	genetrack	.	4461	4481	1110	+	.	stddev=0.0
+chr1	genetrack	.	4500	4520	125	-	.	stddev=0.0
+chr1	genetrack	.	4620	4640	147	+	.	stddev=0.0
+chr1	genetrack	.	4826	4846	1761	+	.	stddev=4.82408982772
+chr1	genetrack	.	4902	4922	710	+	.	stddev=0.0
+chr1	genetrack	.	5110	5130	828	+	.	stddev=0.0
+chr1	genetrack	.	5402	5422	282	-	.	stddev=0.0
+chr1	genetrack	.	5501	5521	75	+	.	stddev=0.0
+chr1	genetrack	.	5707	5727	2	+	.	stddev=0.0
+chr1	genetrack	.	5717	5737	737	-	.	stddev=0.36608362591
+chr1	genetrack	.	6086	6106	646	+	.	stddev=0.039314009595
+chr1	genetrack	.	6098	6118	230	-	.	stddev=0.0657945476105
+chr1	genetrack	.	6187	6207	329	-	.	stddev=0.0
+chr1	genetrack	.	6290	6310	5	+	.	stddev=0.0
+chr1	genetrack	.	6356	6376	285	+	.	stddev=0.0
+chr1	genetrack	.	6380	6400	34	-	.	stddev=0.0
+chr1	genetrack	.	6401	6421	1587	+	.	stddev=5.61831543503
+chr1	genetrack	.	6415	6435	953	-	.	stddev=3.52372902021
+chr1	genetrack	.	6432	6452	742	+	.	stddev=0.0
+chr1	genetrack	.	6496	6516	691	+	.	stddev=0.0
+chr1	genetrack	.	6506	6526	61	-	.	stddev=1.5137105198
+chr1	genetrack	.	6843	6863	28	+	.	stddev=0.0
+chr1	genetrack	.	7058	7078	518	-	.	stddev=0.0
+chr1	genetrack	.	7124	7144	654	+	.	stddev=0.0
+chr1	genetrack	.	7765	7785	714	+	.	stddev=0.0
+chr1	genetrack	.	7847	7867	3	+	.	stddev=0.0
+chr1	genetrack	.	8209	8229	17	+	.	stddev=0.0
+chr1	genetrack	.	8272	8292	2	-	.	stddev=0.0
+chr1	genetrack	.	8459	8479	10	+	.	stddev=0.0
+chr1	genetrack	.	8471	8491	5	-	.	stddev=0.0
+chr1	genetrack	.	8715	8735	5	+	.	stddev=0.0
+chr1	genetrack	.	8834	8854	332	+	.	stddev=0.0
+chr1	genetrack	.	8839	8859	593	-	.	stddev=0.0
+chr1	genetrack	.	9034	9054	24	+	.	stddev=0.0
+chr1	genetrack	.	9058	9078	4	+	.	stddev=0.0
+chr1	genetrack	.	9485	9505	36	+	.	stddev=0.0
+chr1	genetrack	.	9710	9730	480	+	.	stddev=0.0
+chr1	genetrack	.	9923	9943	606	-	.	stddev=0.0
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/genetrack_input3.scidx	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,100 @@
+chrom	index	forward	reverse	value
+chr01	7	2	0	2
+chr01	8	3	0	3
+chr01	10	2	0	2
+chr01	14	1	0	1
+chr01	15	1	0	1
+chr01	16	1	0	1
+chr01	31	1	0	1
+chr01	32	1	0	1
+chr01	35	1	0	1
+chr01	37	2	0	2
+chr01	38	65	0	65
+chr01	39	8	4	12
+chr01	40	3	0	3
+chr01	41	3	2	5
+chr01	42	3	1	4
+chr01	43	2	2	4
+chr01	44	0	1	1
+chr01	46	3	0	3
+chr01	52	0	1	1
+chr01	54	7	0	7
+chr01	55	2	0	2
+chr01	56	1	0	1
+chr01	58	0	1	1
+chr01	61	0	1	1
+chr01	62	1	3	4
+chr01	63	0	1	1
+chr01	64	2	4	6
+chr01	66	0	13	13
+chr01	67	0	4	4
+chr01	72	0	5	5
+chr01	75	4	0	4
+chr01	77	0	6	6
+chr01	81	13	0	13
+chr01	83	1	0	1
+chr01	85	0	4	4
+chr01	86	0	1	1
+chr01	89	4	0	4
+chr01	90	1	0	1
+chr01	91	2	0	2
+chr01	92	4	1	5
+chr01	93	4	0	4
+chr01	96	1	0	1
+chr01	97	1	0	1
+chr01	98	4	0	4
+chr01	99	8	0	8
+chr01	102	1	1	2
+chr01	105	0	1	1
+chr01	107	0	3	3
+chr01	109	0	4	4
+chr01	110	2	1	3
+chr01	111	6	0	6
+chr01	112	0	3	3
+chr01	113	24	0	24
+chr01	114	0	32	32
+chr01	115	5	4	9
+chr01	116	36	5	41
+chr01	117	194	5	199
+chr01	118	100	2	102
+chr01	119	98	19	117
+chr01	120	92	7	99
+chr01	121	112	11	123
+chr01	122	36	3	39
+chr01	123	18	16	34
+chr01	124	9	14	23
+chr01	125	14	12	26
+chr01	126	20	1	21
+chr01	127	4	2	6
+chr01	128	66	6	72
+chr01	129	3	6	9
+chr01	130	2	10	12
+chr01	131	62	0	62
+chr01	132	20	1	21
+chr01	133	64	1	65
+chr01	134	58	0	58
+chr01	135	83	4	87
+chr01	136	92	1	93
+chr01	137	12	0	12
+chr01	138	0	3	3
+chr01	139	0	2	2
+chr01	142	0	48	48
+chr01	143	0	15	15
+chr01	144	0	74	74
+chr01	145	23	11	34
+chr01	146	1	120	121
+chr01	147	0	574	574
+chr01	148	0	11	11
+chr01	150	1	0	1
+chr01	151	3	4	7
+chr01	152	0	2	2
+chr01	153	1	3	4
+chr01	155	8	11	19
+chr01	156	0	14	14
+chr01	158	0	3	3
+chr01	159	0	1	1
+chr01	160	0	1	1
+chr01	161	0	5	5
+chr01	163	0	2	2
+chr01	164	0	29	29
+chr01	165	0	1	1
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/genetrack_input_unsorted4.gff	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,100 @@
+chr1	genetrack	.	1	19	450	+	.	stddev=4.77582440689
+chr1	genetrack	.	31	51	583	+	.	stddev=5.04887975977
+chr1	genetrack	.	34	54	88	-	.	stddev=2.896077962
+chr1	genetrack	.	55	75	689	-	.	stddev=4.30660216862
+chr1	genetrack	.	65	85	688	+	.	stddev=5.42600427981
+chr1	genetrack	.	85	105	676	-	.	stddev=4.49818285888
+chr1	genetrack	.	111	131	5134	+	.	stddev=6.21477198515
+chr1	genetrack	.	135	155	4059	-	.	stddev=4.75030601387
+chr1	genetrack	.	142	162	818	+	.	stddev=4.61530674023
+chr1	genetrack	.	162	182	509	+	.	stddev=5.54002911575
+chr1	genetrack	.	173	193	1844	-	.	stddev=4.57155440482
+chr1	genetrack	.	197	217	513	+	.	stddev=6.03335915224
+chr1	genetrack	.	201	221	889	-	.	stddev=4.25878183562
+chr1	genetrack	.	218	238	1213	-	.	stddev=4.50854580513
+chr1	genetrack	.	237	257	585	+	.	stddev=4.99399478661
+chr1	genetrack	.	245	265	1605	-	.	stddev=5.11855734302
+chr1	genetrack	.	263	283	1914	-	.	stddev=5.30293908946
+chr1	genetrack	.	281	301	144	+	.	stddev=5.80110142799
+chr1	genetrack	.	292	312	919	-	.	stddev=4.83235586146
+chr1	genetrack	.	317	337	227	+	.	stddev=5.39529430562
+chr1	genetrack	.	341	361	399	-	.	stddev=4.78763133608
+chr1	genetrack	.	363	383	494	+	.	stddev=5.06265712248
+chr1	genetrack	.	369	389	643	-	.	stddev=5.18741859391
+chr1	genetrack	.	394	414	548	-	.	stddev=5.35724250762
+chr1	genetrack	.	395	415	286	+	.	stddev=6.40676329916
+chr1	genetrack	.	423	443	462	-	.	stddev=4.14819371578
+chr1	genetrack	.	424	444	225	+	.	stddev=5.43077649451
+chr1	genetrack	.	442	462	325	+	.	stddev=4.88036859603
+chr1	genetrack	.	452	472	446	-	.	stddev=5.2172377816
+chr1	genetrack	.	464	484	288	+	.	stddev=5.3671065623
+chr1	genetrack	.	467	487	899	-	.	stddev=3.96802149366
+chr1	genetrack	.	486	506	216	-	.	stddev=4.58116999076
+chr1	genetrack	.	489	509	14	+	.	stddev=2.56845065662
+chr1	genetrack	.	521	541	212	-	.	stddev=1.41926272859
+chr1	genetrack	.	531	551	22	+	.	stddev=4.29226673437
+chr1	genetrack	.	561	581	15	+	.	stddev=2.31516738056
+chr1	genetrack	.	563	583	9	-	.	stddev=3.97523196
+chr1	genetrack	.	579	599	18	+	.	stddev=4.21819983838
+chr1	genetrack	.	602	622	35	-	.	stddev=3.0751140679
+chr1	genetrack	.	603	623	38	+	.	stddev=0.0
+chr1	genetrack	.	634	654	5	-	.	stddev=0.4
+chr1	genetrack	.	638	658	26	+	.	stddev=3.92251127144
+chr1	genetrack	.	650	670	63	-	.	stddev=1.3232803175
+chr1	genetrack	.	669	689	4	+	.	stddev=1.73205080757
+chr1	genetrack	.	689	709	29	-	.	stddev=2.02716392506
+chr1	genetrack	.	761	781	10	+	.	stddev=1.6
+chr1	genetrack	.	813	833	6	+	.	stddev=0.0
+chr1	genetrack	.	856	876	3	-	.	stddev=0.0
+chr1	genetrack	.	881	901	9	+	.	stddev=0.0
+chr1	genetrack	.	908	928	10	+	.	stddev=3.4292856399
+chr1	genetrack	.	923	943	13	+	.	stddev=2.30769230769
+chr1	genetrack	.	928	948	6	-	.	stddev=4.71404520791
+chr1	genetrack	.	948	968	14	+	.	stddev=0.0
+chr1	genetrack	.	956	976	2	-	.	stddev=0.0
+chr1	genetrack	.	986	1006	2	+	.	stddev=3.0
+chr1	genetrack	.	986	1006	5	-	.	stddev=0.0
+chr1	genetrack	.	1009	1029	11	+	.	stddev=0.0
+chr1	genetrack	.	1050	1070	2	-	.	stddev=0.0
+chr1	genetrack	.	1051	1071	6	+	.	stddev=0.0
+chr1	genetrack	.	1065	1085	5	-	.	stddev=0.0
+chr1	genetrack	.	1092	1112	7	-	.	stddev=0.0
+chr1	genetrack	.	1115	1135	9	-	.	stddev=0.99380799
+chr1	genetrack	.	1118	1138	3	+	.	stddev=0.0
+chr1	genetrack	.	1128	1148	10	-	.	stddev=0.0
+chr1	genetrack	.	1149	1169	5	+	.	stddev=0.0
+chr1	genetrack	.	1224	1244	4	-	.	stddev=0.0
+chr1	genetrack	.	1231	1251	2	+	.	stddev=0.0
+chr1	genetrack	.	1256	1276	7	+	.	stddev=0.0
+chr1	genetrack	.	1297	1317	8	+	.	stddev=0.0
+chr1	genetrack	.	1298	1318	4	-	.	stddev=0.866025403784
+chr1	genetrack	.	1323	1343	55	-	.	stddev=1.36980423286
+chr1	genetrack	.	1382	1402	5	-	.	stddev=0.0
+chr1	genetrack	.	1442	1462	4	-	.	stddev=0.433012701892
+chr1	genetrack	.	1460	1480	2	+	.	stddev=0.0
+chr1	genetrack	.	1522	1542	6	-	.	stddev=0.0
+chr1	genetrack	.	1547	1567	8	-	.	stddev=3.30718913883
+chr1	genetrack	.	1557	1577	7	+	.	stddev=0.451753951453
+chr1	genetrack	.	1604	1624	11	-	.	stddev=0.0
+chr1	genetrack	.	1690	1710	2	+	.	stddev=0.0
+chr1	genetrack	.	1714	1734	3	-	.	stddev=0.0
+chr1	genetrack	.	1747	1767	2	-	.	stddev=0.0
+chr1	genetrack	.	1765	1785	5	-	.	stddev=0.0
+chr1	genetrack	.	1834	1854	5	-	.	stddev=0.0
+chr1	genetrack	.	1859	1879	8	-	.	stddev=0.0
+chr1	genetrack	.	1870	1890	2	+	.	stddev=0.0
+chr1	genetrack	.	1915	1935	5	-	.	stddev=0.4
+chr1	genetrack	.	1935	1955	4	-	.	stddev=0.0
+chr1	genetrack	.	2005	2025	6	+	.	stddev=0.0
+chr1	genetrack	.	2056	2076	5	-	.	stddev=0.0
+chr1	genetrack	.	2079	2099	5	+	.	stddev=4.0
+chr1	genetrack	.	2088	2108	2	-	.	stddev=0.0
+chr1	genetrack	.	2120	2140	31	-	.	stddev=0.706738783878
+chr1	genetrack	.	2131	2151	16	+	.	stddev=0.7806247498
+chr1	genetrack	.	2191	2211	4	+	.	stddev=0.0
+chr1	genetrack	.	2227	2247	6	+	.	stddev=0.0
+chr1	genetrack	.	2235	2255	5	-	.	stddev=1.46969384567
+chr1	genetrack	.	2260	2280	6	+	.	stddev=0.0
+chr1	genetrack	.	2278	2298	2	-	.	stddev=0.0
+chr1	genetrack	.	2421	2441	23	-	.	stddev=3.58583645226
+chr1	genetrack	.	2446	2466	6	-	.	stddev=0.0
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/genetrack_output2.gff	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,92 @@
+chr1	genetrack	.	30	50	2060	+	.	stddev=0.0
+chr1	genetrack	.	63	83	397	+	.	stddev=0.0
+chr1	genetrack	.	79	99	521	+	.	stddev=0.0
+chr1	genetrack	.	113	133	5129	+	.	stddev=0.0
+chr1	genetrack	.	182	202	2538	+	.	stddev=0.0
+chr1	genetrack	.	228	248	2496	+	.	stddev=0.0
+chr1	genetrack	.	244	264	1525	+	.	stddev=0.0
+chr1	genetrack	.	271	291	15	+	.	stddev=0.0
+chr1	genetrack	.	298	318	1544	+	.	stddev=0.0
+chr1	genetrack	.	324	344	533	+	.	stddev=0.0
+chr1	genetrack	.	335	355	286	+	.	stddev=0.0
+chr1	genetrack	.	364	384	608	+	.	stddev=0.0
+chr1	genetrack	.	431	451	1393	+	.	stddev=0.0
+chr1	genetrack	.	473	493	58	+	.	stddev=0.0
+chr1	genetrack	.	747	767	23	+	.	stddev=0.0
+chr1	genetrack	.	789	809	607	+	.	stddev=0.0
+chr1	genetrack	.	834	854	665	+	.	stddev=0.0
+chr1	genetrack	.	867	887	468	+	.	stddev=0.0
+chr1	genetrack	.	1082	1102	740	+	.	stddev=0.0
+chr1	genetrack	.	1173	1193	25	+	.	stddev=0.0
+chr1	genetrack	.	1474	1494	584	+	.	stddev=0.0
+chr1	genetrack	.	2065	2085	1181	+	.	stddev=0.0
+chr1	genetrack	.	2092	2112	481	+	.	stddev=0.0
+chr1	genetrack	.	2442	2462	1246	+	.	stddev=0.0
+chr1	genetrack	.	2592	2612	34	+	.	stddev=0.0
+chr1	genetrack	.	2823	2843	1062	+	.	stddev=0.0
+chr1	genetrack	.	3120	3140	17	+	.	stddev=0.0
+chr1	genetrack	.	3659	3679	845	+	.	stddev=0.0
+chr1	genetrack	.	3858	3878	491	+	.	stddev=0.0
+chr1	genetrack	.	4316	4336	482	+	.	stddev=0.0
+chr1	genetrack	.	4451	4471	1110	+	.	stddev=0.0
+chr1	genetrack	.	4610	4630	147	+	.	stddev=0.0
+chr1	genetrack	.	4816	4836	1761	+	.	stddev=0.0
+chr1	genetrack	.	4892	4912	710	+	.	stddev=0.0
+chr1	genetrack	.	5100	5120	828	+	.	stddev=0.0
+chr1	genetrack	.	5491	5511	75	+	.	stddev=0.0
+chr1	genetrack	.	6076	6096	646	+	.	stddev=0.0
+chr1	genetrack	.	6280	6300	5	+	.	stddev=0.0
+chr1	genetrack	.	6346	6366	285	+	.	stddev=0.0
+chr1	genetrack	.	6391	6411	1587	+	.	stddev=0.0
+chr1	genetrack	.	6422	6442	742	+	.	stddev=0.0
+chr1	genetrack	.	6486	6506	691	+	.	stddev=0.0
+chr1	genetrack	.	6833	6853	28	+	.	stddev=0.0
+chr1	genetrack	.	7114	7134	654	+	.	stddev=0.0
+chr1	genetrack	.	7755	7775	714	+	.	stddev=0.0
+chr1	genetrack	.	8199	8219	17	+	.	stddev=0.0
+chr1	genetrack	.	8449	8469	10	+	.	stddev=0.0
+chr1	genetrack	.	8705	8725	5	+	.	stddev=0.0
+chr1	genetrack	.	8824	8844	332	+	.	stddev=0.0
+chr1	genetrack	.	9024	9044	24	+	.	stddev=0.0
+chr1	genetrack	.	9048	9068	4	+	.	stddev=0.0
+chr1	genetrack	.	9475	9495	36	+	.	stddev=0.0
+chr1	genetrack	.	9700	9720	480	+	.	stddev=0.0
+chr1	genetrack	.	21	41	245	-	.	stddev=0.0
+chr1	genetrack	.	52	72	1300	-	.	stddev=0.0
+chr1	genetrack	.	115	135	4659	-	.	stddev=0.0
+chr1	genetrack	.	145	165	897	-	.	stddev=0.0
+chr1	genetrack	.	161	181	956	-	.	stddev=0.0
+chr1	genetrack	.	174	194	494	-	.	stddev=0.0
+chr1	genetrack	.	196	216	2087	-	.	stddev=0.0
+chr1	genetrack	.	232	252	5047	-	.	stddev=0.0
+chr1	genetrack	.	292	312	626	-	.	stddev=0.0
+chr1	genetrack	.	334	354	726	-	.	stddev=0.0
+chr1	genetrack	.	348	368	792	-	.	stddev=0.0
+chr1	genetrack	.	379	399	126	-	.	stddev=0.0
+chr1	genetrack	.	429	449	618	-	.	stddev=0.0
+chr1	genetrack	.	451	471	754	-	.	stddev=0.0
+chr1	genetrack	.	528	548	1015	-	.	stddev=0.0
+chr1	genetrack	.	718	738	39	-	.	stddev=0.0
+chr1	genetrack	.	893	913	107	-	.	stddev=0.0
+chr1	genetrack	.	1117	1137	940	-	.	stddev=0.0
+chr1	genetrack	.	1281	1301	454	-	.	stddev=0.0
+chr1	genetrack	.	1319	1339	207	-	.	stddev=0.0
+chr1	genetrack	.	2115	2135	199	-	.	stddev=0.0
+chr1	genetrack	.	2828	2848	1144	-	.	stddev=0.0
+chr1	genetrack	.	3001	3021	1212	-	.	stddev=0.0
+chr1	genetrack	.	3106	3126	555	-	.	stddev=0.0
+chr1	genetrack	.	3368	3388	525	-	.	stddev=0.0
+chr1	genetrack	.	3775	3795	23	-	.	stddev=0.0
+chr1	genetrack	.	3837	3857	316	-	.	stddev=0.0
+chr1	genetrack	.	4087	4107	536	-	.	stddev=0.0
+chr1	genetrack	.	4490	4510	125	-	.	stddev=0.0
+chr1	genetrack	.	5392	5412	282	-	.	stddev=0.0
+chr1	genetrack	.	5707	5727	737	-	.	stddev=0.0
+chr1	genetrack	.	6088	6108	230	-	.	stddev=0.0
+chr1	genetrack	.	6177	6197	329	-	.	stddev=0.0
+chr1	genetrack	.	6370	6390	34	-	.	stddev=0.0
+chr1	genetrack	.	6405	6425	953	-	.	stddev=0.0
+chr1	genetrack	.	6496	6516	61	-	.	stddev=0.0
+chr1	genetrack	.	7048	7068	518	-	.	stddev=0.0
+chr1	genetrack	.	8461	8481	5	-	.	stddev=0.0
+chr1	genetrack	.	8829	8849	593	-	.	stddev=0.0
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/genetrack_output3.gff	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,10 @@
+chr1	genetrack	.	1	20	10	+	.	stddev=3.25729949498
+chr1	genetrack	.	28	48	92	+	.	stddev=2.09743181697
+chr1	genetrack	.	73	93	33	+	.	stddev=6.09264259969
+chr1	genetrack	.	85	105	30	+	.	stddev=3.96106046407
+chr1	genetrack	.	109	129	839	+	.	stddev=3.6259952943
+chr1	genetrack	.	121	141	675	+	.	stddev=5.59352633853
+chr1	genetrack	.	31	51	10	-	.	stddev=1.84390889146
+chr1	genetrack	.	56	76	32	-	.	stddev=3.28764258246
+chr1	genetrack	.	110	130	159	-	.	stddev=5.21810978727
+chr1	genetrack	.	136	156	893	-	.	stddev=2.27140632107
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/genetrack_output4.gff	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,84 @@
+chr1	genetrack	.	21	41	583	+	.	stddev=0.0
+chr1	genetrack	.	55	75	688	+	.	stddev=0.0
+chr1	genetrack	.	101	121	5134	+	.	stddev=0.0
+chr1	genetrack	.	132	152	818	+	.	stddev=0.0
+chr1	genetrack	.	152	172	509	+	.	stddev=0.0
+chr1	genetrack	.	187	207	513	+	.	stddev=0.0
+chr1	genetrack	.	227	247	585	+	.	stddev=0.0
+chr1	genetrack	.	271	291	144	+	.	stddev=0.0
+chr1	genetrack	.	307	327	227	+	.	stddev=0.0
+chr1	genetrack	.	353	373	494	+	.	stddev=0.0
+chr1	genetrack	.	385	405	286	+	.	stddev=0.0
+chr1	genetrack	.	414	434	225	+	.	stddev=0.0
+chr1	genetrack	.	432	452	325	+	.	stddev=0.0
+chr1	genetrack	.	454	474	288	+	.	stddev=0.0
+chr1	genetrack	.	479	499	14	+	.	stddev=0.0
+chr1	genetrack	.	521	541	22	+	.	stddev=0.0
+chr1	genetrack	.	551	571	15	+	.	stddev=0.0
+chr1	genetrack	.	569	589	18	+	.	stddev=0.0
+chr1	genetrack	.	593	613	38	+	.	stddev=0.0
+chr1	genetrack	.	628	648	26	+	.	stddev=0.0
+chr1	genetrack	.	659	679	4	+	.	stddev=0.0
+chr1	genetrack	.	751	771	10	+	.	stddev=0.0
+chr1	genetrack	.	803	823	6	+	.	stddev=0.0
+chr1	genetrack	.	871	891	9	+	.	stddev=0.0
+chr1	genetrack	.	898	918	10	+	.	stddev=0.0
+chr1	genetrack	.	913	933	13	+	.	stddev=0.0
+chr1	genetrack	.	938	958	14	+	.	stddev=0.0
+chr1	genetrack	.	999	1019	11	+	.	stddev=0.0
+chr1	genetrack	.	1041	1061	6	+	.	stddev=0.0
+chr1	genetrack	.	1139	1159	5	+	.	stddev=0.0
+chr1	genetrack	.	1246	1266	7	+	.	stddev=0.0
+chr1	genetrack	.	1287	1307	8	+	.	stddev=0.0
+chr1	genetrack	.	1547	1567	7	+	.	stddev=0.0
+chr1	genetrack	.	1995	2015	6	+	.	stddev=0.0
+chr1	genetrack	.	2069	2089	5	+	.	stddev=0.0
+chr1	genetrack	.	2121	2141	16	+	.	stddev=0.0
+chr1	genetrack	.	2181	2201	4	+	.	stddev=0.0
+chr1	genetrack	.	2217	2237	6	+	.	stddev=0.0
+chr1	genetrack	.	2250	2270	6	+	.	stddev=0.0
+chr1	genetrack	.	24	44	88	-	.	stddev=0.0
+chr1	genetrack	.	45	65	689	-	.	stddev=0.0
+chr1	genetrack	.	75	95	676	-	.	stddev=0.0
+chr1	genetrack	.	125	145	4059	-	.	stddev=0.0
+chr1	genetrack	.	163	183	1844	-	.	stddev=0.0
+chr1	genetrack	.	191	211	889	-	.	stddev=0.0
+chr1	genetrack	.	208	228	1213	-	.	stddev=0.0
+chr1	genetrack	.	235	255	1605	-	.	stddev=0.0
+chr1	genetrack	.	253	273	1914	-	.	stddev=0.0
+chr1	genetrack	.	282	302	919	-	.	stddev=0.0
+chr1	genetrack	.	331	351	399	-	.	stddev=0.0
+chr1	genetrack	.	359	379	643	-	.	stddev=0.0
+chr1	genetrack	.	384	404	548	-	.	stddev=0.0
+chr1	genetrack	.	413	433	462	-	.	stddev=0.0
+chr1	genetrack	.	442	462	446	-	.	stddev=0.0
+chr1	genetrack	.	457	477	899	-	.	stddev=0.0
+chr1	genetrack	.	476	496	216	-	.	stddev=0.0
+chr1	genetrack	.	511	531	212	-	.	stddev=0.0
+chr1	genetrack	.	553	573	9	-	.	stddev=0.0
+chr1	genetrack	.	592	612	35	-	.	stddev=0.0
+chr1	genetrack	.	640	660	63	-	.	stddev=0.0
+chr1	genetrack	.	679	699	29	-	.	stddev=0.0
+chr1	genetrack	.	918	938	6	-	.	stddev=0.0
+chr1	genetrack	.	976	996	5	-	.	stddev=0.0
+chr1	genetrack	.	1055	1075	5	-	.	stddev=0.0
+chr1	genetrack	.	1082	1102	7	-	.	stddev=0.0
+chr1	genetrack	.	1106	1126	9	-	.	stddev=0.0
+chr1	genetrack	.	1117	1137	10	-	.	stddev=0.0
+chr1	genetrack	.	1214	1234	4	-	.	stddev=0.0
+chr1	genetrack	.	1288	1308	4	-	.	stddev=0.0
+chr1	genetrack	.	1313	1333	55	-	.	stddev=0.0
+chr1	genetrack	.	1372	1392	5	-	.	stddev=0.0
+chr1	genetrack	.	1432	1452	4	-	.	stddev=0.0
+chr1	genetrack	.	1512	1532	6	-	.	stddev=0.0
+chr1	genetrack	.	1537	1557	8	-	.	stddev=0.0
+chr1	genetrack	.	1594	1614	11	-	.	stddev=0.0
+chr1	genetrack	.	1755	1775	5	-	.	stddev=0.0
+chr1	genetrack	.	1824	1844	5	-	.	stddev=0.0
+chr1	genetrack	.	1849	1869	8	-	.	stddev=0.0
+chr1	genetrack	.	1905	1925	5	-	.	stddev=0.0
+chr1	genetrack	.	1925	1945	4	-	.	stddev=0.0
+chr1	genetrack	.	2046	2066	5	-	.	stddev=0.0
+chr1	genetrack	.	2110	2130	31	-	.	stddev=0.0
+chr1	genetrack	.	2225	2245	5	-	.	stddev=0.0
+chr1	genetrack	.	2411	2431	23	-	.	stddev=0.0
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/tool_dependencies.xml	Tue Dec 22 17:03:27 2015 -0500
@@ -0,0 +1,6 @@
+<?xml version="1.0"?>
+<tool_dependency>
+    <package name="anaconda" version="2.3.0">
+        <repository changeset_revision="94d978ebbfd4" name="package_anaconda_2_3_0" owner="iuc" toolshed="https://toolshed.g2.bx.psu.edu" />
+    </package>
+</tool_dependency>