more changes to workflow

Author: dmnfarrell

README.md

@@ -5,9 +5,9 @@
 
 <img align="right" src=snpgenie/logo.png width=180px>
 
-_snpgenie_ is a tool for microbial variant calling and phylogenetic analysis from raw read data. It was primarily written to be used with bacterial isolates of M. Bovis but can be applied to other species. This is in early stages of development. Anyone interested in using the software is encouraged to make sugggestions on improving or adding features.
+_snpgenie_ is a tool for microbial variant calling and phylogenetic analysis from raw read data. It was primarily written to be used with bacterial isolates of M. Bovis but can be applied to other species. Anyone interested in using the software is encouraged to make suggestions on improving or adding features.
 
-This software is written in Python and is developed with the Qt toolkit using PySide2. It was made on Ubuntu linux but is designed to also run on Windows 10 with a standalone application.
+This software is written in Python. It was developed on Ubuntu linux but is designed to also run on Windows 10 with a standalone application. The GUI is made using the Qt toolkit using PySide2.
 
 ## Documentation
 
@@ -19,7 +19,7 @@ Note for Windows users: a standalone installer will be available.
 
 `pip install -e git+https://github.com/dmnfarrell/snpgenie.git#egg=snpgenie`
 
-(You may need to use pip3 on Ubuntu to ensure you use Python 3)
+Notes: You may need to use pip3 on Ubuntu to ensure you use Python 3. Use sudo if installing system-wide.
 
 ## Usage
 
@@ -34,6 +34,7 @@ For Linux installs, you require Python 3 and the following packages. These will
 * matplotlib
 * biopython
 * pyvcf
+* pyfaidx
 * pyside2 (if using GUI)
 
 Other binaries required:
@@ -42,7 +43,15 @@ Other binaries required:
 * samtools
 * bcftools
 
-The binaries can be installed with apt in Ubuntu. They are downloaded automatically in Windows.
+These binaries can be installed with apt in Ubuntu:
+
+`sudo apt install bwa samtools bcftools`
+
+If you want a tree to be built you should install RaXML:
+
+`sudo apt install raxml`
+
+The binaries are downloaded automatically in Windows.
 
 ## Screenshots
 

doc/source/description.rst

@@ -5,16 +5,13 @@ Introduction
 
 This software is written in Python and is developed with the Qt toolkit using PySide2. It was made on Ubuntu linux but is designed to also run on Windows 10 with a standalone application.
 
-The Desktop application
-=======================
+Command line tool
+-----------------
 
-Unlike many other SNP calling pipelines, this tool is designed around a graphical user interface. Though it will also work from the command line and via Python scripts.
-
-.. image:: scr1.png
-     :scale: 45%
+This tool works from the command line and via Python scripts. Unlike many other SNP calling pipelines, it is also designed to have a graphical user interface, which is in development.
 
 Current Features
-================
+----------------
 
 * load multiple fastq files and process together
 * view fastq quality statistics
@@ -27,7 +24,7 @@ Current Features
 * create phylogenetic tree
 
 Links
-=====
+-----
 
 http://dmnfarrell.github.io/snpgenie
 
@@ -37,14 +34,13 @@ Installation
 Linux
 -----
 
-Install dependencies::
+With pip::
 
-  pip install pandas matplotlib biopython pyside2
-  sudo apt install bcftools samtools bwa
+  pip install -e git+https://github.com/dmnfarrell/snpgenie.git#egg=snpgenie
 
-Get the current version from github::
+Install binary dependencies::
 
-  git clone https://github.com/dmnfarrell/btbgenie.git
+  sudo apt install bcftools samtools bwa
 
 Windows
 -------
@@ -54,4 +50,4 @@ A standalone installer will be used to deploy on windows.
 Mac OSX
 -------
 
-Use the Linux instructions.
+Not tested. You can try the Linux instructions possibly with bioconda for the binaries.

doc/source/usage.rst

@@ -12,7 +12,7 @@ This will run the entire process based on a set of options given at the terminal
   -i FILE, --input FILE
                         input folder(s)
   -l FILE, --labels FILE
-                        sample labels file
+                        sample labels file, optional
   -r FILE, --reference FILE
                         reference genome filename
   -w, --overwrite       overwrite intermediate files
@@ -26,6 +26,23 @@ This will run the entire process based on a set of options given at the terminal
   -v, --version         Get version
   -s, --test            Do test run
 
