Standardize ionic SMILES / check if two ionic SMILES are equal including resonance structures

[simonmb]

Is there a way to standardize ionic SMILES? One example follows where the charge is on a different nitrogen in same the ring.

'CCn1cc[n+](C)c1'
'CC[n+]1ccn(C)c1'

In such a way that both are equal?

[simonmb]

I found out myself that by using the ResonanceMolSupplier the follwing test could be done:

def getResonanceStructureSMILES(SMILES, flags = 0):
    mol = Chem.MolFromSmiles(SMILES)
    suppl = Chem.ResonanceMolSupplier(mol, flags)
    return [Chem.MolToSmiles(mol_) for mol_ in suppl]

smiles1 = getResonanceStructureSMILES('CC[n+]1ccn(C)c1')
smiles2 = getResonanceStructureSMILES('CCn1cc[n+](C)c1')

set(smiles1) == set(smiles2) # returns True

Is this the proper way to test if two ionic SMILES with resonance structures are equal?

[simonmb]

If this does not interest anybody, can someone close this discussion?

[greglandrum]

Hi @simonmb, sorry sometimes I overlook things in the discussions tab.

If it's ok that different tautomers also end up grouped together, I would do this using the MolHash functionality and the tautomer normalization hash:

In [3]: m1 = Chem.MolFromSmiles('CC[n+]1ccn(C)c1')

In [4]: m2 = Chem.MolFromSmiles('CCn1cc[n+](C)c1')

In [5]: from rdkit.Chem import rdMolHash

In [7]: rdMolHash.MolHash(m1,rdMolHash.HashFunction.HetAtomTautomer)
Out[7]: 'CCN1[CH][CH]N(C)[CH]1_0_1'

In [8]: rdMolHash.MolHash(m2,rdMolHash.HashFunction.HetAtomTautomer)
Out[8]: 'CCN1[CH][CH]N(C)[CH]1_0_1'

Explanation of the options for the MolHash code can be found here:
https://nextmovesoftware.github.io/molhash/introduction.html

[simonmb]

Thanks! Unfortunately I cannot mark it as an answer.

[simonmb]

One followup question: Is the way I propose in the original question more restrictive then the one you propose? Or do you expect the output of both to be very similar?

[greglandrum]

Generating alternate tautomers (which moves Hs around along with bond orders) is different than looking at resonance structures (which just swaps bond orders, but can consider more bonds than the RDKit tautomer code).
I'm not sure that I would say that one is more restrictive than the other... they're just a bit different.

[simonmb]

Hi @greglandrum, I found something weird. When I run the following on windows 10, conda environment, rdkit version 2021.03.4:

m1 = Chem.MolFromSmiles('CC[n+]1ccn(C)c1')
suppl = Chem.ResonanceMolSupplier(m1, 0)
f = list(suppl)

execution of the program freezes on generating the list in the third line.

[greglandrum]

I can't reproduce that with the 2021.09.4 rdkit build (though I tried on linux). Can you try updating your RDKit install?

[simonmb]

I can confirm that updating to 2021.09.4 solves my problem.

[simonmb]

Sorry, @greglandrum, I found an example where also in 2021.09.4, Win 10, conda-forge the execution of the program freezes.

m1 = Chem.MolFromSmiles('Nc1ccc2nc3ccc(N)cc3[s+]c2c1')
suppl = Chem.ResonanceMolSupplier(m1, 0)
f = list(suppl)

[ptosco]

@simonmb This looks very much like #4884, I will try to find some time to look into this.

[simonmb]

Hi @greglandrum and @ptosco, since this issue #4491 has not been fixed yet, I am considering using the molhash solution instead. Is either one of the approaches a superset of the other? Meaning would resonance structures be a superset of tautomers, or the other way around, or would there be an intersection of the structures but some might be in the resonance structures set and others in the tautomers set?

My background is testing whether two ionic structures are the same with respect to quantum chemistry calculations, which don't care about where the charge is. But I don't want tautomers to be considered equal.