Re: [Rdkit-discuss] Filtering Out Highly Strained Ring Systems

2017-04-17 Thread Greg Landrum
Hi James,

There's not currently any way to get the contributions to the overall
energy from a force field.
Here's one approach to approximate what you're trying to do:

from rdkit import Chem
from rdkit.Chem import AllChem

def fragmentAndGetEnergies(mol,sma='[r]-[!r;!#1]'):
  pattern = Chem.MolFromSmarts(sma)
  m.GetSubstructMatches(pattern)
  bonds = [m.GetBondBetweenAtoms(x[0],x[1]).GetIdx() for x in matches]
  splitmol = AllChem.FragmentOnBonds(m,bonds)
  frags = Chem.GetMolFrags(splitmol,asMols=True)
  res = []
  for frag in frags:
   ff = AllChem.UFFGetMoleculeForceField(frag)
   res.append((Chem.MolToSmiles(frag,True),ff.CalcEnergy()))
  return res

m = Chem.AddHs(Chem.MolFromSmiles('c1cccn1CCC1=CC1'))
AllChem.EmbedMolecule(m,AllChem.ETKDG())
vs = fragmentAndGetEnergies(m)
for smi,energy in vs:
print(smi,energy)


The function splits molecules into pieces by breaking single bonds between
rings and non-ring atoms and then calculates the energies of those pieces.
This generates warnings from the UFF code since it doesn't have parameters
for the dummy atoms that FragmentOnBonds() inserts. I think this is
probably ok, but you may want to investigate changing those dummies to
hydrogens (and adjusting the corresponding bond lengths).

Best,
-greg



On Sat, Apr 15, 2017 at 5:11 PM, James Johnson 
wrote:

> Dear All,
>
> I'm generating analogs from a parent compound. An unfortunate side effect
> is the generation of various highly strained rings among other compounds.
> I'm developing filters to deal with any unwanted generated compounds, but
> I'm having difficulty with strained ring systems.
>
> Attached are some examples in mol2 format: https://ufile.io/bdmcp
>
> For this problem I'm thinking of using an energy/bond approach:
> ff = AllChem.UFFGetMoleculeForceField(mol)
> energy_per_bond = ff.CalcEnergy()/len(mol.GetBonds())
>
> If energy/bond is above a certain value it gets filtered out. Do you think
> this is the best solution? Is there a way of viewing each bonds energy
> (with increasing amount of bonds an unstable bond might not be filtered)?
> Do you have any suggestions on alternate ways of doing this?
>
> Thank you for the help, it is very much appreciated! :)
>
> - James
>
> 
> --
> Check out the vibrant tech community on one of the world's most
> engaging tech sites, Slashdot.org! http://sdm.link/slashdot
> ___
> Rdkit-discuss mailing list
> Rdkit-discuss@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
>
>
--
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot___
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss


[Rdkit-discuss] Filtering Out Highly Strained Ring Systems

2017-04-15 Thread James Johnson
Dear All,

I'm generating analogs from a parent compound. An unfortunate side effect
is the generation of various highly strained rings among other compounds.
I'm developing filters to deal with any unwanted generated compounds, but
I'm having difficulty with strained ring systems.

Attached are some examples in mol2 format: https://ufile.io/bdmcp

For this problem I'm thinking of using an energy/bond approach:
ff = AllChem.UFFGetMoleculeForceField(mol)
energy_per_bond = ff.CalcEnergy()/len(mol.GetBonds())

If energy/bond is above a certain value it gets filtered out. Do you think
this is the best solution? Is there a way of viewing each bonds energy
(with increasing amount of bonds an unstable bond might not be filtered)?
Do you have any suggestions on alternate ways of doing this?

Thank you for the help, it is very much appreciated! :)

- James
--
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot___
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss