Search code examples
pythonmachine-learningchemistryrdkitcheminformatics

How to save RDKit conformer object into a sdf file?


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!


Solution

  • 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)