6
$\begingroup$

For a substructure search I would like to search for structures containing unfused benzyl. The idea was to explicitly add hydrogen. But apparently this does not give the expected results.

So, I assume either my understanding or expectations are wrong or I'm using RDKit not properly. I just learnt (Substructure search with RDKit) that I better use SMARTS expressions with Chem.MolToSmarts() where I get the expected results.

However, if I have much more complex search structures, how can I get from the SMILES to the SMARTS? Do I have always to figure out the SMARTS code myself? Is there maybe a way just drawing a molecule and add explicit hydrogen and search with this somehow? I guess it should be pretty clear (at least to a human) what I want to find (or not find) with structure 2. Are there maybe any tools or code which can assist here?

enter image description here

Code:

from rdkit import Chem

smiles_strings = '''CC1=CC=CC=C1
CC1=C([H])C([H])=C([H])C([H])=C1[H]
C1(C2=CC=CC=C2)=CC=CC=C1
C12=C(CC3=CC=CC=C3C2)C=CC=C1
C12=CC=CC=C1C=C3C(C=CC=C3)=C2
C12=C(C(C(C=C3)=CC=C4)=C4C=C2)C3=CC=C1
'''
smiles_list = smiles_strings.splitlines()

def search_structure(pattern):
    found = []
    for idx,smiles in enumerate(smiles_list):
        m = Chem.MolFromSmiles(smiles)
        if m.HasSubstructMatch(pattern):
            found.append(idx+1)
    print("Structures found: {}".format(found))

smiles_1  = smiles_list[0]
pattern_1 = Chem.MolFromSmiles(smiles_1)
smiles_2a  = smiles_list[1]
pattern_2a = Chem.MolFromSmiles(smiles_2a)
smiles_2b  = Chem.MolToSmiles(pattern_2a)
pattern_2b = Chem.MolFromSmiles(smiles_2b)
smiles_2c = smiles_2b
pattern_2c = Chem.MolFromSmarts(smiles_2c)
smarts_3 = '[cR1]1[cR1][cR1][cR1][cR1][cR1]1-[c,C]'
pattern_3 = Chem.MolFromSmarts(smarts_3)

print("\nSMILES 1 : {}\n# my expectation: [1, 2, 3, 4, 5, 6]".format(smiles_1))
search_structure(pattern_1)
print("\nSMILES 2a: {}\n# my expectation: [1, 2, 3]".format(smiles_2a))
search_structure(pattern_2a)
print("\nSMILES 2b: {} (Mol from SMILES)\n# my expectation: [1, 2, 3]".format(smiles_2b))
search_structure(pattern_2b)
print("\nSMILES 2c: {} (Mol from SMARTS)\n# my expectation: [1, 2, 3]".format(smiles_2c))
search_structure(pattern_2c)
print("\nSMARTS 3 : {}\n# my expectation: [1, 2, 3]".format(smarts_3))
search_structure(pattern_3)

Result:

SMILES 1 : CC1=CC=CC=C1
# my expectation: [1, 2, 3, 4, 5, 6]
Structures found: [1, 2, 3, 4]

SMILES 2a: CC1=C([H])C([H])=C([H])C([H])=C1[H]
# my expectation: [1, 2, 3]
Structures found: [1, 2, 3, 4]

SMILES 2b: Cc1ccccc1 (Mol from SMILES)
# my expectation: [1, 2, 3]
Structures found: [1, 2, 3, 4]

SMILES 2c: Cc1ccccc1 (Mol from SMARTS)
# my expectation: [1, 2, 3]
Structures found: [1, 2, 4]

SMARTS 3 : [cR1]1[cR1][cR1][cR1][cR1][cR1]1-[c,C]
# my expectation: [1, 2, 3]
Structures found: [1, 2, 3]
$\endgroup$
5
$\begingroup$

If you work in a Jupyter Notebook you can visualize your substructure patterns and the search results, so you can see how the patterns work.

from rdkit import Chem
from rdkit.Chem import Draw

