Re: [Rdkit-discuss] Conformers and Energies

2019-02-14 Thread Greg Landrum
Hi,

ETKDG is the default conformation generator in versions 2018.09.1 and
later, so you don't need to worry about specifying that yourself.

-greg



On Thu, Feb 14, 2019 at 11:52 AM Stamatia Zavitsanou <
stamatia.zavitsa...@oriel.ox.ac.uk> wrote:

> Thank you greg,
>
>
> That was very helpful.
>
>
> I will follow your way. I am just interested in getting rough energies. So
> I am not really sure if should minimise or not(?)(I was using it thinking
> it was a necessary step for DG conformers).
>
>
> Also no I don't want to remove the Hs from my molecule, so I can totally
> ignore the MolFromMolBlock step.
>
>
> Also what if I use the ETKDG :
>
>
> mol2 = Chem.MolFromSmiles(line2)
>
> m = Chem.AddHs(mol2)
>
>
> cids=AllChem.EmbedMultipleConfs(m,numConfs,Chem.rdDistGeom.ETKDG())
>
> ps = AllChem.MMFFGetMoleculeProperties(m)
>
> for k in cids:
>
>  ff = AllChem.MMFFGetMoleculeForceField(m,ps,confId=cid)
>
>  ff.Initialize()
>
>  print(ff.CalcEnergy())
>
>
> would that be correct?
>
>
> Thank you,
>
> Stamatia Zavitsanou
> --
> *From:* Greg Landrum 
> *Sent:* 14 February 2019 10:32:23
> *To:* Stamatia Zavitsanou
> *Cc:* rdkit-discuss@lists.sourceforge.net
> *Subject:* Re: [Rdkit-discuss] Conformers and Energies
>
> Hi Stamatia,
>
> Assuming that you are interested in getting the energies of non-minimized
> conformers, here's a relatively easy way to do what you want. Notice that
> I'm using the MMFF94 force field, since it's probably a better choice for
> energies. You could do the same thing with UFF if you'd prefer.
>
> In [4]: m = Chem.AddHs(Chem.MolFromSmiles('OCCCN'))
>
> In [5]: cids = AllChem.EmbedMultipleConfs(m,10)
>
> In [6]: ps = AllChem.MMFFGetMoleculeProperties(m)
>
> In [8]: for cid in cids:
>...: ff = AllChem.MMFFGetMoleculeForceField(m,ps,confId=cid)
>...: ff.Initialize()
>...: print(ff.CalcEnergy())
>...:
>...:
> 12.233357032097315
> 4.7643146024595335
> 6.163107627780452
> 5.476964847773123
> 4.412833444203187
> 5.665311226271907
> 7.536117466355753
> 3.421494301728207
> 8.97909661823493
> 10.863891606189476
>
>
> Your first approach is inefficient and just a bit wrong: the call to
> MolFromMolBlock() will remove the Hs from your molecule. If you want to go
> that route anyway, you should do:
> b.append(Chem.MolFromMolBlock(a,removeHs=False))
>
> Note that if you also want to minimize the conformers, you can do this:
>
> In [9]: AllChem.MMFFOptimizeMoleculeConfs(m)
> Out[9]:
> [(0, 1.0966514168737902),
>  (0, -1.512072483231206),
>  (0, 0.6847373872800766),
>  (0, -0.6067431133406469),
>  (0, -1.5142123177302143),
>  (0, -2.08974298611677),
>  (0, -2.089742985761273),
>  (0, -2.053178209868828),
>  (0, -3.9515977883892077),
>  (0, -1.0767726284668093)]
>
> There you get a 2-tuple for each conformer. The first element lets you
> know whether or not the minimization converged, and the second is the
> energy of the conformer.
>
> I hope this helps,
> -greg
>
>
> On Tue, Feb 12, 2019 at 5:05 AM Stamatia Zavitsanou <
> stamatia.zavitsa...@oriel.ox.ac.uk> wrote:
>
> Hello everyone,
>
>
> I am a bit confused. I want to make conformers and calculate their
> energies, but I am not very sure which is the correct way and at which
> point I need to embed my molecule.
>
>
>
> I am doing this :
>
>
> for line2 in infile2:
>
>
> mol2 = Chem.MolFromSmiles(line2)
>
> m2 = Chem.AddHs(mol2)
>
>
> cids = AllChem.EmbedMultipleConfs(m2, numConfs)
>
> for k in cids:
>
> a = Chem.MolToMolBlock(m2, confId=k)
>
> b.append(Chem.MolFromMolBlock(a))
>
>
> fff = AllChem.UFFGetMoleculeForceField(b[k])
>
> fff.Initialize()
>
> fff.Minimize()
>
>
> energy_value_U = fff.CalcEnergy()
>
> print("for each conformer")
>
> print(energy_value_U)
>
> energies.append(float(energy_value_U))
>
>
> Should I be doing something else? like:
>
>
> for line2 in infile2:
>
>
> mol2 = Chem.MolFromSmiles(line2)
>
> m2 = Chem.AddHs(mol2)
>
>
> cids = AllChem.EmbedMultipleConfs(m2, numConfs)
>
> for k in cids:
>
> a = Chem.MolToMolBlock(m2, confId=k)
>
> b.append(Chem.MolFromMolBlock(a))
>
>
> AllChem.EmbedMolecule(b[k])
>
> fff = AllChem.UFFGetMoleculeForceField(b[k])
>
> fff.Initialize()
>
> fff.Minimize()
>
> energy_value_U = fff.CalcEnergy()
>
> print("for each conformer")
>
> print(energy_value_U)
>
> energies.append(float(energy_value_U))
>
>
> or like this:
>
>
> for line2 in infile2:
>
> print(line2)
>
> print(i)
>
> mol2 = Chem.MolFromSmiles(line2)
>
> m2 = Chem.AddHs(mol2)
>
> AllChem.EmbedMolecule(m2)
>
> cids = AllChem.EmbedMultipleConfs(m2, numConfs)
>
> for k in cids:
>
> a = 

