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