smiles = ['CC1=CC=CC=C1','CC1=C([H])C([H])=C([H])C([H])=C1[H]',
          'C1(C2=CC=CC=C2)=CC=CC=C1','C12=C(CC3=CC=CC=C3C2)C=CC=C1',
          'C12=CC=CC=C1C=C3C(C=CC=C3)=C2','C12=C(C(C(C=C3)=CC=C4)=C4C=C2)C3=CC=C1']

params = Chem.SmilesParserParams()
params.removeHs=False # draw and work with explicit Hs
mols = [Chem.MolFromSmiles(s, params) for s in smiles]
Draw.MolsToGridImage(mols, molsPerRow=3)

This will give you the compounds and especially the one with the explicit Hs.

Compounds

Now we get the patterns and especially the one with the explicit Hs.

patt1 = Chem.MolFromSmiles('CC1=CC=CC=C1')
patt2 = Chem.MolFromSmiles('CC1=C([H])C([H])=C([H])C([H])=C1[H]', params) # search with explicit Hs
patt3 = Chem.MolFromSmiles('Cc1ccccc1')
patt4 = Chem.MolFromSmarts('Cc1ccccc1')
patt5 = Chem.MolFromSmarts('[cR1]1[cR1][cR1][cR1][cR1][cR1]1-[c,C]')

patts = [patt1,patt2,patt3,patt4,patt5]
leg = ['CC1=CC=CC=C1','CC1=C([H])C([H])=C([H])C([H])=C1[H]','Cc1ccccc1 smiles',
       'Cc1ccccc1 smarts','[cR1]1[cR1][cR1][cR1][cR1][cR1]1-[c,C]']
Draw.MolsToGridImage(patts, molsPerRow=3, legends=leg)

patterns

Search and highlight the result.

ms = []
patt = []
allsubs = []
for m in range(len(mols)):
    for p in range(len(patts)):
        sub = mols[m].GetSubstructMatches(patts[p])
        if len(sub) > 0:
            ms.append(mols[m])
            patt.append(leg[p])
            allsubs.append(sub[0]) # substructures could be find multiple time - just take the first

Draw.MolsToGridImage(ms, molsPerRow=3, legends=patt, highlightAtomLists=allsubs)

Compounds 5 and 6 have no nonaromatic bonds, so there is nothing to find.

Results

$\endgroup$
5
  • $\begingroup$ thank you very much for your code and explanations. It looks like params.removeHs=False is the bit I was missing. I can only test next week and will let you know. $\endgroup$ – theozh Mar 7 '20 at 5:21
  • $\begingroup$ If I'm using params.removeHs=False for searching with structure 2, then HasSubstructMatch(pattern) only finds structure 2 (as in your example). Why are structures 1 and 3 not found as well? Maybe because of the explicit non-aromatic extra carbon in search structure 2? What do I have to change that by drawing structure 2 with explicit hydrogen and unspecified (is it or not?) extra carbon will also find structures 1 and 3? $\endgroup$ – theozh Mar 10 '20 at 6:45
  • $\begingroup$ @theozh It is not the extra carbon, it is just that 1 and 3 do not have explicit Hs. If I add explicit Hs to 3 it is found. I thought C[c;H,h]1[c;H,h][c;H,h][c;H,h][c;H,h][c;H,h]1 (see smartsview.zbh.uni-hamburg.de) will work, but it doesn't. $\endgroup$ – rapelpy Mar 10 '20 at 18:58
  • $\begingroup$ thank you, yes, today I also found out that if you add Hs to all the structures you're searching in, via mol = Chem.AddHs(mol), you will find what is expected. However, there are other cases which I still do not understand. I need more testing. $\endgroup$ – theozh Mar 10 '20 at 19:16
  • $\begingroup$ It took me a while, but you can find 1, 2 and 3 with *[c]1[c;H,h][c;H,h][c;H,h][c;H,h][c;H,h]1 and 1,2,3,4,5 with *[c]1[c;H,h,*][c;H,h][c;H,h][c;H,h][c;H,h]1 $\endgroup$ – rapelpy Mar 12 '20 at 19:48

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.