Re: [Rdkit-discuss] Scaffold Tree implementation

2019-02-14 Thread Francois Berenger

On 14/02/2019 19:44, Colin Bournez wrote:

Dear all,

I would like to know if there is the possibility to use the Scaffold
Tree algorithm in RDKit?
In the documentation, there is an existing page :
http://rdkit.org/docs/source/rdkit.Chem.Scaffolds.ScaffoldTree.html
[1]

But it is empty...


Maybe this open source software can do what you need:

http://scaffoldhunter.sourceforge.net/

Regards,
F.


I also saw message in rdkit list from 2015 but since there nothing.
Are there any news?

Cheers,

--
 Signature
*

 BOURNEZ COLIN [2]  [3]
Chemoinformatics PhD Student

  INSTITUTE OF ORGANIC AND ANALYTICAL CHEMISTRY (ICOA UMR7311)
 Université d'Orléans - Pôle de Chimie
 Rue de Chartres - BP 6759
 45067 Orléans Cedex 2 - France
 +33 (0)2 38 49 45 77 [4]
 SBC Tool Platform [5] - SBC Team [6] [7]

 [8]

  [9]

 [10]



Links:
--
[1] http://rdkit.org/docs/source/rdkit.Chem.Scaffolds.ScaffoldTree.html
[2] https://www.linkedin.com/in/colin-bournez-b9a1b2b7/
[3] http://www.icoa.fr/
[4] tel:+33%202%2038%2049%2045%2077
[5] http://sbc.icoa.fr/
[6] http://www.icoa.fr/bonnet
[7] http://www.icoa.fr/fr/rss.xml
[8]
https://www.facebook.com/pages/Institut-de-Chimie-Organique-et-Analytique-ICOA-umr7311/222060911297163
[9] https://twitter.com/ICOA_UMR7311
[10]
https://www.linkedin.com/company/institut-de-chimie-organique-et-analytique---icoa-umr7311/

___
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss




___
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss


Re: [Rdkit-discuss] Scaffold Tree implementation

2019-02-14 Thread Greg Landrum
Hi Colin,

It would be great to have one, but there is, to the best of my knowledge,
no scaffold tree implementation currently available for the RDKit.

-greg


On Thu, Feb 14, 2019 at 11:45 AM Colin Bournez <
colin.bour...@univ-orleans.fr> wrote:

> Dear all,
>
> I would like to know if there is the possibility to use the Scaffold Tree
> algorithm in RDKit?
> In the documentation, there is an existing page :
> http://rdkit.org/docs/source/rdkit.Chem.Scaffolds.ScaffoldTree.html
>
> But it is empty...
>
> I also saw message in rdkit list from 2015 but since there nothing. Are
> there any news?
>
> Cheers,
>
> --
> *Bournez Colin *  
> 
>
> *Chemoinformatics PhD Student *
> * Institute of Organic and Analytical Chemistry (ICOA UMR7311)*
>  Université d'Orléans - Pôle de Chimie
>  Rue de Chartres - BP 6759
>  45067 Orléans Cedex 2 - France
>  +33 (0)2 38 49 45 77 <+33%202%2038%2049%2045%2077>
>  SBC Tool Platform  - SBC Team
>   
>
> 
>  
>
> 
>
>
> ___
> Rdkit-discuss mailing list
> Rdkit-discuss@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
>
___
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss


[Rdkit-discuss] Scaffold Tree implementation

2019-02-14 Thread Colin Bournez

Dear all,

I would like to know if there is the possibility to use the Scaffold 
Tree algorithm in RDKit?

In the documentation, there is an existing page :
http://rdkit.org/docs/source/rdkit.Chem.Scaffolds.ScaffoldTree.html

But it is empty...

I also saw message in rdkit list from 2015 but since there nothing. Are 
there any news?


Cheers,

--
Signature
*Bournez Colin *  
 

/Chemoinformatics PhD Student
/
* Institute of Organic and Analytical Chemistry (ICOA UMR7311)*
 Université d'Orléans - Pôle de Chimie
 Rue de Chartres - BP 6759
 45067 Orléans Cedex 2 - France
 +33 (0)2 38 49 45 77 
SBC Tool Platform  - SBC Team 
 	 
 


 
 






___
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss


Re: [Rdkit-discuss] Conformers and Energies

2019-02-14 Thread Greg Landrum
Hi Stamatia,

Assuming that you are interested in getting the energies of non-minimized
conformers, here's a relatively easy way to do what you want. Notice that
I'm using the MMFF94 force field, since it's probably a better choice for
energies. You could do the same thing with UFF if you'd prefer.

In [4]: m = Chem.AddHs(Chem.MolFromSmiles('OCCCN'))

In [5]: cids = AllChem.EmbedMultipleConfs(m,10)

In [6]: ps = AllChem.MMFFGetMoleculeProperties(m)

In [8]: for cid in cids:
   ...: ff = AllChem.MMFFGetMoleculeForceField(m,ps,confId=cid)
   ...: ff.Initialize()
   ...: print(ff.CalcEnergy())
   ...:
   ...:
12.233357032097315
4.7643146024595335
6.163107627780452
5.476964847773123
4.412833444203187
5.665311226271907
7.536117466355753
3.421494301728207
8.97909661823493
10.863891606189476


Your first approach is inefficient and just a bit wrong: the call to
MolFromMolBlock() will remove the Hs from your molecule. If you want to go
that route anyway, you should do:
b.append(Chem.MolFromMolBlock(a,removeHs=False))

Note that if you also want to minimize the conformers, you can do this:

In [9]: AllChem.MMFFOptimizeMoleculeConfs(m)
Out[9]:
[(0, 1.0966514168737902),
 (0, -1.512072483231206),
 (0, 0.6847373872800766),
 (0, -0.6067431133406469),
 (0, -1.5142123177302143),
 (0, -2.08974298611677),
 (0, -2.089742985761273),
 (0, -2.053178209868828),
 (0, -3.9515977883892077),
 (0, -1.0767726284668093)]

There you get a 2-tuple for each conformer. The first element lets you know
whether or not the minimization converged, and the second is the energy of
the conformer.

I hope this helps,
-greg


On Tue, Feb 12, 2019 at 5:05 AM Stamatia Zavitsanou <
stamatia.zavitsa...@oriel.ox.ac.uk> wrote:

> Hello everyone,
>
>
> I am a bit confused. I want to make conformers and calculate their
> energies, but I am not very sure which is the correct way and at which
> point I need to embed my molecule.
>
>
>
> I am doing this :
>
>
> for line2 in infile2:
>
>
> mol2 = Chem.MolFromSmiles(line2)
>
> m2 = Chem.AddHs(mol2)
>
>
> cids = AllChem.EmbedMultipleConfs(m2, numConfs)
>
> for k in cids:
>
> a = Chem.MolToMolBlock(m2, confId=k)
>
> b.append(Chem.MolFromMolBlock(a))
>
>
> fff = AllChem.UFFGetMoleculeForceField(b[k])
>
> fff.Initialize()
>
> fff.Minimize()
>
>
> energy_value_U = fff.CalcEnergy()
>
> print("for each conformer")
>
> print(energy_value_U)
>
> energies.append(float(energy_value_U))
>
>
> Should I be doing something else? like:
>
>
> for line2 in infile2:
>
>
> mol2 = Chem.MolFromSmiles(line2)
>
> m2 = Chem.AddHs(mol2)
>
>
> cids = AllChem.EmbedMultipleConfs(m2, numConfs)
>
> for k in cids:
>
> a = Chem.MolToMolBlock(m2, confId=k)
>
> b.append(Chem.MolFromMolBlock(a))
>
>
> AllChem.EmbedMolecule(b[k])
>
> fff = AllChem.UFFGetMoleculeForceField(b[k])
>
> fff.Initialize()
>
> fff.Minimize()
>
> energy_value_U = fff.CalcEnergy()
>
> print("for each conformer")
>
> print(energy_value_U)
>
> energies.append(float(energy_value_U))
>
>
> or like this:
>
>
> for line2 in infile2:
>
> print(line2)
>
> print(i)
>
> mol2 = Chem.MolFromSmiles(line2)
>
> m2 = Chem.AddHs(mol2)
>
> AllChem.EmbedMolecule(m2)
>
> cids = AllChem.EmbedMultipleConfs(m2, numConfs)
>
> for k in cids:
>
> a = Chem.MolToMolBlock(m2, confId=k)
>
> b.append(Chem.MolFromMolBlock(a))
>
>
> AllChem.EmbedMolecule(b[k])
>
> fff = AllChem.UFFGetMoleculeForceField(b[k])
>
> fff.Initialize()
>
> fff.Minimize()
>
> energy_value_U = fff.CalcEnergy()
>
> print("for each conformer")
>
> print(energy_value_U)
>
> energies.append(float(energy_value_U))
>
>
> or maybe non of the above??
>
>
> Thanks a lot,
>
> Stamatia Zavitsanou
> ___
> Rdkit-discuss mailing list
> Rdkit-discuss@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
>
___
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss


Re: [Rdkit-discuss] Double bond in ring

2019-02-14 Thread Greg Landrum
Hi Jean-Marc,

This seems to be a bug in the SMILES parser. Here's a more compact
demonstration of the problem:

In [2]: inchi =
'InChI=1/C20H32/c1-16(2)20-14-12-18(4)10-6-8-17(3)9-7-11-19(5)13-15-20/h8,10-12,14,16,20H,6-7,9,13,15H2,1-5
   ...: H3/b14-12+,17-8+,18-10-,19-11-/t20-/m0/s1'
   ...: m1 = Chem.MolFromInchi(inchi)

In [3]: for bond in m1.GetBonds():
   ...: if bond.GetStereo() != Chem.BondStereo.STEREONONE:
   ...: print(bond.GetStereo())
STEREOE
STEREOZ
STEREOZ
STEREOE

In [4]: smi = Chem.MolToSmiles(m1)

In [5]: print(smi)
C/C1=C/C/C=C(\C)CC/C=C(/C)CC[C@@H](C(C)C)/C=C/1

In [6]: m2 = Chem.MolFromSmiles(smi)

In [8]: for bond in m2.GetBonds():
   ...: if bond.GetStereo() != Chem.BondStereo.STEREONONE:
   ...:
 
print(bond.GetIdx(),bond.GetStereo(),bond.GetBeginAtomIdx(),bond.GetEndAtomIdx())
4 STEREOE 4 5
9 STEREOZ 9 10
18 STEREOE 18 19


It may be that a fix for this is already in progress (as part of the work
I'm doing on this issue: https://github.com/rdkit/rdkit/issues/2148), but I
will do a bit of looking there.[1]

-greg
[1] and do the PR for that branch, which I seem to have forgotten


On Thu, Feb 14, 2019 at 6:28 AM Jean-Marc Nuzillard <
jm.nuzill...@univ-reims.fr> wrote:

> Dear all,
>
> looking at the attached file, you will see that the Mol -> SMILES  -> Mol
> conversions did not return the original molecule.
> Could there be a particular problem here either in SMILES generation
> or in SMILES reading?
> Reading complains with
> # [14:01:45] Conflicting single bond directions around double bond at
> index 1.
> # [14:01:45]   BondStereo set to STEREONONE and single bond directions set
> to NONE.
> However, as far as I can judge, smi1 looks good.
>
> Many thanks in advance,
>
> Jean-Marc
> --
>
> Dr. Jean-Marc Nuzillard
> Institute of Molecular Chemistry, CNRS UMR 7312
> Faculté des Sciences Exactes et Naturelles, Bâtiment 18
> BP 1039
> 51687 REIMS Cedex 2
> France
>
> Tel : 33 3 26 91 82 10
> Fax : 33 3 26 91 31 
> 66http://www.univ-reims.fr/ICMRhttp://eos.univ-reims.fr/LSD/CSNteam.html
> http://www.univ-reims.fr/LSD/http://www.univ-reims.fr/LSD/JmnSoft/
>
> ___
> Rdkit-discuss mailing list
> Rdkit-discuss@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/rdkit-discuss
>
___
Rdkit-discuss mailing list
Rdkit-discuss@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss