Dear Noel,
thanks again, this is of great help.
I think the new code is somehow faulty and for some reason it returns
1.79769313486e+308 for any molecule. The updated code is posted at the
end of the message.
The old code worked for what I thought it was molecules with different
labelling. By adding this to the old code at the end:
for mols in mol:
gs = ob.OBGraphSym(mols.OBMol)
gs.GetSymmetry(symclasses)
ob.CanonicalLabels(mols.OBMol, symclasses, canonlabels)
print list(canonlabels)
I could see where the alignment did not work, e.g. for a working alignment:
0.124327784126
[5, 9, 12, 14, 11, 8, 4, 6, 10, 15, 20, 18, 21, 24, 25, 23, 2, 1, 19,
22, 3, 17, 7, 16, 13]
[5, 9, 12, 14, 11, 8, 4, 6, 10, 15, 20, 18, 21, 24, 25, 23, 2, 1, 19,
22, 3, 17, 7, 16, 13]
and for a non-working alignment:
3.93476583134
[5, 9, 12, 14, 11, 8, 4, 6, 10, 15, 20, 18, 21, 24, 25, 23, 2, 1, 19,
22, 3, 17, 7, 16, 13]
[5, 9, 12, 14, 11, 8, 4, 6, 10, 15, 20, 18, 21, 24, 25, 2, 1, 19, 22, 3,
17, 7, 16, 13, 23]
So I guess, after all, even with the new method (which does not involve
deleting atoms from the OBMol), I should not get a proper alignment?
This is the updated code:
import pybel
ob = pybel.ob
mol = [mola,molb]
smarts = [pybel.Smarts("[OH]C=O"),pybel.Smarts("C=C(C#N)C(=O)O")]
removeanything = True
targets = [[],[]]
c = 0
if removeanything:
for mols in mol:
found = False
for ligands in smarts:
smartsearch = ligands.findall(mols)
smartsearch = [item for sublist in smartsearch for item in
sublist]
if not smartsearch == [] and not found:
for atom in mols:
if atom.OBAtom.GetIdx() not in smartsearch:
targets[c].append(atom.vector)
c += 1
found = True
if not found:
for atom in mols:
targets[c].append(atom.vector)
c += 1
vtarget = ob.vectorVector3(targets[0])
vreference = ob.vectorVector3(targets[1])
align = ob.OBAlign()
align.SetRef(vtarget)
align.SetTarget(vreference)
align.Align()
print align.GetRMSD()
Thank you.
Best,
Pepe
On 03/06/14 15:56, Noel O'Boyle wrote:
Hi Giulio,
Problem 1: check out vectorVector3. Initialized with a list of
vector3, e.g. my_vv = ob.vectorVector3([a, b])
Problem 3: What I meant is that the RMSD only takes into account the
aligned coordinates. If you want to include the RMSD of the other
coordinates, you need to work it out yourself.
Problem 2: Your pseduocode looks fine to me. Maybe use a set for
atomlist for a miniscule speedup. I hope all of your molecules are the
same though, and the atoms are in the same order, as otherwise it
doesn't make sense as the target and reference will involve things
which don't correspond.
- Noel
On 3 June 2014 15:25, Giulio Pepe <gp...@cam.ac.uk
<mailto:gp...@cam.ac.uk>> wrote:
Dear Noel,
just a few problems I am stuck with, using this procedure:
- SetRef() and SetTarget() need a vector of vector3 objects as
argument. While vector3 is an OpenBabel class, vector is not (I
assume it is the C++ vector?). How do I feed them a vector? Is a
python list going to be fine?
- The match to the SMARTS is to exclude from the alignment, so did
I get it right if the one below is the procedure to do it, or is
there a SMARTer (sorry, the pun was too tempting) way? Although
this seems similar to the old (non-working) way.
ob = pybel.ob
smarts = pybel.Smarts("[OH]C=O")
atomlist = []
for i in smarts.findall(mol):
for j in i:
atomlist.append(j),
for atom in mol:
if atom.OBAtom.GetIdx() not in atomlist:
reference.append(atom.vector)
ob.OBAlign().SetTarget(reference)
...
- Why ob.OBAlign().GetRMSD() should not work after using
ob.OBAlign().Align() with this procedure? I am not sure why shall
I apply the rotation matrix.
Thanks a lot for your help!
Best,
Pepe
On 30/05/14 19:46, Noel O'Boyle wrote:
Try using SetRef() and SetTarget() instead with just the
coordinates for the atoms you want to align (in the same order
for both, as given by the match to the Smarts). Afterwards, if
you want the RMSD for the whole molecule, you can manually
calculate the whole molecule RMSD by translating one of the
molecules to its centroid and applying the rotation matrix. (I
can explain more if this is confusing.)
- Noel
On 24 May 2014 13:14, Giulio Pepe <gp...@cam.ac.uk
<mailto:gp...@cam.ac.uk>> wrote:
Hi there,
I am building a script in python that removes a certain
moiety from two
molecules and calculates the RMSD between the remaining
structure.
I am doing this by asking for one or more SMARTS strings to
find in the
molecules, labelling these atoms and deleting them using
OBMol.DeleteAtom(). Then, I align the molecules and find the
RMSD using
the OBAlign() class.
This works fine for some SMARTS string, but for others it
gives a very
high and unrealistic RMSD value. I investigated the reason
for this and
I found out that this is caused by a bad alignment of the two
molecules.
I suspect that this is caused by a disordering of the atom
labels caused
by removing certain SMARTS strings. I tested this by printing the
canonical labels.
I tried to renumber the atoms using OBMol.RenumberAtoms(),
unsuccessfully. Symmetry should not be a problem in these
molecules.
Did anyone come across this, or have a better solution? I am
attaching a
simplified version of my code below:
import pybel
ob = pybel.ob
mol = [mola,molb]
smarts = [pybel.Smarts("[OH]C=O"),pybel.Smarts("C=C(C#N)C(=O)O")]
for mols in mol:
mols.OBMol.StripSalts()
mols.removeh()
if removemoiety:
for ligands in smarts:
i = 0
todelete = []
for mols in mol:
todelete.append([])
for atom in mols:
for l in ligands.findall(mol[i]):
if atom.OBAtom.GetIdx() in l:
todelete[i].append(atom.OBAtom)
i += 1
j = 0
for mols in todelete:
for atom in mols:
mol[j].OBMol.DeleteAtom(atom)
j += 1
align = ob.OBAlign()
align.SetRefMol(mol[0].OBMol)
align.SetTargetMol(mol[1].OBMol)
align.Align()
print align.GetRMSD()
Thank you very much for your help!
Best,
Pepe
------------------------------------------------------------------------------
"Accelerate Dev Cycles with Automated Cross-Browser Testing -
For FREE
Instantly run your Selenium tests across 300+ browser/OS combos.
Get unparalleled scalability from the best Selenium testing
platform available
Simple to use. Nothing to install. Get started now for free."
http://p.sf.net/sfu/SauceLabs
_______________________________________________
OpenBabel-discuss mailing list
OpenBabel-discuss@lists.sourceforge.net
<mailto:OpenBabel-discuss@lists.sourceforge.net>
https://lists.sourceforge.net/lists/listinfo/openbabel-discuss
--
PhD Student
Structure and Dynamics
University of Cambridge
Department of Physics
Cavendish Laboratory
J J Thomson Avenue
CB3 0HE, UK
Tel: +44 (0) 1223 337011 (Room 425)
Fax: +44 (0) 1223 350266
e-mail:gp...@cam.ac.uk <mailto:gp...@cam.ac.uk>
web:www.giuliopepe.com <http://www.giuliopepe.com>
------------------------------------------------------------------------------
Learn Graph Databases - Download FREE O'Reilly Book
"Graph Databases" is the definitive new guide to graph databases and their
applications. Written by three acclaimed leaders in the field,
this first edition is now available. Download your free book today!
http://p.sf.net/sfu/NeoTech
_______________________________________________
OpenBabel-discuss mailing list
OpenBabel-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/openbabel-discuss
