Hello I’m trying to calculate RMSDs between a reference molecule from an X-ray structure and a series of similar but somewhat different molecules with a common scaffold that I have docked to the protein. I don’t want to align the structures, I just want to measure how different their docked poses are to the reference structure. All of the molecules have the same core structure, but have different functional groups. This functionality is available in the recently developed web server LigRMSD (https://ligrmsd.appsbio.utalca.cl), but I need to make some modifications to the procedure.
Following this thread (https://sourceforge.net/p/rdkit/mailman/message/35678724/) and this blog post (http://rdkit.blogspot.com/2013/12/using-allchemconstrainedembed.html), I’ve written the following code: from rdkit import Chem from rdkit.Chem import AllChem from rdkit.Chem.Draw import IPythonConsole from rdkit.Chem.Draw.MolDrawing import MolDrawing, DrawingOptions from rdkit.Chem import rdFMCS import glob # Load a set of mol files. mols = [] for m in glob.glob('*.mol'): tmpmol = Chem.MolFromMolFile(m) name_without_suffix = m[:-4] tmpmol.SetProp('_Name', name_without_suffix) mols.append(tmpmol) # Draw the molecules. Chem.Draw.MolsToGridImage(mols, molsPerRow=2, subImgSize=(350, 350), legends=[mol.GetProp('_Name') for mol in mols]) # Find the MCS. mcs = rdFMCS.FindMCS(mols,threshold=1.0, ringMatchesRingOnly=True, completeRingsOnly=True) # Align everything to the reference molecule. patt = Chem.MolFromSmarts(mcs.smartsString) refMol = mols[0] refMatch = refMol.GetSubstructMatch(patt) rmsVs = [] for probeMol in mols[1:]: mv = probeMol.GetSubstructMatch(patt) rms = AllChem.AlignMol(probeMol, refMol, atomMap=list(zip(mv, refMatch)), maxIters=0) print(probeMol.GetProp('_Name'), f'{rms:.2f}') rmsVs.append(rms) However, even though I set maxIters=0 in the call to AlignMol, the molecules have all been translated to a common origin (which is not what I want). I tried replacing that line of code with the following call to CalcRMS (using the same atom mapping): rms = AllChem.CalcRMS(probeMol, refMol, map=list(zip(mv, refMatch))) but then I got the following error: --------------------------------------------------------------------------- ValueError Traceback (most recent call last) <ipython-input-7-0faef2f852d9> in <module>() 8 mv = probeMol.GetSubstructMatch(patt) 9 #rms = AllChem.AlignMol(probeMol, refMol, atomMap=list(zip(mv, refMatch)), maxIters=0) ---> 10 rms = AllChem.CalcRMS(probeMol, refMol, map=list(zip(mv, refMatch))) 11 print(probeMol.GetProp('_Name'), f'{rms:.2f}') 12 rmsVs.append(rms) ValueError: sequence does not support length query Any suggestions for fixing the code would be appreciated. Thank you, Jerry Parks
_______________________________________________ Rdkit-discuss mailing list Rdkit-discuss@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
