Thanks Andrew, the SMILES approach seemed to have quite a few edge cases so
I wrote something to work directly on a molecule.
#!/usr/bin/env python
import sys
from rdkit import Chem
from collections import defaultdict
from rdkit.Chem.rdchem import EditableMol
# Thanks to steeveslab-blog for example of how to edit RDKit molecules
# http://asteeves.github.io/blog/2015/01/14/editing-in-rdkit/
# Thanks to Andrew Dalke for the function name
def weld_r_groups(input_mol):
# First pass loop over atoms and find the atoms with an AtomMapNum
join_dict = defaultdict(list)
for atom in input_mol.GetAtoms():
map_num = atom.GetAtomMapNum()
if map_num > 0:
join_dict[map_num].append(atom)
# Second pass, transfer the atom maps to the neighbor atoms
for idx, atom_list in join_dict.items():
if len(atom_list) == 2:
atm_1, atm_2 = atom_list
nbr_1 = [x.GetOtherAtom(atm_1) for x in atm_1.GetBonds()][0]
nbr_1.SetAtomMapNum(idx)
nbr_2 = [x.GetOtherAtom(atm_2) for x in atm_2.GetBonds()][0]
nbr_2.SetAtomMapNum(idx)
# Nuke all of the dummy atoms
new_mol = Chem.DeleteSubstructs(input_mol, Chem.MolFromSmarts('[#0]'))
# Third pass - arrange the atoms with AtomMapNum, these will be
connected
bond_join_dict = defaultdict(list)
for atom in new_mol.GetAtoms():
map_num = atom.GetAtomMapNum()
if map_num > 0:
bond_join_dict[map_num].append(atom.GetIdx())
# Make an editable molecule and add bonds between atoms with
correspoing AtomMapNum
em = EditableMol(new_mol)
for idx, atom_list in bond_join_dict.items():
if len(atom_list) == 2:
start_atm, end_atm = atom_list
em.AddBond(start_atm, end_atm,
order=Chem.rdchem.BondType.SINGLE)
final_mol = em.GetMol()
# remove the AtomMapNum values
for atom in final_mol.GetAtoms():
atom.SetAtomMapNum(0)
final_mol = Chem.RemoveHs(final_mol)
return final_mol
if __name__ == "__main__":
mol_to_weld = Chem.MolFromSmiles(
"CN(C)CC(Br)c1cc([*:2])c([*:1])cn1.[H]C([H])([H])[*:1].[H][*:2]")
welded_mol = weld_r_groups(mol_to_weld)
print(Chem.MolToSmiles(welded_mol))
Best,
Pat
On Sun, Apr 15, 2018 at 12:16 PM, Patrick Walters <wpwalt...@gmail.com>
wrote:
> Hi All,
>
> I was about to write a function to reassemble a molecule from a core +
> R-groups, but I thought I'd check and see if such a function already
> exists. This is work with the output of rdRGroupDecomposition
>
> Gvien a core:
> CN(C)CC(Br)c1cc([*:2])c([*:1])cn1
>
> Plus a set of R-groups
> [H]C([H])([H])[*:1]
> [H][*:2]
>
> Reconnect the pieces to generate a molecule
> CN(C)CC(Br)c1ccc(C)cn1
>
> Thanks,
>
> Pat
>
>
>
>
------------------------------------------------------------------------------
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot
_______________________________________________
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
