Sorry for the repeated posting -- I had asked a similar question for
2D some time ago
(http://www.mail-archive.com/[email protected]/msg01916.html),
but now I need the 3D case and I am stuck.
Using RDKit 2011_09_01.
I am simply protonating a molecule via some simple SMART replacements.
I would like the final output molecule to have the same 3D co
ordinates as the input one. I have code which works for the 2D case,
(basically the below but using
GenerateDepictionMatching2DStructure(...) instead.
A copy and paste example:
#!/usr/bin/env python
"""
Protonates the input molecules...
@author: JP
"""
from rdkit import Chem
from rdkit.Chem import AllChem
orig_mol = """LIGAND
Insight 3D 0
47 48 0 0 0 0 0
16.1500 22.3350 42.8110 C 0 0 0 0 0 0
15.3160 21.8930 41.9800 O 0 0 0 0 0 0
16.9980 23.2070 42.5060 O 0 0 0 0 0 0
16.1410 21.7710 44.2030 C 0 0 0 0 0 0
15.1160 22.5540 45.0230 C 0 0 0 0 0 0
15.0619 22.1799 46.0577 H 0 0 0 0 0 0
15.3649 23.6269 45.0656 H 0 0 0 0 0 0
14.1065 22.4630 44.5883 H 0 0 0 0 0 0
17.5310 21.9290 44.8610 C 0 0 0 0 0 0
17.5541 21.4127 45.8334 H 0 0 0 0 0 0
18.3473 21.5168 44.2477 H 0 0 0 0 0 0
17.7642 22.9903 45.0477 H 0 0 0 0 0 0
15.7170 20.3840 44.1390 O 0 0 0 0 0 0
16.5230 19.2430 43.6250 C 0 0 0 0 0 0
16.2460 17.9700 44.1170 C 0 0 0 0 0 0
15.4668 17.8367 44.8675 H 0 0 0 0 0 0
16.9520 16.8770 43.6600 C 0 0 0 0 0 0
16.7273 15.8853 44.0527 H 0 0 0 0 0 0
17.9400 17.0290 42.7100 C 0 0 0 0 0 0
17.5250 19.4050 42.6630 C 0 0 0 0 0 0
17.7525 20.3959 42.2701 H 0 0 0 0 0 0
18.2350 18.2860 42.2070 C 0 0 0 0 0 0
19.0179 18.4043 41.4579 H 0 0 0 0 0 0
18.6450 15.8760 42.2400 C 0 0 0 0 0 0
18.9367 15.2468 43.1044 H 0 0 0 0 0 0
19.6099 16.1765 41.7887 H 0 0 0 0 0 0
17.8550 15.0410 41.2140 C 0 0 0 0 0 0
16.7450 15.4130 40.8360 O 0 0 0 0 0 0
18.2860 13.8040 41.0020 N 0 0 0 0 0 0
19.1093 13.4965 41.5391 H 0 0 0 0 0 0
17.7310 12.8920 40.1270 C 0 0 0 0 0 0
16.3370 12.7960 39.9850 C 0 0 0 0 0 0
15.6892 13.5478 40.4357 H 0 0 0 0 0 0
15.7870 11.7440 39.2710 C 0 0 0 0 0 0
16.6210 10.7910 38.7010 C 0 0 0 0 0 0
16.1929 9.9591 38.1417 H 0 0 0 0 0 0
18.5490 11.9350 39.5470 C 0 0 0 0 0 0
19.6312 12.0105 39.6525 H 0 0 0 0 0 0
18.0040 10.8950 38.8420 C 0 0 0 0 0 0
14.3030 11.6400 39.1100 C 0 0 0 0 0 0
13.8285 12.6314 39.1672 H 0 0 0 0 0 0
14.0383 11.1938 38.1399 H 0 0 0 0 0 0
13.8615 11.0120 39.9008 H 0 0 0 0 0 0
18.9090 9.8760 38.2190 C 0 0 0 0 0 0
18.4728 9.4727 37.2928 H 0 0 0 0 0 0
19.8880 10.3114 37.9693 H 0 0 0 0 0 0
19.0811 9.0290 38.9030 H 0 0 0 0 0 0
1 2 1 0 0 0
1 3 2 0 0 0
1 4 1 0 0 0
4 5 1 0 0 0
4 9 1 0 0 0
4 13 1 0 0 0
5 6 1 0 0 0
5 7 1 0 0 0
5 8 1 0 0 0
9 10 1 0 0 0
9 11 1 0 0 0
9 12 1 0 0 0
13 14 1 0 0 0
14 15 2 0 0 0
14 20 1 0 0 0
15 17 1 0 0 0
15 16 1 0 0 0
17 19 2 0 0 0
17 18 1 0 0 0
19 22 1 0 0 0
19 24 1 0 0 0
20 22 2 0 0 0
20 21 1 0 0 0
22 23 1 0 0 0
24 27 1 0 0 0
24 25 1 0 0 0
24 26 1 0 0 0
27 28 2 0 0 0
27 29 1 0 0 0
29 31 1 0 0 0
29 30 1 0 0 0
31 32 1 0 0 0
31 37 2 0 0 0
32 34 2 0 0 0
32 33 1 0 0 0
34 35 1 0 0 0
34 40 1 0 0 0
35 39 2 0 0 0
35 36 1 0 0 0
37 39 1 0 0 0
37 38 1 0 0 0
39 44 1 0 0 0
40 41 1 0 0 0
40 42 1 0 0 0
40 43 1 0 0 0
44 45 1 0 0 0
44 46 1 0 0 0
44 47 1 0 0 0
M END"""
rules=(
# Carboxylic Acid ionization
('[O-]','[$([OH]C(=O))]'),
)
def protonate(mol):
Chem.SanitizeMol(mol)
orig = Chem.MolToSmiles(mol)
for repl, patt in rules:
repl = Chem.MolFromSmiles(repl, sanitize=False)
patt = Chem.MolFromSmarts(patt)
ps = AllChem.ReplaceSubstructs(mol, patt, repl, replaceAll=True)
ps[0].UpdatePropertyCache(False)
Chem.GetSymmSSSR(ps[0])
new = Chem.MolToSmiles(ps[0])
modified = Chem.MolToSmiles(mol) != Chem.MolToSmiles(ps[0])
print "Input: %s, Output: %s, Deprotonated: %s" % (orig,
new, modified)
if modified:
# emebed these as
# the following keeps the co ordinates in mol, and just
adds co ordinates
# for the new (replaced) atoms
AllChem.GenerateDepictionMatching3DStructure(ps[0], mol)
mol = ps[0]
orig = Chem.MolToSmiles(mol)
# Finally, clean up and write smiles and sdf versions
Chem.SanitizeMol(mol)
return mol
mol = Chem.MolFromMolBlock(orig_mol)
protonated_mol = protonate(mol)
# WHERE ARE MY ORIG CO-ORDINATES !?
print Chem.MolToMolBlock(protonated_mol)
Can someone enlighten me what is going on? I can see that the final
line in AllChem.GenerateDepictionMatching3DStructure(...) is
Compute2DCoordsMimicDistmat(...) which doesn't bode well for me.
Is there some other way to do this?
-
Jean-Paul Ebejer
Early Stage Researcher
------------------------------------------------------------------------------
All the data continuously generated in your IT infrastructure
contains a definitive record of customers, application performance,
security threats, fraudulent activity, and more. Splunk takes this
data and makes sense of it. IT sense. And common sense.
http://p.sf.net/sfu/splunk-novd2d
_______________________________________________
Rdkit-discuss mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss