Hi guys,
I am trying to get RMSD between embedded conformers and reference 3D
structure (no Hs) in the input file. Please see example below:
from rdkit import Chem
from rdkit.Chem import AllChem
suppl = Chem.SDMolSupplier("mol.sdf")
for mol in suppl:
# Adding Hs
mh = Chem.AddHs(mol, addCoords=True)
# embedding
AllChem.EmbedMultipleConfs(mh, numConfs=100, numThreads=0)
cids = [conf.GetId() for conf in mh.GetConformers()]
# optimization
for cid in cids:
mp = AllChem.MMFFGetMoleculeProperties(mh, mmffVariant='MMFF94s')
ff = AllChem.MMFFGetMoleculeForceField(mh, mp, confId=cid)
ff.Initialize()
ff.Minimize(maxIts=1000)
# rms to reference structure in input file
m = Chem.RemoveHs(mh)
bestRmsList = []
for i in range(len(cids)):
bestRmsList.append(AllChem.GetBestRMS(m, mol, prbId=cids[i]))
print(min(bestRmsList))
As you can see, Hs are added for embedding and optimization, and then
removed for RMSD calculation (only compare heavy atoms). I pick "m" as
"probe" and "mol" as "reference" (not quite sure about this part).
Is this the right or proper way to do it? Any suggestions? I am still
learning RDKit. Many thanks!
Best,
Leon
_______________________________________________
Rdkit-discuss mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