+Example::
+
+  snpgenie -r reference.fa -g reference.gff -i data_files -t 8 -o results
+
+From Python
+-----------
+
+You can run a workflow from within Python::
+
+  from sngenie import app
+  args = {'threads':8, 'outdir': 'results', 'labelsep':'-',
+          'input':['/my/folder/',
+                   '/my/other/folder'],
+          'reference': None, 'overwrite':False}
+  W = app.WorkFlow(**args)
+  st = W.setup()
+  W.run()
 
 Desktop Application
 -------------------

notebooks/albania.ipynb

@@ -12,7 +12,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 2,
+   "execution_count": 1,
    "metadata": {},
    "outputs": [],
    "source": [
@@ -31,7 +31,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 3,
+   "execution_count": 2,
    "metadata": {},
    "outputs": [],
    "source": [
@@ -54,6 +54,8 @@
       "reference None\n",
       "overwrite False\n",
       "quality 25\n",
+      "filters 'QUAL>=40 && INFO/DP>=10 && MQ>40'\n",
+      "gff_file None\n",
       "14 samples were loaded:\n",
       "-------------\n",
       "                   name      sample                                           filename  pair  read_length\n",
@@ -108,11 +110,6 @@
       "14 samples\n",
       "using 4020 sites\n",
       "\n",
-      "building tree\n",
-      "-------------\n",
-      "raxmlHPC-PTHREADS -f a -N 10 -T 8 -m GTRCAT -V -p 77901206 -x 18311491 -n variants -w /home/damien/gitprojects/snpgenie/test_results -s ../test_results/core.fa\n",
-      "/home/damien/gitprojects/snpgenie/test_results/RAxML_bipartitions.variants\n",
-      "\n",
       "Done. Sample summary:\n",
       "---------------------\n",
       "        sample                 name                               bam_file  read_length\n",
@@ -129,7 +126,12 @@
       "10  SRR5486071         SRR5486071_1  ../test_results/mapped/SRR5486071.bam          227\n",
       "11  SRR8063654         SRR8063654_1  ../test_results/mapped/SRR8063654.bam          217\n",
       "12  SRR8063665         SRR8063665_1  ../test_results/mapped/SRR8063665.bam          231\n",
-      "13  SRR8065079         SRR8065079_1  ../test_results/mapped/SRR8065079.bam          233\n"
+      "13  SRR8065079         SRR8065079_1  ../test_results/mapped/SRR8065079.bam          233\n",
+      "\n",
+      "building tree\n",
+      "-------------\n",
+      "raxmlHPC-PTHREADS -f a -N 10 -T 8 -m GTRCAT -V -p 72177716 -x 81946655 -n variants -w /home/damien/gitprojects/snpgenie/test_results -s ../test_results/core.fa\n",
+      "/home/damien/gitprojects/snpgenie/test_results/RAxML_bipartitions.variants\n"
      ]
     }
    ],

setup.py

@@ -22,6 +22,7 @@
                       'matplotlib>=3.0',
                       'biopython>=1.5',
                       'pyvcf>=0.6',
