In the course of our work towards an automated Structure Activity Relation (SAR) analysis system that’s intuitive for medicinal chemists to use, we needed to depict R-Groups on our web interface in an aligned manner. We use #OpenEyeScientific chemistry toolkits for our production servers as the speed, quality, and reliability is high. The challenge is the alignment is to just two atoms ([R1]-*), one of which is a wildcard. So how do you do this with @OpenEyeSoftware?
We include a simple Python script (v3.7) to illustrate how this is done – this is a modification from OpenEye’s docs. In production we depict each R-group separately and store them in a repository; from our matched-pair analysis, we know the most frequent observed R-groups, and so we can render these and stored them in advance for retrieval at speed. Here we have produced a 3 by 3 grid of the depictions, with default bond width and scaling.
In the OpenEye world a zero atomic number atom can be mapped to denote an R-group [R1],[R2]…[Rn], and so are the equivalent of [*:1], [*:2]..[*:n]. The alignment SMARTS [R1]-*, is parsed into a graph, and then a custom function prepares co-ordinates for the depiction. The custom function NonCanonicalPrepareDepiction, performs the critical operations to prepare the graph for depiction, but as the name suggests, does not canonicalize the SMARTS. If this step is included the depictions will be orientated with [R1] on the right hand side, not the left. Swap out this function for OEPrepareDepiction, and you will see this occur: this might actually be desirable and could be a feature?
A list of OEMolGraph of the R-groups is generated ready for iteration over and placement into the Image Grid. After some scaling work we are ready to align the structures. The key OE method (OEPrepareAlignedDepiction) takes the Rgroup graph, the align Rgroup graph and adjusts the 2D co-ordinates to this reference. An OEMatch object is required to inform which atoms and bonds are the points of alignment. In order for this to work, the custom function _find_Rgroup_match takes the two graphs and finds [R1] in both molecules, captures that pair of atoms, the connecting bond pair and matching neighbour atoms, which, of course, do not have actually be the same atom types. These pairs are added to an OEMatch object – very nice API we must say. The display mol is prepared, and the image rendered into the grid.
We would like to thank OE support, and Luke in particular, for solving the [R1] group left / right problem.
# MedChemica Limited
# Janurary 2022
from openeye.oechem import (
OECreateSmiString,
OEGraphMol,
OEIsRGroup,
OEMatch,
OESmilesToMol,
OESetDimensionFromCoords,
OEPerceiveChiral,
OE3DToAtomStereo,
OE3DToBondStereo,
OESuppressHydrogens,
OEMDLPerceiveBondStereo,
)
from openeye.oedepict import (
OE2DMolDisplay,
OE2DMolDisplayOptions,
OEImage,
OEImageGrid,
OEGetMoleculeScale,
OEPrepareDepiction,
OEPrepareAlignedDepiction,
OERenderMolecule,
OEScale_AutoScale,
OEWriteImage,
OEAddDepictionHydrogens,
OEDepictCoordinates)
def NonCanonicalPrepareDepiction(mol, clearcoords=False, suppressH=True):
OESetDimensionFromCoords(mol)
OEPerceiveChiral(mol)
if mol.GetDimension() != 2 or clearcoords:
if mol.GetDimension() == 3:
OE3DToBondStereo(mol)
OE3DToAtomStereo(mol)
if suppressH:
OESuppressHydrogens(mol)
OEAddDepictionHydrogens(mol)
OEDepictCoordinates(mol)
OEMDLPerceiveBondStereo(mol)
mol.SetDimension(2)
return True
def _find_R1bond_and_neighbour(mol):
'''
Specifically find [R1] atom, its bond and neighbour atom
'''
_align_atoms = []
_align_bonds = []
for _atom in mol.GetAtoms(OEIsRGroup()):
if _atom.GetMapIdx() == 1:
_align_atoms.append(_atom)
for _bond in _atom.GetBonds():
_align_bonds.append(_bond)
_align_atoms.append(_bond.GetNbr(_atom))
return _align_atoms, _align_bonds
def _find_Rgroup_match(align_mol, match_mol):
_align_atoms, _align_bonds = _find_R1bond_and_neighbour(align_mol)
_match_atoms, _match_bonds = _find_R1bond_and_neighbour(match_mol)
match = OEMatch()
for p1, p2 in zip(_align_atoms, _match_atoms):
match.AddPair(p1, p2)
for p1, p2 in zip(_align_bonds, _match_bonds):
match.AddPair(p1, p2)
return match
# Align Rgroups to what SMILES
MCPAIRS_OE_RGROUP_ALIGN_TO = '[R1]-*'
rgroups_to_depict = [
'[R1]CCCCC R1-n-propyl',
'[R1]C(C)C R1-iso-propyl',
'[R1]OC(C)C R1-O-isopropyl',
'[R1]c1ccccc1[R2] R1-benzene-R2',
'[R1]NC(C)C R1-N-isopropyl',
'[R1]NC(=O)[R2] R1-amide-R2',
'[R1]N(C)C(=O)[R2] R1-N(CH3)amide',
'[R1]C1CCC([R2])CC1 R1-cyclohexyl-R2',
'[R1]Cc1ccccc1 R1-benzyl',
]
_align_rgroup_mol = OEGraphMol()
OESmilesToMol(_align_rgroup_mol, MCPAIRS_OE_RGROUP_ALIGN_TO)
NonCanonicalPrepareDepiction(_align_rgroup_mol)
_rgroup_mols = []
for smi in rgroups_to_depict:
mol = OEGraphMol()
OESmilesToMol(mol, smi)
OEPrepareDepiction(mol)
_rgroup_mols.append(OEGraphMol(mol))
image = OEImage(6400, 3600)
rows, cols = 3, 3
grid = OEImageGrid(image, rows, cols)
opts2 = OE2DMolDisplayOptions(
grid.GetCellWidth(), grid.GetCellHeight(), OEScale_AutoScale)
minscale = float("inf")
for _rgroup_mol in _rgroup_mols:
minscale = min(minscale, OEGetMoleculeScale(_rgroup_mol, opts2))
opts2.SetScale(minscale)
for _rgroup_mol, cell in zip(_rgroup_mols, grid.GetCells()):
_matches = _find_Rgroup_match(_align_rgroup_mol, _rgroup_mol)
alignres = OEPrepareAlignedDepiction(_rgroup_mol, _align_rgroup_mol, _matches)
disp = OE2DMolDisplay(_rgroup_mol, opts2)
if not alignres:
print(f'depiction rgroup alignment failed : '
f'Rgroup {OECreateSmiString(_rgroup_mol)} : align_to {MCPAIRS_OE_RGROUP_ALIGN_TO}')
OERenderMolecule(cell, disp)
OEWriteImage("RgroupsAlignedMatrixExample.png", image)
@OpenEyeSoftware
#OpenEyeScientific
Footnote: Can this be done with RDkit? Yes we have a method for that too…..