Fix: Only write auspice.json when there is sequence information

.gitignore

@@ -96,4 +96,7 @@ Icon?
 ehthumbs.db
 Thumbs.db
 
-test/treetime_examples
+test/treetime_examples
+
+# Output directories produced by bash test/command_line_tests.sh
+test/20*/

changelog.md

@@ -7,6 +7,7 @@ Information on outliers is saved in a pandas DataFrame `self.outliers` of `TreeT
 
 ### Other fixes
  * error when rate estimate is negative during the rate susceptibility calculation. Give hint in error message to specify the rate and its uncertainty explicitly.
+ * Fix bug [issue #250](https://github.com/neherlab/treetime/issues/250) introduced in 0.10.0 where treetime fails in absence of an alignment when trying to create an auspice json file. [PR #251](https://github.com/neherlab/treetime/pull/251)
 
 # 0.10.1: bug fix release
 

test.sh

@@ -3,10 +3,21 @@
 set -euo pipefail
 
 cd test
+
+# Remove treetime_examples in case it exists to not fail
+rm -rf treetime_examples
 git clone https://github.com/neherlab/treetime_examples.git
+
 bash command_line_tests.sh
 OUT=$?
 if [ "$OUT" != 0 ]; then
   exit 1
 fi
+
 pytest test_treetime.py
+if [ "$OUT" != 0 ]; then
+  exit 1
+fi
+
+# Clean up, the 202* is to remove auto-generated output dirs
+rm -rf treetime_examples __pycache__ 202*

test/command_line_tests.sh

@@ -63,6 +63,16 @@ else
 	echo "skyline approximation failed $retval"
 fi
 
+# From https://github.com/neherlab/treetime/issues/250 
+treetime --tree treetime_examples/data/ebola/ebola.nwk --dates treetime_examples/data/ebola/ebola.metadata.csv --sequence-length 1000
+retval="$?"
+if [ "$retval" == 0 ]; then
+	echo "sequence length only ok"
+else
+	((all_tests++))
+	echo "sequence length only failed $retval"
+fi
+
 if [ "$all_tests" == 0 ];then
 	echo "All tests passed"
 	exit 0

treetime/CLI_io.py

@@ -182,7 +182,8 @@ def export_sequences_and_tree(tt, basename, is_vcf=False, zero_based=False,
         Phylo.write(tt.tree, outtree_name, 'nexus', format_branch_length=fmt_bl)
     print("--- tree saved in nexus format as  \n\t %s\n"%outtree_name)
 
-    auspice = create_auspice_json(tt, timetree=timetree, confidence=confidence)
+    # Only create auspice json if there is sequence information
+    auspice = create_auspice_json(tt, timetree=timetree, confidence=confidence, seq_info=seq_info)
     outtree_name_json = basename + f'auspice_tree{tree_suffix}.json'
     with open(outtree_name_json, 'w') as fh:
         import json
@@ -229,7 +230,7 @@ def print_save_plot_skyline(tt, n_std=2.0, screen=True, save='', plot='', gen=50
 
 
 
-def create_auspice_json(tt, timetree=False, confidence=False):
+def create_auspice_json(tt, timetree=False, confidence=False, seq_info=False):
     # mock up meta data for auspice json
     from datetime import datetime
     meta = {
@@ -280,11 +281,12 @@ def node_to_json(n, pdiv=0.0):
         j["node_attrs"]["div"] = float(pdiv + n.mutation_length)
         j["node_attrs"]["bad_branch"] = {"value": "Yes" if n.bad_branch else "No"}
 
-        j["branch_attrs"]["mutations"] = {"nuc": [f"{a}{pos+1}{d}" for a,pos,d in n.mutations if d in "ACGT-"]}
-        # generate bootstrap confidence substitute via the negative exponential of the number of mutations
-        # this is the bootstrap confidence for iid mutations (only ACGT mutations)
-        j["node_attrs"]["confidence"] = {"value":round(1-np.exp(-len([pos for a,pos,d in n.mutations if d in "ACGT"])),3)
-                                          if not n.is_terminal() else 1.0}
+        if seq_info: # only add mutations to the json if run with sequence data (fasta or vcf)
+            j["branch_attrs"]["mutations"] = {"nuc": [f"{a}{pos+1}{d}" for a,pos,d in n.mutations if d in "ACGT-"]}
+            # generate bootstrap confidence substitute via the negative exponential of the number of mutations
+            # this is the bootstrap confidence for iid mutations (only ACGT mutations)
+            j["node_attrs"]["confidence"] = {"value":round(1-np.exp(-len([pos for a,pos,d in n.mutations if d in "ACGT"])),3)
+                                            if not n.is_terminal() else 1.0}
         return j
 
     # create the tree data structure from the Biopython tree