-
Notifications
You must be signed in to change notification settings - Fork 1
/
recap.sql.in
155 lines (139 loc) · 5.64 KB
/
recap.sql.in
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
Drop Table If Exists _OCSCHEMA_.recap;
Create Table _OCSCHEMA_.recap (id integer, description text, smarts text);
-- each smarts must contain only two atoms representing the bond to be broken.
-- of course, each atom may be a complex atom smarts, ala [$(whatever)]
Copy _OCSCHEMA_.recap (id, description, smarts) From STDIN;
1 amide [$([C;!$(C([#7])[#7])](=!@[O]))]!@[$([#7;+0;!D1])]
2 ester [$(C=!@O)]!@[$([O;+0])]
3 amine [$([N;!D1;+0;!$(N-C=[#7,#8,#15,#16])](-!@[*]))]-!@[$([*])]
4 urea [$(C(=!@O)([#7;+0;D2,D3])!@[#7;+0;D2,D3])]!@[$([#7;+0;D2,D3])]
5 ether [$([O;+0](-!@[#6!$(C=O)])-!@[#6!$(C=O)])]-!@[$([#6!$(C=O)])]
6 olefin C=!@C
7 quaternaryN [N;+1;D4]!@[#6]
8 aromaticN-aliphaticC [$([n;+0])]-!@C
9 lactamN-aromaticC [$([O]=[C]-@[N;+0])]-!@[$([C])]
10 aromaticC-aromaticC c-!@c
11 sulphonamide [$([#7;+0;D2,D3])]-!@[$([S](=[O])=[O])]
\.
Grant Select On _OCSCHEMA_.recap To Public;
Create Or Replace Function _OCSCHEMA_.recap_set(smiles Text, label boolean) Returns SetOf Text As $EOPY$
import openbabel
class recap:
# RECAP-Retrosynthetic Combinatorial Analysis Procedure
# J. Chem. Inf. Comput. Sci. 1998, 38, 511-522
def __init__(self, mol, minsize=5):
self.mol = mol;
# minimum allowed size (atom count) of fragment
self.minsize = minsize;
# bonded atom pairs populated by the apply method,
# subsequently used by split and add_star
self.atom_pairs = list()
def apply(self, pat, patnum):
if pat.Match(self.mol):
# find all atom pairs that match
for p in pat.GetUMapList():
i = 0
atoms = list()
for a in openbabel.OBMolAtomIter(self.mol):
i += 1
if i in p:
atoms.append(a)
if self.small_fragment(atoms[0], atoms[1]):
#print True
pass
else:
atoms.append(patnum)
self.atom_pairs.append(atoms)
#print False
return True
else:
return False
def split(self, label):
for a in self.atom_pairs:
if label:
a[0].SetIsotope(a[2])
a[1].SetIsotope(a[2])
bond = a[0].GetBond(a[1])
# bond could be null if already deleted when smarts matched multiple times
if bond: self.mol.DeleteBond(bond)
def add_star(self):
for pair in self.atom_pairs:
self.mol.AddBond(pair[0].GetIdx(),self.mol.NewAtom().GetIdx(),1)
self.mol.AddBond(pair[1].GetIdx(),self.mol.NewAtom().GetIdx(),1)
def decide_multiples(self):
# some smarts (e.g. ether, amine) allow multiple bonds to the
# central atom to be broken. Yet it appears the central atom
# needs to be retained in one of the multiple fragments.
# If multiple fragments, let it stay with the smallest fragment.
# If tied, pick the first fragment.
multiples = [0]*self.mol.NumAtoms()
for pair in self.atom_pairs:
multiples[pair[0].GetIdx()] += 1
multiples[pair[1].GetIdx()] += 1
#print multiples
currsize = -1
currpair = None
for pair in self.atom_pairs:
a = pair[0]
b = pair[1]
if multiples[a.GetIdx()] > 1 or multiples[b.GetIdx()] > 1:
# remove larger fragment(s) if a-b were broken
#print a.GetIdx(),b.GetIdx(),
fsize = self.fragment_size(a,b)
if currpair == None:
currpair = pair
currsize = fsize
else:
if fsize < currsize:
self.atom_pairs.remove(pair)
else:
self.atom_pairs.remove(currpair)
currpair = pair
currsize = fsize
def fragment_size(self, a, b):
# size of fragment b if a-b were broken
c1 = openbabel.vectorInt()
self.mol.FindChildren(c1,a.GetIdx(),b.GetIdx())
#for atom in c1:
# if self.mol.GetAtom(atom).GetValence() == 1:
return 1+len(c1)
def small_fragment(self, a, b):
# if we were to break the bond between a and b,
# would either fragment be too small?
#print a.GetIdx(), b.GetIdx(),
if self.fragment_size(a,b) < self.minsize: return True
if self.fragment_size(b,a) < self.minsize: return True
return False
import obchord
oc=obchord.obchord(GD)
# make copy in order to alter mol
amol = openbabel.OBMol(oc.parse_smi(smiles))
amol.SetTitle("")
Recap = recap(amol,4)
rx = plpy.execute("Select * from _OCSCHEMA_.recap")
for row in rx:
pat = oc.parse_sma(row["smarts"])
Recap.apply(pat, row["id"])
Recap.decide_multiples()
Recap.split(label)
#Recap.add_star()
#return [oc.writestring(amol,"can")]
return oc.writestring(amol,"can").split('.')
$EOPY$ Language plpythonu Immutable Returns Null On Null Input;
Comment On FUNCTION _OCSCHEMA_.recap_set(text, boolean)
Is 'return recap fragments for input smiles with optional isotope labels';
Create Or Replace Function _OCSCHEMA_.recap_set(smiles Text) Returns SetOf Text As $EOSQL$
Select _OCSCHEMA_.recap_set($1, false);
$EOSQL$ Language SQL Immutable Returns Null On Null Input;
Comment On FUNCTION _OCSCHEMA_.recap_set(text)
Is 'return recap fragments for input smiles without isotope labels';
Create Or Replace Function _OCSCHEMA_.recap(smiles Text, label boolean) Returns Text As $EOSQL$
Select string_agg(recap_set,'.') From (Select _OCSCHEMA_.recap_set($1,$2)) as atmp;
$EOSQL$ Language SQL Immutable Returns Null On Null Input;
Comment On FUNCTION _OCSCHEMA_.recap(text, boolean)
Is 'return recap fragments as dot-separated smiles for input smiles with optional isotope labels';
Create Or Replace Function _OCSCHEMA_.recap(smiles Text) Returns Text As $EOSQL$
Select string_agg(recap_set,'.') From (Select _OCSCHEMA_.recap_set($1, false)) As atmp;
$EOSQL$ Language SQL Immutable Returns Null On Null Input;
Comment On FUNCTION _OCSCHEMA_.recap(text)
Is 'return recap fragments as dot-separated smiles for input smiles without isotope labels';