[Rdkit-discuss] gromacs force field and topology from SMILES

2016-10-14 Thread Daniel Seeliger
Hi Rafal,

as far as I know the MMFF is not compatible with any contemporary MD
forcefield. The usual way to generate topologies for ligands is the
generalised amber forcefield (GAFF). What you need are ambertools (
http://ambermd.org/) and acpype (https://github.com/t-/acpype). As input
you need a mol2 file of your ligand which you can generate with the rdkit
from smiles. Make sure to add hydrogens and the protonation state you want.
Partial charges can either be computed with AM1-BCC or Gaussian (RESP).
Acpype can generate topologies for Amber or Gromacs. In either case you
should use an amber forcefield for the protein (I guess you want to run a
simulation of a complex).

I recommend to check the structure after minimisation, especially if you’re
attempting to generate a topology for a molecule which is not common (e.g.
a transition state) since GAFF sometimes assigns incorrect angle and
dihedral parameters.

Best,
Daniel


On 14 Oct 2016, at 15:12, Paolo Tosco > wrote:

Dear Rafal,

answers below.

On 10/14/16 13:53, Rafal Roszak wrote:

Dear all,

1. I want to generate gromacs topology with force field parameters
from SMILES. RDKit has two FFs so this should be easy task at least for
this two force fields. Bond and angle parameters can be simply take
from Get*BondStretchParams and Get*StretchBendParams, but vdW and
torsian parameters are not so obvious.
-) How to assign vdW params to atoms? How can I extract this from rdkit?

The RDKit currently returns has two Get*VdWParams() functions which
retrieve the vdW parameters for a specific *pair* of atoms, i.e. after
the {MMFF, UFF} combination rules have been applied. I can add getters
for the raw per-atom vdW parameters; I'll let you know once they are
implemented.

-) How to convert result of Get*TorsionParams (Inversion and OopBend) to
something which can be used by gromacs? (I know, it is not strictly
question about RDKit)

I am not myself a Gromacs user and I don't know about the Gromacs
topology format, so I am afraid I can't help with this one.


2. It think I found bug in MMFF implementation:

acetona=Chem.AddHs(Chem.MolFromSmiles('[CH2-]C(=O)C') )
AllChem.MMFFHasAllMoleculeParams(acetona)

True

acetmm=AllChem.MMFFGetMoleculeProperties(acetona)
sumch=0
for i in range( len(acetona.GetAtoms()) ):

...  sumch+=acetmm.GetMMFFPartialCharge(i)
...

round(sumch,2)

0.0

for i in range( len(acetona.GetAtoms()) ):

...  sumch+=acetmm.GetMMFFFormalCharge(i)
...

round(sumch,2)

0.0

The molecule is anion but sum of formal or partial charge is zero. In
my opinion it is not currect result :) The same behaviour I also
observed for other exotic C-anion.

The MMFF force field does not have parameters for carbanions nor for
carbocations, so AllChem.MMFFHasAllMoleculeParams() should return False.
Thanks for reporting this; I'll file a bug and fix it as soon as I have
a bit of time.

Best,
Paolo

Regards,

Rafał


--
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
--
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


Re: [Rdkit-discuss] gromacs force field and topology from SMILES

2016-10-14 Thread Daniel Seeliger
Hi Rafal,

as far as I know the MMFF is not compatible with any contemporary MD 
forcefield. The usual way to generate topologies for ligands is the generalised 
amber forcefield (GAFF). What you need are ambertools (http://ambermd.org/ 
) and acpype (https://github.com/t-/acpype 
). As input you need a mol2 file of your ligand 
which you can generate with the rdkit from smiles. Make sure to add hydrogens 
and the protonation state you want. Partial charges can either be computed with 
AM1-BCC or Gaussian (RESP). 
Acpype can generate topologies for Amber or Gromacs. In either case you should 
use an amber forcefield for the protein (I guess you want to run a simulation 
of a complex).

I recommend to check the structure after minimisation, especially if you’re 
attempting to generate a topology for a molecule which is not common (e.g. a 
transition state) since GAFF sometimes assigns incorrect angle and dihedral 
parameters.

Best,
Daniel


> On 14 Oct 2016, at 15:12, Paolo Tosco  wrote:
> 
> Dear Rafal,
> 
> answers below.
> 
> On 10/14/16 13:53, Rafal Roszak wrote:
>> Dear all,
>> 
>> 1. I want to generate gromacs topology with force field parameters
>> from SMILES. RDKit has two FFs so this should be easy task at least for
>> this two force fields. Bond and angle parameters can be simply take
>> from Get*BondStretchParams and Get*StretchBendParams, but vdW and
>> torsian parameters are not so obvious.
>> -) How to assign vdW params to atoms? How can I extract this from rdkit?
> The RDKit currently returns has two Get*VdWParams() functions which 
> retrieve the vdW parameters for a specific *pair* of atoms, i.e. after 
> the {MMFF, UFF} combination rules have been applied. I can add getters 
> for the raw per-atom vdW parameters; I'll let you know once they are 
> implemented.
>> -) How to convert result of Get*TorsionParams (Inversion and OopBend) to
>> something which can be used by gromacs? (I know, it is not strictly
>> question about RDKit)
> I am not myself a Gromacs user and I don't know about the Gromacs 
> topology format, so I am afraid I can't help with this one.
>> 
>> 2. It think I found bug in MMFF implementation:
>> 
> acetona=Chem.AddHs(Chem.MolFromSmiles('[CH2-]C(=O)C') )
> AllChem.MMFFHasAllMoleculeParams(acetona)
>> True
> acetmm=AllChem.MMFFGetMoleculeProperties(acetona)
> sumch=0
> for i in range( len(acetona.GetAtoms()) ):
>> ...  sumch+=acetmm.GetMMFFPartialCharge(i)
>> ...
> round(sumch,2)
>> 0.0
> for i in range( len(acetona.GetAtoms()) ):
>> ...  sumch+=acetmm.GetMMFFFormalCharge(i)
>> ...
> round(sumch,2)
>> 0.0
>> 
>> The molecule is anion but sum of formal or partial charge is zero. In
>> my opinion it is not currect result :) The same behaviour I also
>> observed for other exotic C-anion.
>> 
> The MMFF force field does not have parameters for carbanions nor for 
> carbocations, so AllChem.MMFFHasAllMoleculeParams() should return False. 
> Thanks for reporting this; I'll file a bug and fix it as soon as I have 
> a bit of time.
> 
> Best,
> Paolo
>> Regards,
>> 
>> Rafał
>> 
>> --
>> 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 
> 
--
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


Re: [Rdkit-discuss] gromacs force field and topology from SMILES

2016-10-14 Thread Paolo Tosco
Dear Rafal,

answers below.

On 10/14/16 13:53, Rafal Roszak wrote:
> Dear all,
>
> 1. I want to generate gromacs topology with force field parameters
> from SMILES. RDKit has two FFs so this should be easy task at least for
> this two force fields. Bond and angle parameters can be simply take
> from Get*BondStretchParams and Get*StretchBendParams, but vdW and
> torsian parameters are not so obvious.
> -) How to assign vdW params to atoms? How can I extract this from rdkit?
The RDKit currently returns has two Get*VdWParams() functions which 
retrieve the vdW parameters for a specific *pair* of atoms, i.e. after 
the {MMFF, UFF} combination rules have been applied. I can add getters 
for the raw per-atom vdW parameters; I'll let you know once they are 
implemented.
> -) How to convert result of Get*TorsionParams (Inversion and OopBend) to
> something which can be used by gromacs? (I know, it is not strictly
> question about RDKit)
I am not myself a Gromacs user and I don't know about the Gromacs 
topology format, so I am afraid I can't help with this one.
>
> 2. It think I found bug in MMFF implementation:
>
 acetona=Chem.AddHs(Chem.MolFromSmiles('[CH2-]C(=O)C') )
 AllChem.MMFFHasAllMoleculeParams(acetona)
> True
 acetmm=AllChem.MMFFGetMoleculeProperties(acetona)
 sumch=0
 for i in range( len(acetona.GetAtoms()) ):
> ...  sumch+=acetmm.GetMMFFPartialCharge(i)
> ...
 round(sumch,2)
> 0.0
 for i in range( len(acetona.GetAtoms()) ):
> ...  sumch+=acetmm.GetMMFFFormalCharge(i)
> ...
 round(sumch,2)
> 0.0
>
> The molecule is anion but sum of formal or partial charge is zero. In
> my opinion it is not currect result :) The same behaviour I also
> observed for other exotic C-anion.
>
The MMFF force field does not have parameters for carbanions nor for 
carbocations, so AllChem.MMFFHasAllMoleculeParams() should return False. 
Thanks for reporting this; I'll file a bug and fix it as soon as I have 
a bit of time.

Best,
Paolo
> Regards,
>
> Rafał
>
> --
> 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] gromacs force field and topology from SMILES

2016-10-14 Thread Rafal Roszak
Dear all,

1. I want to generate gromacs topology with force field parameters
from SMILES. RDKit has two FFs so this should be easy task at least for
this two force fields. Bond and angle parameters can be simply take
from Get*BondStretchParams and Get*StretchBendParams, but vdW and
torsian parameters are not so obvious.
-) How to assign vdW params to atoms? How can I extract this from rdkit?
-) How to convert result of Get*TorsionParams (Inversion and OopBend) to
something which can be used by gromacs? (I know, it is not strictly
question about RDKit)

2. It think I found bug in MMFF implementation:

>>> acetona=Chem.AddHs(Chem.MolFromSmiles('[CH2-]C(=O)C') )
>>> AllChem.MMFFHasAllMoleculeParams(acetona)
True
>>> acetmm=AllChem.MMFFGetMoleculeProperties(acetona)
>>> sumch=0
>>> for i in range( len(acetona.GetAtoms()) ):
...  sumch+=acetmm.GetMMFFPartialCharge(i)
... 
>>> round(sumch,2)
0.0
>>> for i in range( len(acetona.GetAtoms()) ):
...  sumch+=acetmm.GetMMFFFormalCharge(i)
... 
>>> round(sumch,2)
0.0

The molecule is anion but sum of formal or partial charge is zero. In
my opinion it is not currect result :) The same behaviour I also
observed for other exotic C-anion.


Regards,

Rafał

--
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