Merge pull request #142 from jaleezyy/patch_1.6.4

Patch 1.6.4 including modified output in place of error and breseq parameter update

Snakefile

@@ -45,7 +45,8 @@ validate(config, 'resources/config.schema.yaml')
 samples = pd.read_table(config['samples'], sep=',')
 validate(samples, 'resources/sample.schema.yaml')
 
-# manual assignment of breseq reference
+# manual assignment of breseq parameters
+# breseq reference, minor variant thresholds
 try:
     if os.path.exists(config['breseq_reference']):
         breseq_ref = config['breseq_reference']
@@ -54,6 +55,20 @@ try:
 except TypeError:
     breseq_ref = ""
 
+if "polymorphism_variant_coverage" in config.keys():
+    breseq_cov = config['polymorphism_variant_coverage']
+    if breseq_cov == "":
+        breseq_cov = 2
+else:
+    breseq_cov = 2
+
+if "polymorphism_frequency" in config.keys():
+    breseq_freq = config['polymorphism_frequency']
+    if breseq_freq == "":
+        breseq_freq = 0.05
+else:
+    breseq_freq = 0.05
+
 # set output directory
 exec_dir = os.getcwd()
 workdir: os.path.abspath(config['result_dir'])
@@ -588,10 +603,12 @@ rule run_breseq:
         "{sn}/benchmarks/{sn}_run_breseq.benchmark.tsv"
     params:
         ref = os.path.join(exec_dir, breseq_ref),
+        cov = breseq_cov,
+        freq = breseq_freq,
         outdir = '{sn}/breseq'
     shell:
         """
-        breseq --reference {params.ref} --num-processors {threads} --polymorphism-prediction --brief-html-output --output {params.outdir} {input} > {log} 2>&1 || touch {output}
+        breseq --reference {params.ref} --num-processors {threads} --polymorphism-prediction --polymorphism-minimum-variant-coverage-each-strand {params.cov} --polymorphism-frequency-cutoff {params.freq} --brief-html-output --output {params.outdir} {input} > {log} 2>&1 || touch {output}
         """
 
 ################## Based on https://github.com/jts/ncov2019-artic-nf/blob/be26baedcc6876a798a599071bb25e0973261861/modules/illumina.nf ##################
@@ -731,7 +748,6 @@ rule run_kraken2:
 
 ##################################  Based on scripts/quast.sh   ####################################
 
-
 rule run_quast:
     threads: 1
     conda: 'conda_envs/assembly_qc.yaml'

scripts/signal_postprocess.py

@@ -17,7 +17,7 @@
 
 assert long_git_id.startswith('$Id: ')
 #short_git_id = long_git_id[5:12]
-short_git_id = "v1.6.3"
+short_git_id = "v1.6.4"
 
 # Suppresses matplotlib warning (https://github.com/jaleezyy/covid-19-signal/issues/59)
 # Creates a small memory leak, but it's nontrivial to fix, and won't be a practical concern!
@@ -1343,18 +1343,41 @@ def __init__(self, name, ivarlin, fblin):
 
 
         if ivarlin['lineage'] != fblin['lineage'] and fblin['lineage'] is not None:
-            assert ivarlin['pangolin_ver'] == fblin['pangolin_ver']
-            assert ivarlin['pangodata_ver'] == fblin['pangodata_ver']
+           # pangolin version check (should match, but will report discrepency
+            try:
+                assert ivarlin['pangolin_ver'] == fblin['pangolin_ver']
+                pangolin_version = f"{ivarlin['pangolin_ver']}"
+            except AssertionError:
+                pangolin_version = f"{ivarlin['pangolin_ver']} (FB: {fblin['pangolin_ver']})"
+
+            # pangodata version (should match, but will report discrepency)
+            try:
+                assert ivarlin['pangodata_ver'] == fblin['pangodata_ver']
+                pangodata_version = f"{ivarlin['pangodata_ver']}"
+            except AssertionError:
+                pangodata_version = f"{ivarlin['pangodata_ver']} (FB: {fblin['pangodata_ver']})"
+
+            # report assigned clades between iVar and FreeBayes
             if ivarlin['clade'] == fblin['clade']:
-                 self.lineage = { 'lineage': str(ivarlin['lineage'] + " (FB: %s)" %(fblin['lineage'])),
-                                 'pangolin_ver': ivarlin['pangolin_ver'],
-                                 'pangodata_ver': ivarlin['pangodata_ver'],
-                                 'clade': ivarlin['clade'] }
+                assigned_clade = f"{ivarlin['clade']}"
             else:
-                 self.lineage = { 'lineage': str(ivarlin['lineage'] + " (FB: %s)" %(fblin['lineage'])),
-                                 'pangolin_ver': ivarlin['pangolin_ver'],
-                                 'pangodata_ver': ivarlin['pangodata_ver'],
-                                 'clade': str(ivarlin['clade'] + " (FB: %s)" %(fblin['clade'])) }
+                assigned_clade = f"{ivarlin['clade']} (FB: {fblin['clade']})"
+
+            self.lineage = { 'lineage': str(ivarlin['lineage'] + " (FB: %s)" %(fblin['lineage'])),
+                            'pangolin_ver': pangolin_version,
+                            'pangodata_ver': pangodata_version,
+                             'clade': assigned_clade }
+
+            # if ivarlin['clade'] == fblin['clade']:
+                 # self.lineage = { 'lineage': str(ivarlin['lineage'] + " (FB: %s)" %(fblin['lineage'])),
+                                 # 'pangolin_ver': ivarlin['pangolin_ver'],
+                                 # 'pangodata_ver': ivarlin['pangodata_ver'],
+                                 # 'clade': ivarlin['clade'] }
+            # else:
+                 # self.lineage = { 'lineage': str(ivarlin['lineage'] + " (FB: %s)" %(fblin['lineage'])),
+                                 # 'pangolin_ver': ivarlin['pangolin_ver'],
+                                 # 'pangodata_ver': ivarlin['pangodata_ver'],
+                                 # 'clade': str(ivarlin['clade'] + " (FB: %s)" %(fblin['clade'])) }
         else:
             self.lineage = ivarlin
 

signalexe.py

@@ -184,6 +184,10 @@ def write_config_file(run_name, config_file, data_directory, opt_tasks):
 run_breseq: {opt_tasks[0]}
 # Used as --reference argument to 'breseq'
 breseq_reference: "{data_directory}/MN908947.3.gbk"
+# Used as --polymorphism-minimum-variant-coverage-each-strand, --polymorphism-frequency-cutoff arguments
+# Parameters needed to determine thresholds for minor variant detection
+polymorphism_variant_coverage: 2
+polymorphism_frequency: 0.05
 
 # run freebayes for variant and consensus calling (as well as ivar)
 run_freebayes: {opt_tasks[1]}
@@ -274,7 +278,7 @@ def install_signal(frontend, data='data', unlock=False):
 	# note: add root_dir to determine the root directory of SIGNAL
 	script_path = os.path.join(os.path.abspath(sys.argv[0]).rsplit("/",1)[0])
 	args, allowed = create_parser()
-	version = 'v1.6.3'
+	version = 'v1.6.4'
 	alt_options = []
 
 	if args.version: