getting ready for release

Author: Dom Laetsch

.gitignore

@@ -1,7 +1,3 @@
-testing/
-testing_fly/
-*.json
+nodesDB.txt
 *.pyc
-*.txt
-*.png
-*.tsv
+*.png

lib/BtCore.py bloblib/BtCore.py

@@ -1,10 +1,15 @@
 #!/usr/bin/env python
 # -*- coding: utf-8 -*-
 
+"""
+File        : BtCore.py
+Author      : Dominik R. Laetsch, dominik.laetsch at gmail dot com
+"""
+
 from __future__ import division
-import lib.BtLog as BtLog
-import lib.BtIO as BtIO
-import lib.BtTax as BtTax
+import bloblib.BtLog as BtLog
+import bloblib.BtIO as BtIO
+import bloblib.BtTax as BtTax
 from os.path import abspath, isfile, basename, isdir, join
 from os import getcwd, mkdir
 import json
@@ -39,7 +44,6 @@ def view(self, **kwargs):
         seqs = kwargs['seqs']
         cov_libs = kwargs['cov_libs']
         progress_bar = kwargs['progressbar']
-
         # Default sequences if no subset
         if not (seqs):
             seqs = self.order_of_blobs
@@ -219,9 +223,6 @@ def getPlotData(self, rank, min_length, hide_nohits, taxrule, c_index, catcolour
         read_cov_dict = {}
         max_cov = 0.0
 
-        if taxrule not in self.taxrules:
-            BtLog.error('11', taxrule, self.taxrules)
-
         cov_lib_dict = self.covLibs
         cov_lib_names_l = self.covLibs.keys() # does not include cov_sum
         if len(cov_lib_names_l) > 1:
@@ -236,8 +237,13 @@ def getPlotData(self, rank, min_length, hide_nohits, taxrule, c_index, catcolour
             if (catcolour_dict): # annotation with categories specified in catcolour
                 group = str(catcolour_dict[name])
             elif (c_index): # annotation with c_index instead of taxonomic group
-                group = str(blob['taxonomy'][taxrule][rank]['c_index'])
+                if taxrule not in self.taxrules:
+                    BtLog.error('11', taxrule, self.taxrules)
+                else:
+                    group = str(blob['taxonomy'][taxrule][rank]['c_index'])
             else: # annotation with taxonomic group
+                if taxrule not in self.taxrules:
+                    BtLog.warn_d['9'] % (taxrule, self.taxrules)
                 group = str(blob['taxonomy'][taxrule][rank]['tax'])
             if not group in data_dict:
                 data_dict[group] = {
@@ -360,6 +366,9 @@ def parseCoverage(self, **kwargs):
                     covLib.cov_sum += cov
                     self.dict_of_blobs[name].addCov(covLib.name, cov)
                     self.dict_of_blobs[name].addReadCov(covLib.name, read_cov_dict[name])
+                out_f = BtIO.getOutFile(covLib.f, None, None)
+                covView = ViewObj(name="covlib", out_f=out_f, suffix="cov", header="", body=[])
+                self.view(viewObjs=[covView], ranks=None, taxrule=None, hits_flag=None, seqs=None, cov_libs=[covLib.name], progressbar=False)
 
             elif covLib.fmt == 'cov':
                 base_cov_dict, covLib.reads_total, covLib.reads_mapped, covLib.reads_unmapped, read_cov_dict = BtIO.parseCov(covLib.f, set(self.dict_of_blobs))

lib/BtIO.py bloblib/BtIO.py

@@ -3,10 +3,7 @@
 
 """
 File        : BtIO.py
-Version     : 0.1
 Author      : Dominik R. Laetsch, dominik.laetsch at gmail dot com
-Bugs        : ?
-To do       : ?
 """
 
 from __future__ import division
@@ -16,7 +13,7 @@
 import shutil
 import os
 import sys
-import lib.BtLog as BtLog
+import bloblib.BtLog as BtLog
 from collections import deque
 
 

lib/BtLog.py bloblib/BtLog.py

@@ -3,11 +3,9 @@
 
 """
 File        : BtLog.py
-Version     : 0.1
 Author      : Dominik R. Laetsch, dominik.laetsch at gmail dot com
-Bugs        : ?
-To do       : ?
 """
+
 from __future__ import division
 import sys
 from os.path import basename
@@ -60,13 +58,13 @@ def progress(iteration, steps, max_value, no_limit=False):
     '22' : '[ERROR:22]\t: Tax file %s seems to have no taxids.',
     '23' : '[ERROR:23]\t: Catcolour file %s does not seem to have the right format.',
     '24' : '[ERROR:24]\t: Catcolour file incompatible with c-index colouring.',
-    '25' : '[ERROR:25]\t: Cov file %s does not seem to have the right format.',
+    '25' : '[ERROR:25]\t: COV file %s does not seem to have the right format.',
     '26' : '[ERROR:26]\t: Please specify type similarity search result (--rnacentral FILE or --diamond FILE) or a single TAXID (--taxid INT)',
     '27' : '[ERROR:27]\t: nodesDB ("--db") %s could not be read.',
     '28' : '[ERROR:28]\t: Please specify "--names" and "--nodes", or "--db"',
     '29' : '[ERROR:29]\t: No mapping reads found in %s',
     '30' : '[ERROR:30]\t: The module docopt is not installed. Please install it to run blobtools\n\tpip install docopt',
-    '31' : '[ERROR:31]\t: Please specify a read mapping file',
+    '31' : '[ERROR:31]\t: Please specify a read mapping file (BAM/SAM/CAS)',
     '32' : '[ERROR:32]\t: Choose either --cumulative or --multiplot',
     '33' : '[ERROR:33] : CovLib(s) not found. The available covlibs are: \n%s',
     '34' : '[ERROR:34] : Invalid plot type : %s',
@@ -75,6 +73,7 @@ def progress(iteration, steps, max_value, no_limit=False):
     '37' : '[ERROR:37] : %s does not seem to be a valid %s file',
     '38' : '[ERROR:38] : %s is not an integer',
     '39' : '[ERROR:39] : Please specify a taxid file (mapping subjects to taxids)',
+    '40' : '[ERROR:40] : CovLib \'%s\' not specified in refcov file'
 }
 
 warn_d = {
@@ -86,7 +85,9 @@ def progress(iteration, steps, max_value, no_limit=False):
     '5' : '[WARN]\t\t: Hit for sequence %s in tax file %s has multiple taxIds, only first one is used.',
     '6' : '[WARN]\t\t: Sum of coverage in cov lib %s is 0.0. Please ignore this warning if "--no_base_cov" was specified.',
     '7' : '[WARN]\t\t: No taxonomy information found.',
-    '8' : '[WARN]\t\t: Duplicated sequences found :\n\t\t\t%s'
+    '8' : '[WARN]\t\t: Duplicated sequences found :\n\t\t\t%s',
+    '9' : '[WARN]\t\t: Taxrule "%s" was not computed for this BlobDb. Available taxrule(s) : %s. Will proceed without taxonomic annotation ...'
+
 }
 status_d = {
     '0': '[STATUS]\t: Nothing to be done. %s',
@@ -98,15 +99,15 @@ def progress(iteration, steps, max_value, no_limit=False):
     '6' : '[STATUS]\t: Computing taxonomy using taxrule(s) %s',
     '7' : '[STATUS]\t: Generating BlobDB and writing to file %s',
     '8' : '[STATUS]\t: Plotting %s',
-    '9' : '[STATUS]\t: Reading BlobDb %s',
+    '9' : '[STATUS]\t: Reading BlobDB %s',
     '10': '[STATUS]\t: \tChecking with \'samtools flagstat\'',
     '11': '[STATUS]\t: \tMapping reads = %s, total reads = %s (mapping rate = %s)',
     '12': '[STATUS]\t: \tChecking with \'clc_mapping_info\'',
     '13': '[STATUS]\t: \tWriting %s',
     '14': '[STATUS]\t: Preparing view(s) ...',
-    '15': '[STATUS]\t: \tLoading BlobDb into memory ...',
-    '16': '[STATUS]\t: \tSerialising BlobDb (using \'%s\' module) (this may take a while) ...',
-    '17': '[STATUS]\t: \tSerialising BlobDb (using \'%s\' module) (this may take a while, consider installing the \'ujson\' module) ...',
+    '15': '[STATUS]\t: \tLoading BlobDB into memory ...',
+    '16': '[STATUS]\t: \tSerialising BlobDB (using \'%s\' module) (this may take a while) ...',
+    '17': '[STATUS]\t: \tSerialising BlobDB (using \'%s\' module) (this may take a while, consider installing the \'ujson\' module) ...',
     '18': '[STATUS]\t: Extracting data for plots ...',
     '19': '[STATUS]\t: Writing output ...',
     '20': '[STATUS]\t: \tFinished in %ss',

lib/BtPlot.py bloblib/BtPlot.py

@@ -3,17 +3,15 @@
 
 """
 File        : BtPlot.py
-Version     : 0.1
 Author      : Dominik R. Laetsch, dominik.laetsch at gmail dot com
-Bugs        : ?
-To do       : ?
 """
+
 from __future__ import division
 from numpy import array, arange, logspace, mean, std
 import math
-import lib.BtLog as BtLog
-import lib.BtIO as BtIO
-import lib.BtTax as BtTax
+import bloblib.BtLog as BtLog
+import bloblib.BtIO as BtIO
+import bloblib.BtTax as BtTax
 import matplotlib as mat
 from matplotlib import cm
 from matplotlib.ticker import NullFormatter, MultipleLocator, AutoMinorLocator
@@ -91,17 +89,13 @@ def set_format_scatterplot(axScatter, **kwargs):
         min_y, max_y = 0.005, kwargs['max_cov']+1000
         axScatter.set_yscale('log')
         axScatter.set_xscale('linear')
-        axScatter.set_xlabel("GC proportion")
-        axScatter.set_ylabel("Coverage")
         axScatter.xaxis.set_major_locator(major_xticks)
         axScatter.xaxis.set_minor_locator(minor_xticks)
     elif kwargs['plot'] == 'covplot':
         min_x, max_x = 0.005, kwargs['max_cov']+1000
         min_y, max_y = 0.005, kwargs['max_cov']+1000
         axScatter.set_yscale('log')
         axScatter.set_xscale('log')
-        axScatter.set_xlabel(kwargs['xlabel'])
-        axScatter.set_ylabel(kwargs['ylabel'])
     else:
         BtLog.error('34' % kwargs['plot'])
     axScatter.set_xlim( (min_x, max_x) )
@@ -211,12 +205,14 @@ def __init__(self, data_dict, cov_lib_dict, cov_lib_selection, plot_type):
         self.labels = {'all'}
         self.plot = plot_type # type of plot
         self.group_labels = {}
+        self.cov_lib_dict = cov_lib_dict
         self.cov_libs = self.subselect_cov_libs(cov_lib_dict, cov_lib_selection)
         self.cov_libs_total_reads_dict = self.get_cov_libs_total_reads_dict(cov_lib_dict)
         self.cov_libs_mapped_reads_dict = self.get_cov_libs_mapped_reads_dict(cov_lib_dict)
         self.data_dict = data_dict
         self.stats = {}
         self.exclude_groups = []
+        self.version = None
         self.colours = {}
         self.group_order = []
         self.plot_order = []
@@ -230,10 +226,8 @@ def __init__(self, data_dict, cov_lib_dict, cov_lib_selection, plot_type):
         self.cumulative_flag = ''
 
         self.cov_y_dict = {}
-        self.xlabel = ''
-        self.ylabel = ''
-        self.xmax = ''
-        self.ymax = ''
+        self.xlabel = None
+        self.ylabel = None
 
         self.refcov_dict = {}
 
@@ -298,18 +292,33 @@ def write_stats(self, out_f):
                 for g, labels in self.group_labels.items():
                     if label in labels:
                         stats.append(self.get_stats_for_group(g))
-
+        output = []
+        output.append('## %s' % self.version)
+        for cov_lib, cov_lib_dict in self.cov_lib_dict.items():
+           if cov_lib in self.cov_libs:
+                output.append("## %s=%s" % (cov_lib, cov_lib_dict['f']))
+        fields = ['name', 'colour', 'count_visible', 'count_visible_perc', 'span_visible','span_visible_perc', 'n50', 'gc_mean', 'gc_std']
+        header = [field for field in fields]
+        for cov_lib in sorted(self.cov_libs):
+            header.append('%s_mean' % cov_lib)
+            header.append('%s_std' % cov_lib)
+            header.append('%s_read_map' % cov_lib)
+            header.append('%s_read_map_p' % cov_lib)
+        output.append('# %s' % "\t".join(header))
+        for stat in stats:
+            line = []
+            for field in fields:
+                line.append(stat[field])
+            for cov_lib in sorted(self.cov_libs):
+                line.append(stat['cov_mean'][cov_lib])
+                line.append(stat['cov_std'][cov_lib])
+                line.append(stat['reads_mapped'][cov_lib])
+                line.append(stat['reads_mapped_perc'][cov_lib])
+            output.append("%s" % "\t".join(line))
         out_f = "%s.stats.txt" % out_f
         with open(out_f, 'w') as fh:
-            for cov_lib in sorted(self.cov_libs):
-                fh.write("# %s - %s\n" % (self.out_f, cov_lib))
-                fh.write("%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\n" % ('Group', 'colour', 'count', 'visible-%', 'span', 'visible-%', 'n50', 'GC-mean', 'GC-std', 'Cov-mean', 'Cov-std', 'Mapped-reads', 'Mapped-reads-%'))
-                for stat in stats:
-                    fh.write("%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\t%s\n" % (\
-                            stat['name'], stat['colour'], stat['count_visible'], stat['count_visible_perc'], stat['span_visible'], \
-                            stat['span_visible_perc'], stat['n50'], stat['gc_mean'], stat['gc_std'], stat['cov_mean'][cov_lib], \
-                            stat['cov_std'][cov_lib], stat['reads_mapped'][cov_lib], stat['reads_mapped_perc'][cov_lib]))
-
+            print BtLog.status_d['24'] % ("%s" % out_f)
+            fh.write("\n".join(output))
 
     def compute_stats(self):
         stats = {}
@@ -334,6 +343,7 @@ def compute_stats(self):
                             'count_hidden' : 0,
                             'span_hidden' : 0
                             }
+
         for group, labels in self.group_labels.items():
             for label in labels:
                 stats[label]['name'] = stats[label]['name'] + self.data_dict[group]['name']
@@ -372,10 +382,8 @@ def relabel_and_colour(self, colour_dict, user_labels):
             colour_groups = [group if not (group in user_labels) else user_labels[group] for group in groups]
             colour_dict = generateColourDict(colour_groups)
         for idx, group in enumerate(self.group_order):
-            if (self.exclude_groups):
-                if group in self.exclude_groups:
-                    self.group_labels[group].add('other')
-                    self.colours[group] = WHITE
+            if group in self.exclude_groups:
+                pass
             elif group in user_labels:
                 label = user_labels[group]
                 self.group_labels[group].add(label)
@@ -408,7 +416,9 @@ def setupPlot(self, plot):
             # Setting up plots and axes
             fig = plt.figure(1, figsize=(35,35), dpi=400)
             axScatter = plt.axes(rect_scatter, axisbg=BGGREY)
-            axScatter = set_format_scatterplot(axScatter, max_cov=self.max_cov, plot=plot, xlabel=self.xlabel, ylabel=self.ylabel)
+            axScatter = set_format_scatterplot(axScatter, max_cov=self.max_cov, plot=plot)
+            axScatter.set_xlabel(self.xlabel)
+            axScatter.set_ylabel(self.ylabel)
             axHistx = plt.axes(rect_histx, axisbg=BGGREY)
             axHistx = set_format_hist_x(axHistx, axScatter)
             axHisty = plt.axes(rect_histy, axisbg=BGGREY)
@@ -475,16 +485,20 @@ def plotBar(self, cov_lib, out_f):
         ax_main_data['labels'].append('Mapped (assembly)')
         ax_main_data['values'].append(reads_mapped/reads_total)
         ax_main_data['colours'].append(DGREY)
-        if cov_lib in self.refcov_dict:
-            reads_total_ref = refcov_dict[cov_lib]['reads_total']
-            reads_mapped_ref = refcov_dict[cov_lib]['reads_mapped']
-            reads_unmapped_ref = reads_total_ref - reads_mapped_ref
-            ax_main_data['labels'].append('Unmapped (ref)')
-            ax_main_data['values'].append(reads_unmapped_ref/reads_total_ref)
-            ax_main_data['colours'].append(DGREY)
-            ax_main_data['labels'].append('Mapped (ref)')
-            ax_main_data['values'].append(reads_mapped_ref/reads_total_ref)
-            ax_main_data['colours'].append(DGREY)
+        if (self.refcov_dict):
+            if cov_lib in self.refcov_dict:
+                reads_total_ref = self.refcov_dict[cov_lib]['reads_total']
+                reads_mapped_ref = self.refcov_dict[cov_lib]['reads_mapped']
+                reads_unmapped_ref = reads_total_ref - reads_mapped_ref
+                ax_main_data['labels'].append('Unmapped (ref)')
+                ax_main_data['values'].append(reads_unmapped_ref/reads_total_ref)
+                ax_main_data['colours'].append(DGREY)
+                ax_main_data['labels'].append('Mapped (ref)')
+                ax_main_data['values'].append(reads_mapped_ref/reads_total_ref)
+                ax_main_data['colours'].append(DGREY)
+            else:
+                BtLog.error('40', cov_lib)
+
         # mapped plotted groups
         for group in self.plot_order:
            ax_group_data['labels'].append(group)
@@ -524,7 +538,7 @@ def plotScatter(self, cov_lib, info_flag, out_f):
             if group == 'no-hit':
                 alpha = 0.5
             group_length_array = array(self.stats[group]['length'])
-            if len(group_length_array) > 0:
+            if len(group_length_array) > 0 and group not in self.exclude_groups:
                 colour = self.colours[group]
                 group_x_array = ''
                 group_y_array = ''
@@ -567,7 +581,7 @@ def plotScatter(self, cov_lib, info_flag, out_f):
                 axScatter.scatter(group_x_array, group_y_array, color = colour, s = marker_size_array, lw = lw, alpha=alpha, edgecolor=BLACK, label=label)
                 axLegend.axis('off')
                 if (self.multiplot):
-                    fig_m, axScatter_m, axHistx_m, axHisty_m, axLegend_m = self.setupPlot()
+                    fig_m, axScatter_m, axHistx_m, axHisty_m, axLegend_m, top_bins, right_bins = self.setupPlot(self.plot)
                     legend_handles_m = []
                     legend_labels_m = []
                     legend_handles_m.append(Line2D([0], [0], linewidth = 0.5, linestyle="none", marker="o", alpha=1, markersize=24, markerfacecolor=colour))

lib/BtTax.py bloblib/BtTax.py

@@ -3,11 +3,9 @@
 
 """
 File        : BtTax.py
-Version     : 0.1
 Author      : Dominik R. Laetsch, dominik.laetsch at gmail dot com
-Bugs        : ?
-To do       : ?
 """
+
 from __future__ import division
 RANKS = ['species', 'genus', 'family', 'order', 'phylum', 'superkingdom']
 TAXRULES = ['bestsum', 'bestsumorder'] # this should be re-named colour rules at one point