Dear JP, On Tue, Nov 29, 2011 at 5:47 PM, JP <jeanpaul.ebe...@inhibox.com> wrote: > > I am not sure that is what I need/I explained myself well. > > My question is: If you have a 3D molecule, and you replace a part of > it with a new group based on some SMARTS expression (with no 3D > coorinates) - what is the best way to keep all the 3D info for the > atoms in the molecule which have not changed while creating some > "sensible" co-cordinates (in the current frame) for the newly appended > atoms to the molecule.
Ok, got it. > I do not see how AddHs will help me when replacing [OH] with [O-], but > that might just be me being thick. > My problem is that ps = AllChem.ReplaceSubstructs(mol, patt, repl, > replaceAll=True) in my original attached code throws away the 3D > coordinates (ps[0] atoms are all on 0,0,0) - and I don't know of a way > to regain the original co-ordinates back. You need two things here: 1) the coordinates of the input molecule should be copied to the result molecule 2) some way to generate sensible coordinates for the new atoms. The first is already supposed to be happening; the fact that it's not is a bug which I will fix (https://sourceforge.net/tracker/?func=detail&aid=3453144&group_id=160139&atid=814650). The second piece is doable using the function AllChem.ConstrainedEmbed() assuming that you can define a core from your original molecule. Here's an example of using AllChem.ConstrainedEmbed() to update the coordinates of part of a molecule: Start by getting a molecule with coordinates: In [16]: m = Chem.MolFromSmiles('CCO') In [17]: AllChem.EmbedMolecule(m) Out[17]: 0 Add an atom to it: In [19]: em = Chem.EditableMol(m) In [21]: em.AddAtom(Chem.Atom(6)) In [22]: em.AddBond(2,3,Chem.BondType.SINGLE) In [23]: m2 = em.GetMol() In [24]: Chem.SanitizeMol(m2) Notice that the atom I just added has zero coords: In [25]: print Chem.MolToMolBlock(m2) RDKit 3D 4 3 0 0 0 0 0 0 0 0999 V2000 1.2251 -0.2656 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0 -0.0350 0.5554 -0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0 -1.1901 -0.2898 -0.0000 O 0 0 0 0 0 0 0 0 0 0 0 0 0.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0 1 2 1 0 2 3 1 0 3 4 1 0 M END Here's how ConstrainedEmbed is used: In [32]: m3=AllChem.ConstrainedEmbed(m2,m) Check that it did what it's supposed to... original coordinates: In [33]: print Chem.MolToMolBlock(m) RDKit 3D 3 2 0 0 0 0 0 0 0 0999 V2000 1.2251 -0.2656 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0 -0.0350 0.5554 -0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0 -1.1901 -0.2898 -0.0000 O 0 0 0 0 0 0 0 0 0 0 0 0 1 2 1 0 2 3 1 0 M END And new coordinates: In [34]: print Chem.MolToMolBlock(m3) RDKit 3D 4 3 0 0 0 0 0 0 0 0999 V2000 1.2253 -0.2708 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0 -0.0370 0.5648 -0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0 -1.1882 -0.2939 0.0000 O 0 0 0 0 0 0 0 0 0 0 0 0 -2.3089 0.6044 -0.0001 C 0 0 0 0 0 0 0 0 0 0 0 0 1 2 1 0 2 3 1 0 3 4 1 0 M END If you take a look at how the tethering in ConstrainedEmbed() is implemented, you could probably skip the distance geometry step and just move right to the force field minimization with the atoms from your "core" tethered. It probably would be a good idea to assign random coordinates to the new atoms (instead of the zeroes they start with) to avoid numeric instabilities. -greg ------------------------------------------------------------------------------ Cloud Services Checklist: Pricing and Packaging Optimization This white paper is intended to serve as a reference, checklist and point of discussion for anyone considering optimizing the pricing and packaging model of a cloud services business. Read Now! http://www.accelacomm.com/jaw/sfnl/114/51491232/ _______________________________________________ Rdkit-discuss mailing list Rdkit-discuss@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
