Package rdkit :: Package Chem :: Package Scaffolds :: Module MurckoScaffold
[hide private]
[frames] | no frames]

Source Code for Module rdkit.Chem.Scaffolds.MurckoScaffold

  1  # $Id: MurckoScaffold.py 3672 2010-06-14 17:10:00Z landrgr1 $ 
  2  # 
  3  # Created by Peter Gedeck, September 2008 
  4  # 
  5  """ 
  6    Generation of Murcko scaffolds from a molecule 
  7  """ 
  8   
  9  from rdkit import Chem 
 10  from rdkit.Chem import AllChem 
 11   
 12  murckoTransforms=[AllChem.ReactionFromSmarts('[*:1]-[!#1;D1]>>[*:1][H]'), 
 13                    AllChem.ReactionFromSmarts('[*:1]-[!#1;D2]#[AD1]>>[*:1][H]'), 
 14                    AllChem.ReactionFromSmarts('[*:1]-[!#1;D2]=[AD1]>>[*:1][H]'), 
 15                    AllChem.ReactionFromSmarts('[*:1]-[!#1;D3](=[AD1])=[AD1]>>[*:1][H]')] 
 16   
 17   
18 -def MakeScaffoldGeneric(mol):
19 """ Makes a Murcko scaffold generic (i.e. all atom types->C and all bonds ->single 20 21 >>> Chem.MolToSmiles(MakeScaffoldGeneric(Chem.MolFromSmiles('c1ccccc1'))) 22 'C1CCCCC1' 23 >>> Chem.MolToSmiles(MakeScaffoldGeneric(Chem.MolFromSmiles('c1ncccc1'))) 24 'C1CCCCC1' 25 26 The following were associated with sf.net issue 246 27 >>> Chem.MolToSmiles(MakeScaffoldGeneric(Chem.MolFromSmiles('c1[nH]ccc1'))) 28 'C1CCCC1' 29 >>> Chem.MolToSmiles(MakeScaffoldGeneric(Chem.MolFromSmiles('C1[NH2+]C1'))) 30 'C1CC1' 31 >>> Chem.MolToSmiles(MakeScaffoldGeneric(Chem.MolFromSmiles('C1[C@](Cl)(F)O1'))) 32 'CC1(C)CC1' 33 34 """ 35 res = Chem.Mol(mol) 36 for atom in res.GetAtoms(): 37 if atom.GetAtomicNum()!=1: 38 atom.SetAtomicNum(6) 39 atom.SetIsAromatic(False) 40 atom.SetFormalCharge(0) 41 atom.SetChiralTag(Chem.ChiralType.CHI_UNSPECIFIED) 42 atom.SetNoImplicit(0) 43 atom.SetNumExplicitHs(0) 44 for bond in res.GetBonds(): 45 bond.SetBondType(Chem.BondType.SINGLE) 46 bond.SetIsAromatic(False) 47 return Chem.RemoveHs(res)
48 49 50 murckoPatts = ['[!#1;D3;$([D3]-[!#1])](=[AD1])=[AD1]', 51 '[!#1;D2;$([D2]-[!#1])]=,#[AD1]', 52 '[!#1;D1;$([D1]-[!#1;!n])]'] 53 murckoQ = '['+','.join(['$(%s)'%x for x in murckoPatts])+']' 54 murckoQ = Chem.MolFromSmarts(murckoQ) 55 murckoPatts = [Chem.MolFromSmarts(x) for x in murckoPatts] 56 aromaticNTransform = AllChem.ReactionFromSmarts('[n:1]-[D1]>>[nH:1]')
57 -def GetScaffoldForMol(mol):
58 """ Return molecule object containing scaffold of mol 59 60 >>> m = Chem.MolFromSmiles('Cc1ccccc1') 61 >>> GetScaffoldForMol(m) 62 <rdkit.Chem.rdchem.Mol object at 0x...> 63 >>> Chem.MolToSmiles(GetScaffoldForMol(m)) 64 'c1ccccc1' 65 66 >>> m = Chem.MolFromSmiles('Cc1cc(Oc2nccc(CCC)c2)ccc1') 67 >>> Chem.MolToSmiles(GetScaffoldForMol(m)) 68 'c1ccc(Oc2ccccn2)cc1' 69 70 """ 71 if 1: 72 res=Chem.MurckoDecompose(mol) 73 res.ClearComputedProps() 74 res.UpdatePropertyCache() 75 Chem.GetSymmSSSR(res) 76 else: 77 res=_pyGetScaffoldForMol(mol) 78 return res
79
80 -def _pyGetScaffoldForMol(mol):
81 while mol.HasSubstructMatch(murckoQ): 82 for patt in murckoPatts: 83 mol = Chem.DeleteSubstructs(mol,patt) 84 for atom in mol.GetAtoms(): 85 if atom.GetAtomicNum()==6 and atom.GetNoImplicit() and atom.GetExplicitValence()<4: 86 atom.SetNoImplicit(False) 87 h = Chem.MolFromSmiles('[H]') 88 mol = Chem.ReplaceSubstructs(mol,Chem.MolFromSmarts('[D1;$([D1]-n)]'),h,True)[0]; 89 mol=Chem.RemoveHs(mol) 90 #while 1: 91 # ps = aromaticNTransform.RunReactants([mol]) 92 # if ps: 93 # mol = ps[0][0] 94 # else: 95 # break 96 return mol 97
98 -def MurckoScaffoldSmilesFromSmiles(smiles,includeChirality=False):
99 """ Returns MurckScaffold Smiles from smiles 100 101 >>> MurckoScaffoldSmilesFromSmiles('Cc1cc(Oc2nccc(CCC)c2)ccc1') 102 'c1ccc(Oc2ccccn2)cc1' 103 104 """ 105 mol = Chem.MolFromSmiles(smiles) 106 scaffold = GetScaffoldForMol(mol) 107 if not scaffold: return None 108 return Chem.MolToSmiles(scaffold,includeChirality)
109
110 -def MurckoScaffoldSmiles(smiles=None, mol=None, includeChirality=False):
111 """ Returns MurckScaffold Smiles from smiles 112 113 >>> MurckoScaffoldSmiles('Cc1cc(Oc2nccc(CCC)c2)ccc1') 114 'c1ccc(Oc2ccccn2)cc1' 115 116 >>> MurckoScaffoldSmiles(mol=Chem.MolFromSmiles('Cc1cc(Oc2nccc(CCC)c2)ccc1')) 117 'c1ccc(Oc2ccccn2)cc1' 118 119 """ 120 if smiles: 121 mol = Chem.MolFromSmiles(smiles) 122 else: 123 mol = mol 124 if mol is None: 125 raise ValueError('No molecule provided') 126 scaffold = GetScaffoldForMol(mol) 127 if not scaffold: return None 128 return Chem.MolToSmiles(scaffold,includeChirality)
129 130 #------------------------------------ 131 # 132 # doctest boilerplate 133 #
134 -def _test():
135 import doctest,sys 136 return doctest.testmod(sys.modules["__main__"],optionflags=doctest.ELLIPSIS)
137 138 if __name__ == '__main__': 139 import sys 140 failed,tried = _test() 141 sys.exit(failed) 142