...

epdb/models.py

@@ -1573,11 +1573,9 @@ class Pathway(EnviPathModel, AliasMixin, ScenarioMixin):
         while len(queue):
             current = queue.pop()
             processed.add(current)
-
             nodes.append(current.d3_json())
 
-            for e in self.edges:
-                if current in e.start_nodes.all():
+            for e in self.edges.filter(start_nodes=current).distinct():
                 for prod in e.end_nodes.all():
                     if prod not in queue and prod not in processed:
                         queue.append(prod)

epdb/tasks.py

@@ -172,6 +172,7 @@ def predict(
 
     except Exception as e:
         pw.kv.update({"status": "failed"})
+        pw.kv.update(**{"error": str(e)})
         pw.save()
 
         if JobLog.objects.filter(task_id=self.request.id).exists():

utilities/chem.py

@@ -2,12 +2,13 @@ import logging
 import re
 from abc import ABC
 from collections import defaultdict
-from typing import List, Optional, Dict, TYPE_CHECKING
+from typing import List, Optional, Dict, TYPE_CHECKING, Union
 
 from indigo import Indigo, IndigoException, IndigoObject
 from indigo.renderer import IndigoRenderer
 from rdkit import Chem, rdBase
 from rdkit.Chem import MACCSkeys, Descriptors
+from rdkit.Chem import rdchem
 from rdkit.Chem import rdChemReactions
 from rdkit.Chem.Draw import rdMolDraw2D
 from rdkit.Chem.MolStandardize import rdMolStandardize
@@ -90,8 +91,15 @@ class FormatConverter(object):
         return Chem.MolToSmiles(mol, canonical=canonical)
 
     @staticmethod
-    def InChIKey(smiles):
-        return Chem.MolToInchiKey(FormatConverter.from_smiles(smiles))
+    def InChIKey(mol_or_smiles: Union[rdchem.Mol | str]):
+        if isinstance(mol_or_smiles, str):
+            mol_or_smiles = mol_or_smiles.replace("~", "")
+            mol_or_smiles = FormatConverter.from_smiles(mol_or_smiles)
+
+        if mol_or_smiles is None:
+            return None
+
+        return Chem.MolToInchiKey(mol_or_smiles)
 
     @staticmethod
     def InChI(smiles):
@@ -288,7 +296,8 @@ class FormatConverter(object):
                             product = GetMolFrags(product, asMols=True)
                             for p in product:
                                 p = FormatConverter.standardize(
-                                    Chem.MolToSmiles(p), remove_stereo=remove_stereo
+                                    Chem.MolToSmiles(p).replace("~", ""),
+                                    remove_stereo=remove_stereo,
                                 )
                                 prods.append(p)
 

utilities/misc.py

@@ -35,6 +35,7 @@ from epdb.models import (
     RuleBasedRelativeReasoning,
     Scenario,
     SequentialRule,
+    Setting,
     SimpleAmbitRule,
     SimpleRDKitRule,
     SimpleRule,
@@ -1258,3 +1259,46 @@ class PathwayUtils:
                 res[edge.url] = rule_chain
 
         return res
+
+    def _find_intermediates(self, data_pathway, pred_pathway):
+        pass
+
+    def engineer(self, setting: "Setting"):
+        from epdb.logic import SPathway
+
+        # get a fresh copy
+        pw = Pathway.objects.get(id=self.pathway.pk)
+
+        root_nodes = [n.default_node_label.smiles for n in pw.root_nodes]
+
+        if len(root_nodes) != 1:
+            logger.warning(f"Pathway {pw.name} has {len(root_nodes)} root nodes")
+            return
+
+        spw = SPathway(root_nodes[0], None, setting)
+
+        level = 0
+        while not spw.done:
+            spw.predict_step(from_depth=level)
+            level += 1
+
+        # Generate Node / SMILES mapping
+        node_mapping = {}
+        from utilities.chem import FormatConverter
+
+        for node in pw.nodes:
+            for snode in spw.smiles_to_node.values():
+                data_smiles = node.default_node_label.smiles
+                pred_smiles = snode.smiles
+
+                data_key = FormatConverter.InChIKey(data_smiles.replace("~", ""))
+                pred_key = FormatConverter.InChIKey(pred_smiles.replace("~", ""))
+
+                if data_key == pred_key:
+                    node_mapping[snode] = node
+
+        print(node_mapping)
+
+        return spw
+
+        pass