Dear Thomas,
if you wish to align two structures by their MCS, O3A is probably not
the tool for you. O3A is meant for unsupervised alignment, and it will
attempt to align two structures matching the most similar pairs of atoms
between the two. Similarity is defined by the closeness of their MMFF94
atom types and charges. Constraints force the addition of certain atom
pairs in the atom map which is built in an unsupervised fashion. Whether
or not such constraints succeed in giving you the desired alignment
depends on their weight compared to the rest of the atom map.
I think you have two options:
1) Use the AlignMol()method after having identified the MCSs/fragments
of interest in the 2 molecules and using them as atomMap
2) Try to carry out the alignment with O3A in two steps:
a) constrained alignment, global search (options = 0): this will give
you an initial
alignment,which may be a bit skewed by the high constraint weight and
thus sub-optimal
b) unconstrained alignment, local-only search(options = 4): this
might succeed in snapping
in place correctly matching atoms, provided that the initial
alignment is
good enough to allow the local optimization to do its job
HTH,
Paolo
On 01/07/14 09:37, Thomas Strunz wrote:
Hi all,
here the code for Open3DALIGN:
cids = generateConformers(mol, numConformers)
prbPyMP = AllChem.MMFFGetMoleculeProperties(mol)
maxScore = 0;
for refCid in refCids:
for cid in cids:
alignment = AllChem.GetO3A(mol, refMol, prbPyMP, refPyMP,
prbCid=cid, refCid=refCid,)
score = alignment.Score()
logger.debug('Score: %.2f', score)
if score > maxScore:
logger.info('New max. Score: %.2f', score)
maxScore = score
refConformerId = refCid
molCid = cid
alignment = AllChem.GetO3A(mol, refMol, prbPyMP, refPyMP,
prbCid=molCid, refCid=refConformerId)
alignment.Align()
# show in PyMol
One question I have is how I can align this to a specific fragment or
the MCS? I tried using constraintMap but for some reason o3a does not
want to align the 2 structure on it. I suspect because then the rest
of the alignment will not be very good in that case.
The code:
constraintMap = []
constraintWeights = []
mcs = MCS.FindMCS([refMol,mol], bondCompare='bondtypes',
ringMatchesRingOnly=False)
if mcs.completed == 1 and mcs.numAtoms > 0:
core = Chem.MolFromSmarts(mcs.smarts) # or use specific
smarts pattern through argument option
logger.info('MCS: %s', Chem.MolToSmiles(core))
refMatch = refMol.GetSubstructMatch(core)
match = mol.GetSubstructMatch(core)
for idx, val in enumerate(match):
constraintMap.append((val, refMatch[idx]))
constraintWeights.append(100.0)
logger.info(constraintMap)
...
alignment = AllChem.GetO3A(mol, refMol, prbPyMP, refPyMP,
prbCid=cid, refCid=refCid, constraintMap=constraintMap,
constraintWeights=constraintWeights)
The alinment is different but it still does not want to align to my
defined fragment. even if I set constraintWeights to very high value.
In my case I especially want to align heteroatoms properly. For that I
tried:
for idx, val in enumerate(match):
constraintMap.append((val, refMatch[idx]))
if mol.GetAtomWithIdx(val).GetAtomicNum() != 6:
constraintWeights.append(1000.0)
logger.info('Set weigth to 1000')
else:
constraintWeights.append(10.0)
But the single hetero atoms in this my test case are still not aligned
on each other. I also tried ConstrainedEmbed (or how it is called)
but that results in either an error or with relaxed parameters in a
completely useless conformation, namely on a ring.
What are my options?
------------------------------------------------------------------------
From: [email protected]
Date: Fri, 27 Jun 2014 08:20:56 +0200
Subject: Re: [Rdkit-discuss] 3D alignment in Python: align conformers
of 2 molecules
To: [email protected]
CC: [email protected]
On Fri, Jun 27, 2014 at 7:57 AM, Thomas Strunz <[email protected]
<mailto:[email protected]>> wrote:
thanks for your quick reply. This helped to improve the alignment.
I'm glad to hear it!
How can I reproduce the alignment done in with the Open3DAlign
Node in Python? Is it possible at all?
But of course. :-)
There's some example code that shows how to do it on page 37 of
Paolo's presentation from the last RDKit UGM:
https://github.com/rdkit/UGM_2013/raw/master/Presentations/Tosco.RDKit_UGM2013.pdf
-greg
------------------------------------------------------------------------------
Open source business process management suite built on Java and Eclipse
Turn processes into business applications with Bonita BPM Community Edition
Quickly connect people, data, and systems into organized workflows
Winner of BOSSIE, CODIE, OW2 and Gartner awards
http://p.sf.net/sfu/Bonitasoft
_______________________________________________
Rdkit-discuss mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
------------------------------------------------------------------------------
Open source business process management suite built on Java and Eclipse
Turn processes into business applications with Bonita BPM Community Edition
Quickly connect people, data, and systems into organized workflows
Winner of BOSSIE, CODIE, OW2 and Gartner awards
http://p.sf.net/sfu/Bonitasoft
_______________________________________________
Rdkit-discuss mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
