9
$\begingroup$

Given a SMILES string for a molecule, I want to know the elements in the molecule, I then send this to basis_set_exchange to get a "custom" basis set for use in Gaussian.

This seems like something that RDKit can do, or PubChem.

I do know that given a molecular formula, like "H2O" molmass can return the elements. So, it would also suffice to have a program that simply takes a SMILES string and returns a molecular formula, then use molmass to get the elements.

However, this seems like such a trivial exercise, it must already be done. I don't want to make my own script for converting a SMILES to its elements just because I don't want to have to worry about corner cases or reinvent a wheel.

$\endgroup$
7
  • $\begingroup$ +1 because, Charlie Crown is back! By the way, did you ever visit the Gaussian chat room? $\endgroup$ – Nike Dattani Mar 27 at 5:09
  • $\begingroup$ Now I did, I liked the historical context in that question. Personally, I am going to try psi4 since it is python friendly $\endgroup$ – B. Kelly Mar 27 at 5:13
  • $\begingroup$ is smiles ID supposed to mean smiles string? I've never heard of a smiles ID $\endgroup$ – Cody Aldaz Mar 27 at 6:06
  • $\begingroup$ Probably, I may be misusing the term ID in this case. "CCCCCO" is the SMILES for pentanol, for instance $\endgroup$ – B. Kelly Mar 27 at 6:09
  • 1
    $\begingroup$ haha, sorry didn't mean to make you feel bad ;) was just confused. $\endgroup$ – Cody Aldaz Mar 27 at 6:13
10
$\begingroup$

RDKit

This is pretty easy to do in RDKit. If you want the molecular formula, you can just use CalcMolFormula():

from rdkit import Chem
from rdkit.Chem.rdMolDescriptors import CalcMolFormula

# e.g. cysteine
mol = Chem.MolFromSmiles("C([C@@H](C(=O)O)N)S")
formula = CalcMolFormula(mol)

It is also very easy to get all of the unique heteroatoms:

from rdkit import Chem

# e.g. cysteine
mol = Chem.MolFromSmiles("C([C@@H](C(=O)O)N)S")
hetatms = {atom.GetSymbol() for atom in mol.GetAtoms()}
$\endgroup$
9
$\begingroup$

Pybel

Pybel is good for this type of task https://open-babel.readthedocs.io/en/latest/UseTheLibrary/Python_PybelAPI.html

import pybel
global __ATOM_LIST__
 __ATOM_LIST__ = \ 
      ['h',  'he',                                                                                                                                                                                                                                                 
      'li', 'be', 'b',  'c',  'n',  'o',  'f',  'ne',
      'na', 'mg', 'al', 'si', 'p',  's',  'cl', 'ar',
      'k',  'ca', 'sc', 'ti', 'v ', 'cr', 'mn', 'fe', 'co', 'ni', 'cu',
      'zn', 'ga', 'ge', 'as', 'se', 'br', 'kr',
      'rb', 'sr', 'y',  'zr', 'nb', 'mo', 'tc', 'ru', 'rh', 'pd', 'ag',
      'cd', 'in', 'sn', 'sb', 'te', 'i',  'xe',
      'cs', 'ba', 'la', 'ce', 'pr', 'nd', 'pm', 'sm', 'eu', 'gd', 'tb', 'dy',
      'ho', 'er', 'tm', 'yb', 'lu', 'hf', 'ta', 'w',  're', 'os', 'ir', 'pt',
      'au', 'hg', 'tl', 'pb', 'bi', 'po', 'at', 'rn',
      'fr', 'ra', 'ac', 'th', 'pa', 'u',  'np', 'pu']

mol = pybel.readstring("smi", smi)
numbers = [atom.atomicnum for atom in mol]

def from_number_to_element(i):
    return __ATOM_LIST__[i-1].capitalize()

names= [ from_number_to_element(i) for i in numbers]
    
$\endgroup$
4
  • 2
    $\begingroup$ Oof, looks good but pybel doesn't have readstring... turns out I need to install pip install openbabel not pip install pybel...turns out openbabel doesn't want to install, turns out it may be python2? python... so much ductape :) will try in the morning, need sleep $\endgroup$ – B. Kelly Mar 27 at 7:14
  • $\begingroup$ You want from openbabel import pybel for OB v3. $\endgroup$ – Geoff Hutchison Mar 28 at 17:32
  • 1
    $\begingroup$ As far as the "other pybel" - there's not much that can be done. I use conda, so conda install -c conda-forge openbabel $\endgroup$ – Geoff Hutchison Mar 28 at 17:34
  • $\begingroup$ @GeoffHutchison yeah it seems a bit tragic, I found myself at the github discussion of it. $\endgroup$ – B. Kelly Mar 30 at 22:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.