+                      'pyfaidx',
                       'pyside2>=5.1',
                       'future'],
     entry_points = {

snpgenie/app.py

@@ -42,19 +42,20 @@
 module_path = os.path.dirname(os.path.abspath(__file__)) #path to module
 datadir = os.path.join(module_path, 'data')
 sequence_path = os.path.join(config_path, 'genome')
+annotation_path = os.path.join(config_path, 'annotation')
 ref_genome = os.path.join(sequence_path, 'Mbovis_AF212297.fa')
 ref_gff = os.path.join(datadir, 'Mbovis_csq_format.gff')
 #windows only path to binaries
 bin_path = os.path.join(config_path, 'binaries')
-default_filter = "'QUAL>=40 && INFO/DP>=10 && MQ>40'"
+default_filter = 'QUAL>=40 && INFO/DP>=10 && MQ>40'
 
 if not os.path.exists(config_path):
     try:
         os.makedirs(config_path, exist_ok=True)
     except:
         os.makedirs(config_path)
 
-defaults = {'threads':4, 'labelsep':'_','quality':25,
+defaults = {'threads':4, 'labelsep':'_','quality':25, 'filters': default_filter,
             'reference': None, 'gff_file': None, 'overwrite':False}
 
 def check_platform():
@@ -285,11 +286,11 @@ def mpileup_gnuparallel(bam_files, ref, outpath, threads=4, callback=None):
     return rawbcf
 
 def variant_calling(bam_files, ref, outpath, relabel=True, threads=4,
-                    callback=None, overwrite=False, filter=None, gff_file=None, **kwargs):
+                    callback=None, overwrite=False, filters=None, gff_file=None, **kwargs):
     """Call variants with bcftools"""
 
-    if filter == None:
-        filter = default_filter
+    if filters == None:
+        filters = default_filter
     rawbcf = os.path.join(outpath,'raw.bcf')
     bcftoolscmd = tools.get_cmd('bcftools')
     if not os.path.exists(rawbcf) or overwrite == True:
@@ -321,11 +322,11 @@ def variant_calling(bam_files, ref, outpath, relabel=True, threads=4,
 
     #filter variants
     final = os.path.join(outpath,'filtered.vcf.gz')
-    cmd = '{bc} filter -i {f} -o {o} -O z {i}'.format(bc=bcftoolscmd,i=vcfout,o=final,f=filter)
+    cmd = '{bc} filter -i "{f}" -o {o} -O z {i}'.format(bc=bcftoolscmd,i=vcfout,o=final,f=filters)
     print (cmd)
     tmp = subprocess.check_output(cmd,shell=True)
     if callback != None:
-        callback(tmp)
+        callback(cmd)
 
     #consequence calling
     if gff_file != None:
@@ -399,6 +400,9 @@ def setup(self):
         self.filenames = get_files_from_paths(self.input)
         self.threads = int(self.threads)
         df = get_samples(self.filenames, sep=self.labelsep)
+        if len(df) == 0:
+            print ('no samples provided')
+            return
         df['read_length'] = df.filename.apply(tools.get_fastq_info)
         self.fastq_table = df
         sample_size = len(df['sample'].unique())
@@ -444,7 +448,7 @@ def run(self):
         print ('----------------')
         bam_files = list(samples.bam_file.unique())
         self.vcf_file = variant_calling(bam_files, self.reference, self.outdir, threads=self.threads,
-                                        gff_file=self.gff_file,
+                                        gff_file=self.gff_file, filters=self.filters,
                                         overwrite=self.overwrite)
         print (self.vcf_file)
         print ()
@@ -468,6 +472,8 @@ def run(self):
             print ('building tree')
             print ('-------------')
             treefile = trees.run_RAXML(outfasta, outpath=self.outdir)
+            if treefile == None:
+                return
             print (treefile)
             #labelmap = dict(zip(sra.filename,sra.geo_loc_name_country))
             t,ts = trees.create_tree(treefile)#, labelmap)
@@ -491,25 +497,25 @@ def main():
 
     import sys, os
     from argparse import ArgumentParser
-    parser = ArgumentParser(description='snpgenie CLI tool. https://github.com/dmnfarrell/btbgenie')
-    #parser.add_argument("-f", "--fasta", dest="filenames",default=[],
-    #                    help="input fasta file", metavar="FILE")
+    parser = ArgumentParser(description='snpgenie CLI tool. https://github.com/dmnfarrell/snpgenie')
     parser.add_argument("-i", "--input", action='append', dest="input", default=[],
                         help="input folder(s)", metavar="FILE")
-    parser.add_argument("-l", "--labels", dest="labels", default=[],
-                        help="sample labels file", metavar="FILE")
+    #parser.add_argument("-l", "--labels", dest="labels", default=[],
+    #                    help="sample labels file, optional", metavar="FILE")
     parser.add_argument("-e", "--labelsep", dest="labelsep", default=',',
-                        help="symbol to split sample labels on")
+                        help="symbol to split the sample labels on")
     parser.add_argument("-r", "--reference", dest="reference", default=None,
                         help="reference genome filename", metavar="FILE")
     parser.add_argument("-g", "--gff", dest="gff_file", default=None,
                         help="reference gff, optional", metavar="FILE")
     parser.add_argument("-w", "--overwrite", dest="overwrite", action="store_true", default=False,
-                        help="overwrite intermediate files" )
+                        help="overwrite intermediate files")
     parser.add_argument("-m", "--trim", dest="trim", action="store_true", default=False,
-                        help="trim fastq files" )
+                        help="whether to trim fastq files" )
     parser.add_argument("-q", "--quality", dest="quality", default=25,
                         help="trim quality" )
+    parser.add_argument("-f", "--filters", dest="filters", default=default_filter,
+                        help="variant calling post-filters" )
     parser.add_argument("-t", "--threads", dest="threads",default=4,
                         help="cpu threads to use")
     parser.add_argument("-o", "--outdir", dest="outdir",
@@ -527,6 +533,11 @@ def main():
         from . import __version__
         print ('snpgenie version %s' %__version__)
         print ('https://github.com/dmnfarrell/btbgenie')
+    elif args['outdir'] == None:
+        print ('No input or output folders provided. These are required.')
+        print ('Example:')
+        print ('snpgenie -r <reference> -i <input folder with fastq.gz files> -o <output folder>')
+        print ('Use -h for more help on options.')
     else:
         W = WorkFlow(**args)
         st = W.setup()