fixes implementation of transmuting_methyl_into_ring_score

openfe/setup/atom_mapping/lomap_scorers.py

@@ -291,10 +291,11 @@ def creates_heterocyle(mol):
 
 def transmuting_methyl_into_ring_score(mapping: LigandAtomMapping,
                                        beta=0.1, penalty=6.0):
-    """Penalises ring forming
+    """Penalises having a non-mapped ring atoms become a non-ring
 
-    Check if any atoms transition to/from rings in the mapping, if so
-    returns a score of::
+    This score would for example penalise R-CH3 to R-Ph where R is the same
+    mapped atom and both CH3 and Ph are unmapped. Does not penalise R-H to R-Ph.
+    If any atoms trigger the rule returns a score of::
 
       1 - exp(-1 * beta * penalty)
 
@@ -312,12 +313,15 @@ def transmuting_methyl_into_ring_score(mapping: LigandAtomMapping,
     molB = mapping.componentB.to_rdkit()
     molA_to_molB = mapping.componentA_to_componentB
 
+
     ringbreak = False
     for i, j in molA_to_molB.items():
         atomA = molA.GetAtomWithIdx(i)
 
         for bA in atomA.GetBonds():
             otherA = bA.GetOtherAtom(atomA)
+            if otherA.GetAtomicNum() == 1:
+                continue
             if otherA.GetIdx() in molA_to_molB:
                 # if other end of bond in core, ignore
                 continue
@@ -328,6 +332,8 @@ def transmuting_methyl_into_ring_score(mapping: LigandAtomMapping,
                 otherB = bB.GetOtherAtom(atomB)
                 if otherB.GetIdx() in molA_to_molB.values():
                     continue
+                if otherB.GetAtomicNum() == 1:
+                    continue
 
                 if otherA.IsInRing() ^ otherB.IsInRing():
                     ringbreak = True

openfe/tests/setup/atom_mapping/test_lomap_scorers.py

@@ -265,3 +265,34 @@ def test_lomap_regression(lomap_basic_test_files_dir,  # in a dir for lomap
         scores[i, i] = 0
 
     assert_allclose(matrix, scores)
+
+
+def test_transmuting_methyl_into_ring_score():
+    """
+    Sets up two mappings:
+      RC_to_RPh = [CCC]C -> [CCC]Ph
+      RH_to_RPh = [CCC]H -> [CCC]Ph
+    Where square brackets show mapped (core) region
+
+    The first mapping should trigger this rule, the second shouldn't
+    """
+    def makemol(smi):
+        m = Chem.MolFromSmiles(smi)
+        m = Chem.AddHs(m)
+        m.Compute2DCoords()
+
+        return openfe.SmallMoleculeComponent(m)
+
+    core = 'CCC{}'
+    RC = makemol(core.format('C'))
+    RPh = makemol(core.format('c1ccccc1'))
+    RH = makemol(core.format('[H]'))
+
+    RC_to_RPh = openfe.LigandAtomMapping(RC, RPh, {i: i for i in range(3)})
+    RH_to_RPh = openfe.LigandAtomMapping(RH, RPh, {i: i for i in range(3)})
+
+    score1 = lomap_scorers.transmuting_methyl_into_ring_score(RC_to_RPh)
+    score2 = lomap_scorers.transmuting_methyl_into_ring_score(RH_to_RPh)
+
+    assert score2 == 0.0
+    assert score1 == pytest.approx(1 - math.exp(-0.1 * 6.0))