Package rdkit :: Package Chem :: Package AtomPairs :: Module Torsions
[hide private]
[frames] | no frames]

Source Code for Module rdkit.Chem.AtomPairs.Torsions

  1  # $Id$ 
  2  # 
  3  #  Copyright (C) 2004-2006 Rational Discovery LLC 
  4  # 
  5  #   @@ All Rights Reserved @@ 
  6  #  This file is part of the RDKit. 
  7  #  The contents are covered by the terms of the BSD license 
  8  #  which is included in the file license.txt, found at the root 
  9  #  of the RDKit source tree. 
 10  # 
 11  """ 
 12  Contains an implementation of Topological-torsion fingerprints, as 
 13  described in: 
 14   
 15  R. Nilakantan, N. Bauman, J. S. Dixon, R. Venkataraghavan; 
 16  "Topological Torsion: A New Molecular Descriptor for SAR Applications. 
 17  Comparison with Other Descriptors" JCICS 27, 82-85 (1987). 
 18   
 19  The fingerprints can be accessed through the following functions: 
 20  - GetTopologicalTorsionFingerprint 
 21  - GetHashedTopologicalTorsionFingerprint 
 22  - GetTopologicalTorsionFingerprintAsIntVect (identical to GetTopologicalTorsionFingerprint) 
 23  - GetTopologicalTorsionFingerprintAsIds 
 24   
 25  """ 
 26  from rdkit.Chem import rdMolDescriptors 
 27  from rdkit.Chem.AtomPairs import Utils 
 28  from rdkit.Chem.rdMolDescriptors import (GetTopologicalTorsionFingerprint, 
 29                                           GetHashedTopologicalTorsionFingerprint) 
 30  GetTopologicalTorsionFingerprintAsIntVect = rdMolDescriptors.GetTopologicalTorsionFingerprint 
 31   
 32   
33 -def pyScorePath(mol, path, size, atomCodes=None):
34 """ Returns a score for an individual path. 35 36 >>> from rdkit import Chem 37 >>> m = Chem.MolFromSmiles('CCCCC') 38 >>> c1 = Utils.GetAtomCode(m.GetAtomWithIdx(0),1) 39 >>> c2 = Utils.GetAtomCode(m.GetAtomWithIdx(1),2) 40 >>> c3 = Utils.GetAtomCode(m.GetAtomWithIdx(2),2) 41 >>> c4 = Utils.GetAtomCode(m.GetAtomWithIdx(3),1) 42 >>> t = c1 | (c2 << rdMolDescriptors.AtomPairsParameters.codeSize) | (c3 << (rdMolDescriptors.AtomPairsParameters.codeSize*2)) | (c4 << (rdMolDescriptors.AtomPairsParameters.codeSize*3)) 43 >>> pyScorePath(m,(0,1,2,3),4)==t 44 1 45 46 The scores are path direction independent: 47 >>> pyScorePath(m,(3,2,1,0),4)==t 48 1 49 50 >>> m = Chem.MolFromSmiles('C=CC(=O)O') 51 >>> c1 = Utils.GetAtomCode(m.GetAtomWithIdx(0),1) 52 >>> c2 = Utils.GetAtomCode(m.GetAtomWithIdx(1),2) 53 >>> c3 = Utils.GetAtomCode(m.GetAtomWithIdx(2),2) 54 >>> c4 = Utils.GetAtomCode(m.GetAtomWithIdx(4),1) 55 >>> t = c1 | (c2 << rdMolDescriptors.AtomPairsParameters.codeSize) | (c3 << (rdMolDescriptors.AtomPairsParameters.codeSize*2)) | (c4 << (rdMolDescriptors.AtomPairsParameters.codeSize*3)) 56 >>> pyScorePath(m,(0,1,2,4),4)==t 57 1 58 59 """ 60 codes = [None] * size 61 for i in range(size): 62 if i == 0 or i == (size - 1): 63 sub = 1 64 else: 65 sub = 2 66 if not atomCodes: 67 codes[i] = Utils.GetAtomCode(mol.GetAtomWithIdx(path[i]), sub) 68 else: 69 base = atomCodes[path[i]] 70 codes[i] = base - sub 71 72 # "canonicalize" the code vector: 73 beg = 0 74 end = len(codes) - 1 75 while (beg < end): 76 if codes[beg] > codes[end]: 77 codes.reverse() 78 break 79 elif codes[beg] == codes[end]: 80 beg += 1 81 end -= 1 82 else: 83 break 84 accum = 0 85 for i in range(size): 86 accum |= codes[i] << (rdMolDescriptors.AtomPairsParameters.codeSize * i) 87 return accum
88 89
90 -def ExplainPathScore(score, size=4):
91 """ 92 93 >>> from rdkit import Chem 94 >>> m = Chem.MolFromSmiles('C=CC') 95 >>> score=pyScorePath(m,(0,1,2),3) 96 >>> ExplainPathScore(score,3) 97 (('C', 1, 0), ('C', 2, 1), ('C', 1, 1)) 98 99 Again, it's order independent: 100 >>> score=pyScorePath(m,(2,1,0),3) 101 >>> ExplainPathScore(score,3) 102 (('C', 1, 0), ('C', 2, 1), ('C', 1, 1)) 103 104 >>> m = Chem.MolFromSmiles('C=CO') 105 >>> score=pyScorePath(m,(0,1,2),3) 106 >>> ExplainPathScore(score,3) 107 (('C', 1, 1), ('C', 2, 1), ('O', 1, 0)) 108 109 >>> m = Chem.MolFromSmiles('OC=CO') 110 >>> score=pyScorePath(m,(0,1,2,3),4) 111 >>> ExplainPathScore(score,4) 112 (('O', 1, 0), ('C', 2, 1), ('C', 2, 1), ('O', 1, 0)) 113 114 >>> m = Chem.MolFromSmiles('CC=CO') 115 >>> score=pyScorePath(m,(0,1,2,3),4) 116 >>> ExplainPathScore(score,4) 117 (('C', 1, 0), ('C', 2, 1), ('C', 2, 1), ('O', 1, 0)) 118 119 120 >>> m = Chem.MolFromSmiles('C=CC(=O)O') 121 >>> score=pyScorePath(m,(0,1,2,3),4) 122 >>> ExplainPathScore(score,4) 123 (('C', 1, 1), ('C', 2, 1), ('C', 3, 1), ('O', 1, 1)) 124 >>> score=pyScorePath(m,(0,1,2,4),4) 125 >>> ExplainPathScore(score,4) 126 (('C', 1, 1), ('C', 2, 1), ('C', 3, 1), ('O', 1, 0)) 127 128 129 >>> m = Chem.MolFromSmiles('OOOO') 130 >>> score=pyScorePath(m,(0,1,2),3) 131 >>> ExplainPathScore(score,3) 132 (('O', 1, 0), ('O', 2, 0), ('O', 2, 0)) 133 >>> score=pyScorePath(m,(0,1,2,3),4) 134 >>> ExplainPathScore(score,4) 135 (('O', 1, 0), ('O', 2, 0), ('O', 2, 0), ('O', 1, 0)) 136 137 """ 138 codeMask = (1 << rdMolDescriptors.AtomPairsParameters.codeSize) - 1 139 res = [None] * size 140 for i in range(size): 141 if i == 0 or i == (size - 1): 142 sub = 1 143 else: 144 sub = 2 145 code = score & codeMask 146 score = score >> rdMolDescriptors.AtomPairsParameters.codeSize 147 symb, nBranch, nPi = Utils.ExplainAtomCode(code) 148 expl = symb, nBranch + sub, nPi 149 res[i] = expl 150 return tuple(res)
151 152
153 -def GetTopologicalTorsionFingerprintAsIds(mol, targetSize=4):
154 iv = GetTopologicalTorsionFingerprint(mol, targetSize) 155 res = [] 156 for k, v in iv.GetNonzeroElements().items(): 157 res.extend([k] * v) 158 res.sort() 159 return res
160 161 162 # ------------------------------------ 163 # 164 # doctest boilerplate 165 #
166 -def _runDoctests(verbose=None): # pragma: nocover
167 import sys 168 import doctest 169 failed, _ = doctest.testmod(optionflags=doctest.ELLIPSIS, verbose=verbose) 170 sys.exit(failed) 171 172 173 if __name__ == '__main__': # pragma: nocover 174 _runDoctests() 175