I believe (Greg can correct me) to align the bemis-murcko scaffold, you
could

(1) extract the original atom pairs and send them to the RMSD algorithm
(2) take a bemis scaffold and convert it to a substructure query for use in
the RMSD algorithm.  In either case the RMSD is the rmsd of the scaffold
atoms, not the rest of the molecule.  Below is a little snippet that I
believes does this.

from rdkit.Chem import AllChem
from rdkit.Chem import rdqueries, rdMolAlign

mol = Chem.MolFromMolBlock( mol block )
print Chem.MolToSmiles(m)

murcko = AllChem.MurckoDecompose(mol)
print "murko",  Chem.MolToSmiles(murcko)

# bemis scaffolds match everything
q = rdqueries.AtomNumGreaterQueryAtom(0)
bemis = Chem.RWMol(murcko)
for atom in bemis.GetAtoms():
   bemis.ReplaceAtom(atom.GetIdx(), q)

rmsd = Chem.rdMolAlign.AlignMol( bemis, m )
print rmsd





On Mon, Feb 20, 2017 at 11:21 AM, Greg Landrum <greg.land...@gmail.com>
wrote:

> HI Thomas,
>
> To be sure we're talking about the same thing: rdMolAlign.GetO3A() is an
> implementation of the Open3DAlign algorithm. This is an unsupervised
> approach that uses a clever algorithm to come up with an atom-atom mapping
> between the two molecules you give it. It's not always going to pick the
> same atoms to align that you would.
>
> To answer the original question: if the two molecules you want to align do
> not share the same scaffold (or at least have a lot in common in the core
> of the molecule), it's unlikely that an MCS-based alignment is going to
> help.
>
> To answer your direct question here, the scaffold finding code in
> rdkit.Chem should be preserving coordinates. Here's a simple demonstration
> of that:
>
> In [3]: m =Chem.AddHs(Chem.MolFromSmiles('CC1CO1'))
>
> In [4]: AllChem.EmbedMolecule(m,AllChem.ETKDG())
> Out[4]: 0
>
> In [5]: nh = Chem.RemoveHs(m)
>
> In [6]: print(Chem.MolToMolBlock(nh))
>
>      RDKit          3D
>
>   4  4  0  0  0  0  0  0  0  0999 V2000
>    -1.1850    0.2738   -0.1814 C   0  0  0  0  0  0  0  0  0  0  0  0
>     0.0992   -0.4987   -0.2391 C   0  0  0  0  0  0  0  0  0  0  0  0
>     1.3229    0.1241    0.2290 C   0  0  0  0  0  0  0  0  0  0  0  0
>     0.6874   -0.7558    1.1165 O   0  0  0  0  0  0  0  0  0  0  0  0
>   1  2  1  0
>   2  3  1  0
>   3  4  1  0
>   4  2  1  0
> M  END
>
>
> In [7]: from rdkit.Chem import MurckoDecompose
>
> In [8]: scaff = Chem.MurckoDecompose(nh)
>
> In [10]: print(Chem.MolToMolBlock(scaff))
>
>      RDKit          3D
>
>   3  3  0  0  0  0  0  0  0  0999 V2000
>     0.0992   -0.4987   -0.2391 C   0  0  0  0  0  0  0  0  0  0  0  0
>     1.3229    0.1241    0.2290 C   0  0  0  0  0  0  0  0  0  0  0  0
>     0.6874   -0.7558    1.1165 O   0  0  0  0  0  0  0  0  0  0  0  0
>   1  2  1  0
>   2  3  1  0
>   3  1  1  0
> M  END
>
>
> You could help me provide a better answer here by providing a couple of
> example SDFs that you'd like to align, ideally together with a bit of
> sample code showing what you have tried that produces alignments you are
> unhappy with.
>
>
> -greg
>
>
>
>
> On Mon, Feb 20, 2017 at 1:54 PM, Thomas Evangelidis <teva...@gmail.com>
> wrote:
>
>> As a follow up question on this topic, I would like to ask if
>> MurckoScaffold.GetScaffoldForMol(mol) returns the scaffold of mol with
>> different coordinates?
>> I am asking this because when I use the transformation matrix of the
>> alignment of the cores of the probe and the reference molecules, in order
>> to align the whole probe to the reference molecule, the two molecules don't
>> seem to be aligned (they are in distance). Basically I do this:
>>
>> qcore = MurckoScaffold.GetScaffoldForMol(qmol)
>> refcore = MurckoScaffold.GetScaffoldForMol(refmol)
>> pyO3A = rdMolAlign.GetO3A(qcore, refcore, prbCid=qconfID, refCid=0,
>> reflect=True)
>> AllChem.TransformMol(qmol, bestRMSDTrans[1], confId=bestconfID,
>> keepConfs=False)
>>
>> and then I write the qmol in an sdf file. But when I visualize it the
>> qmol is far from the refmol!
>>
>>
>>
>>
>>
>> On 20 February 2017 at 02:33, Thomas Evangelidis <teva...@gmail.com>
>> wrote:
>>
>>> Dear all,
>>>
>>> I want to align 250 compounds that binding to the same pocket to one of
>>> the 9 available crystal ligands. I chose the reference ligand based on the
>>> Morgan2 similarity to the probe molecule. Then I align the 2 compounds
>>> using:
>>>
>>> pyO3A = rdMolAlign.GetO3A(qmol, refmol, prbCid=qconfID, refCid=0,
>>> reflect=True)
>>> RMSD = pyO3A.Align()
>>>
>>> ​and keep only the conformer of the probe with the lowest RMSD to the
>>> reference compound. However, the alignment looks terrible when I
>>> visualize it, so I would like to ask if there is any way to align the
>>> maximum common substructure only. I tried to align only the core of both
>>> molecules as defined by MurckoScaffold.GetScaffoldForMol(mol)​, but
>>> still the alignment looks bad. I have seen in the documentation how to find
>>> the maximum common substructure with rdFMCS.FindMCS but before I engage
>>> into programming it I would like to know if there is any automatic way to
>>> find it on the fly while aligning the 2 molecules.
>>>
>>>
>>>
>>> --
>>>
>>> ======================================================================
>>>
>>> Thomas Evangelidis
>>>
>>> Research Specialist
>>> CEITEC - Central European Institute of Technology
>>> Masaryk University
>>> Kamenice 5/A35/1S081,
>>> 62500 Brno, Czech Republic
>>>
>>> email: tev...@pharm.uoa.gr
>>>
>>>           teva...@gmail.com
>>>
>>>
>>> website: https://sites.google.com/site/thomasevangelidishomepage/
>>>
>>>
>>
>>
>> --
>>
>> ======================================================================
>>
>> Thomas Evangelidis
>>
>> Research Specialist
>> CEITEC - Central European Institute of Technology
>> Masaryk University
>> Kamenice 5/A35/1S081,
>> 62500 Brno, Czech Republic
>>
>> email: tev...@pharm.uoa.gr
>>
>>           teva...@gmail.com
>>
>>
>> website: https://sites.google.com/site/thomasevangelidishomepage/
>>
>>
>> ------------------------------------------------------------
>> ------------------
>> 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
>
>
------------------------------------------------------------------------------
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