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

[Daniel S. <see...@go...>]

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 <pao...@un...
<javascript:_e(%7B%7D,'cvml','pao...@un...');>> 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://slashdot.org/>! http://
sdm.link/slashdot
_______________________________________________
Rdkit-discuss mailing list
Rdk...@li...
<javascript:_e(%7B%7D,'cvml','Rdk...@li...');>
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://slashdot.org/>! http://
sdm.link/slashdot
_______________________________________________
Rdkit-discuss mailing list
Rdk...@li...
<javascript:_e(%7B%7D,'cvml','Rdk...@li...');>
https://lists.sourceforge.net/lists/listinfo/rdkit-discuss