3

I generated a bunch of conformers for a molecule. For each conformed, I want to save the coordinates in a SDF file. I tried the following, but the coordinates in the sdf file is different from that of the conformer.

from rdkit import Chem
from rdkit.Chem import AllChem
from rdkit.Chem.EnumerateStereoisomers import EnumerateStereoisomers
from rdkit.Chem.EnumerateStereoisomers import StereoEnumerationOptions


aspirin = 'CC(=O)OC1=CC=CC=C1C(=O)O'
mol = Chem.AddHs(Chem.MolFromSmiles(aspirin))

conformers = AllChem.EmbedMultipleConfs(mol, numConfs=8)

conformer_i = mol.GetConformer(0)  # This is the conformer that I wanted to save in a SDF file.
print(conformer_i.GetPositions())

I get

[[ 3.55910965  0.20840444  0.20174025]
 [ 2.13148034 -0.07274183 -0.11680036]
 [ 1.81233825 -0.09134332 -1.32465968]
 [ 1.16950868 -0.30747936  0.83214171]
 [-0.12944319 -0.55824117  0.41262056]
 [-0.53676538 -1.86517356  0.15656512]
 [-1.82111824 -2.14574436 -0.26154327]
 [-2.69581394 -1.09437287 -0.4203516 ]
 [-2.35275981  0.20578091 -0.18448083]
 [-1.05455197  0.45649503  0.23480361]
 [-0.64543959  1.83612099  0.49878179]
 [ 0.50523818  2.10782743  0.87545738]
 [-1.54126908  2.879518    0.33273924]
 [ 4.15835976 -0.73143248  0.14404327]
 [ 3.62587517  0.63489849  1.2047807 ]
 [ 3.98056046  0.90402053 -0.55673262]
 [ 0.20229068 -2.6478098   0.3014611 ]
 [-2.09903401 -3.17965212 -0.45066362]
 [-3.71036497 -1.30377103 -0.74938855]
 [-3.04003996  1.04851376 -0.30756941]
 [-1.51816102  3.71618233  0.87474425]]

But when I tried to save the molecule with this conformer structure into a SDF file with,

w = Chem.SDWriter('conformer.sdf')
mol.AddConformer(conformer_i)
w.write(mol)
w.close()

I got the following:


     RDKit          3D

 21 21  0  0  0  0  0  0  0  0999 V2000
   -2.8028   -2.2971    0.2528 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.9196   -1.1193    0.3821 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.9380   -0.4309    1.4053 O   0  0  0  0  0  0  0  0  0  0  0  0
   -1.0329   -0.6957   -0.5799 O   0  0  0  0  0  0  0  0  0  0  0  0
   -0.2038    0.3759   -0.5076 C   0  0  0  0  0  0  0  0  0  0  0  0
   -0.5008    1.6592   -0.9114 C   0  0  0  0  0  0  0  0  0  0  0  0
    0.4068    2.6946   -0.8011 C   0  0  0  0  0  0  0  0  0  0  0  0
    1.6647    2.4889   -0.2762 C   0  0  0  0  0  0  0  0  0  0  0  0
    1.9781    1.2093    0.1330 C   0  0  0  0  0  0  0  0  0  0  0  0
    1.0695    0.1765    0.0217 C   0  0  0  0  0  0  0  0  0  0  0  0
    1.4778   -1.1541    0.4786 C   0  0  0  0  0  0  0  0  0  0  0  0
    0.6542   -2.1230    0.3851 O   0  0  0  0  0  0  0  0  0  0  0  0
    2.7241   -1.3647    1.0002 O   0  0  0  0  0  0  0  0  0  0  0  0
   -2.9017   -2.5297   -0.8210 H   0  0  0  0  0  0  0  0  0  0  0  0
   -3.7860   -2.0195    0.6858 H   0  0  0  0  0  0  0  0  0  0  0  0
   -2.4405   -3.1762    0.7986 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.4956    1.7914   -1.3196 H   0  0  0  0  0  0  0  0  0  0  0  0
    0.1087    3.6835   -1.1386 H   0  0  0  0  0  0  0  0  0  0  0  0
    2.3837    3.2960   -0.1858 H   0  0  0  0  0  0  0  0  0  0  0  0
    2.9838    1.0330    0.5554 H   0  0  0  0  0  0  0  0  0  0  0  0
    3.5700   -1.4979    0.4428 H   0  0  0  0  0  0  0  0  0  0  0  0
  1  2  1  0
  2  3  2  0
  2  4  1  0
  4  5  1  0
  5  6  2  0
  6  7  1  0
  7  8  2  0
  8  9  1  0
  9 10  2  0
 10 11  1  0
 11 12  2  0
 11 13  1  0
 10  5  1  0
  1 14  1  0
  1 15  1  0
  1 16  1  0
  6 17  1  0
  7 18  1  0
  8 19  1  0
  9 20  1  0
 13 21  1  0
M  END
$$$$

The coordinates in the molecular sdf file is different from the coordinates in the conformer_i. Does anyone have a insight about this issue? Thank you!

Zhen Liu
  • 57
  • 1
  • 5

1 Answers1

4

When you use SDWriter.write you need to supply the ID of the conformer you wish to write to the file:

writer = Chem.SDWriter('aspirin_confs.sdf')
for cid in range(mol.GetNumConformers()):
    writer.write(mol, confId=cid)

Edit:

If you are only interested in writing this property to the file then why not just overwrite the molecule property each iteration i.e.

writer = Chem.SDWriter('aspirin_confs.sdf')
for cid in range(mol.GetNumConformers()):
    mol.SetProp('ID', f'aspirin_conformer_{cid}')
    writer.write(mol, confId=cid)
Oliver Scott
  • 1,673
  • 8
  • 17
  • Great, this worked! BTW, Do you know how to add an ID field in the SDF file for each conformer? I wanna be able to distinguish them by giving each conformed a name (or so-called ID). I thought to convert each `conformer` into a `rdkit molecule object`, then use `mol.SetProp('ID', 'aspirin_conformer_1')`. But I don't know how to do convert conformer to molecule yet... – Zhen Liu Oct 19 '21 at 14:03
  • No problem, if this answered your question please mark it as the accepted answer. Also, check the edit for an answer to your other question. If you would like to know how to convert a conformed into a molecule then you should open another question so others may also find the answer. – Oliver Scott Oct 19 '21 at 16:28