Another cool application of the distance bounds matrix: Nik Stiefl did a
nice tutorial on generating conformations that satisfy a 3D pharmacophore :
https://github.com/rdkit/UGM_2016/blob/master/Notebooks/Stiefl_RDKitPh4FullPublication.ipynb

On Thu, Mar 1, 2018 at 8:56 AM, Jan Halborg Jensen <jhjen...@chem.ku.dk>
wrote:

> Hi Brian
>
> My interest stems from this paper http://dx.doi.org/10.1039/c5cp04706d
> where they alter the bounds matrix to impose constraints on the embedding.
>
> In principle this could also be done with EmbedMolecule and coordMap but
> that doesn't seem to be working
> https://sourceforge.net/p/rdkit/mailman/message/36139015/
>
> Best regards, Jan
> ------------------------------
> *From:* Bennion, Brian [benni...@llnl.gov]
> *Sent:* Wednesday, February 28, 2018 7:09 PM
> *To:* Greg Landrum
>
> *Cc:* rdkit-discuss@lists.sourceforge.net
> *Subject:* Re: [Rdkit-discuss] problems with EmbedMol
>
> Hello Greg and Jan,
>
>
> This is a real newbie question, but what is the use case for this
> function?  Is it used to generate all possible connections (limited by some
> distance) between 3 or more atoms given in a smiles string?
>
>
> Brian
>
>
> ------------------------------
> *From:* Greg Landrum <greg.land...@gmail.com>
> *Sent:* Wednesday, February 28, 2018 8:53:59 AM
> *To:* Jan Halborg Jensen
> *Cc:* rdkit-discuss@lists.sourceforge.net
> *Subject:* Re: [Rdkit-discuss] problems with EmbedMol
>
> Hi Jan,
>
> It took me much longer than it should have to figure this one out...
>
> The bounds matrix that is returned by GetMoleculeBoundsMatrix() needs to
> have triangle bounds smoothing applied to it before it can be embedded. The
> bounds smoothing process narrows the possible distance ranges between the
> atoms. Here's a quick demo of that.
>
> We start with your example:
>
> In [19]: mol = Chem.MolFromSmiles("CCC")
>     ...: mol = Chem.AddHs(mol)
>     ...: bounds = AllChem.GetMoleculeBoundsMatrix(mol)
>     ...: EmbedLib.EmbedMol(mol,bounds)
>     ...:
> ------------------------------------------------------------
> ---------------
> ValueError                                Traceback (most recent call last)
> <ipython-input-19-60befaacbe48> in <module>()
>       2 mol = Chem.AddHs(mol)
>       3 bounds = AllChem.GetMoleculeBoundsMatrix(mol)
> ----> 4 EmbedLib.EmbedMol(mol,bounds)
>
> c:\Users\glandrum\RDKit_git\rdkit\Chem\Pharm3D\EmbedLib.py in
> EmbedMol(mol, bm, atomMatch, weight, randomSeed, excludedVolumes)
>     183       for i in range(nAts):
>     184         weights.append((i, idx, weight))
> --> 185   coords = DG.EmbedBoundsMatrix(bm, weights=weights,
> numZeroFail=1, randomSeed=randomSeed)
>     186   # for row in coords:
>     187   #  print(', '.join(['%.2f'%x for x in row]))
>
> ValueError: could not embed matrix
>
>
>
> But if we do the triangle bounds smoothing things embed without problems:
>
> In [20]: from rdkit import DistanceGeometry
>
> In [21]: DistanceGeometry.DoTriangleSmoothing(bounds)
> Out[21]: True
>
> In [22]: EmbedLib.EmbedMol(mol,bounds)
>
> In [23]:
>
>
> There is a good argument to be made for GetMoleculeBoundsMatrix()
> returning the smoothed bounds matrix by default. I'll put that on the list
> for the next release.
>
> Best,
> -greg
>
>
>
> On Wed, Feb 28, 2018 at 10:41 AM, Jan Halborg Jensen <jhjen...@chem.ku.dk>
> wrote:
>
> The following code works fine with ethane (CC) but for propane (CCC) or
> anything else I get the following error
> ValueError: could not embed matrix
>
> Any ideas or solutions would be appreciated
>
> Best regards, Jan
>
>
> from rdkit import Chem
> from rdkit.Chem import AllChem
> from rdkit.Chem.Pharm3D import EmbedLib
>
>
> mol = Chem.MolFromSmiles("CCC")
> mol = Chem.AddHs(mol)
> bounds = AllChem.GetMoleculeBoundsMatrix(mol)
>
> EmbedLib.EmbedMol(mol,bounds)
> EmbedLib.OptimizeMol(mol, bounds)
>
> ------------------------------------------------------------
> ------------------
> Check out the vibrant tech community on one of the world's most
> engaging tech sites, Slashdot.org! http://sdm.link/slashdot
> _______________________________________________
> Rdkit-discuss mailing list
> Rdkit-discuss@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
>
>
>
------------------------------------------------------------------------------
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot
_______________________________________________
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss

Reply via email to