Re: [AMBER] Can Parmed convert AMBER files to GROMACS made with the AMBER19SB ? --update

From: ABEL Stephane <Stephane.ABEL.cea.fr>
Date: Thu, 18 Nov 2021 12:23:18 +0000

Dear all,

Below a little update of of my testing of AMBER19SB with GROMACS(2021.3). I have converted prmto/impcrd files (equivalent to prm7 and rst7) of a protein modeled with the AMBER19SB into GROMACS gro/top files using Parmed (AT21, fully patched) and ACEPYPE (last version).
The AMBER files were obtained using tleap (AT21). From the results below with the same gro file

I used the ParmEd commands below from Jason Swails

import parmed as pmd
parm = pmd.load_file("MYTOPFF19SBparm7", "MYTOPFF19SB.rst7")
parm.save("MYTOPFF19SB.top") # .top saves to GROMACS topology file format
parm.save("MYTOPFF19SB.gro")
```
And execute it as `python <file.py>`. Does that work for you?

The results with the same mdp of GROMACS at step0

>From Parmed (Parmed_Conversion_GROMACS_AMBER19SB.py, VERSION 3.4.1) at step0
           Step Time
              0 0.00000

   Energies (kJ/mol)
           Bond Angle Proper Dih. Improper Dih. CMAP Dih.
    1.52227e+03 7.26001e+03 1.15737e+04 2.23252e+02 3.32921e+03
          LJ-14 Coulomb-14 LJ (SR) Coulomb (SR) Potential
    1.02351e+04 9.98003e+04 -1.76602e+04 -1.87376e+05 -7.10925e+04
 Pressure (bar)
    1.77444e-01

>From ACEPYPE at step0
; 3ffh_ShuA_Alone_ff19SB_test_GMX.top created by acpype (v: 2021-02-05T22:15:50CET) on Mon Nov 15 13:06:39 2021

           Step Time
              0 0.00000
   Energies (kJ/mol)
           Bond Angle Proper Dih. Improper Dih. LJ-14
    1.52227e+03 7.25999e+03 1.15737e+04 2.23260e+02 1.02350e+04
     Coulomb-14 LJ (SR) Coulomb (SR) Potential Pressure (bar)
    9.97956e+04 -1.76603e+04 -1.87376e+05 -7.44266e+04 1.77431e-01

>From Sander with coordinates file but not the same

 NSTEP = 0 TIME(PS) = 0.000 TEMP(K) = 0.00 PRESS = 0.0
 Etot = -15720.7550 EKtot = 0.0000 EPtot = -15720.7550
 BOND = 348.4593 ANGLE = 1731.7886 DIHED = 2819.0630
 1-4 NB = 2449.9439 1-4 EEL = 23852.5075 VDWAALS = -4733.5065
 EELEC = -42775.8491 EHBOND = 0.0000 RESTRAINT = 0.0000
 CMAP = 586.8384


SANDER* (AT21) Kcal/mol** (KJ/mol)*** with 1kjoul/mol = 4.1867 Kcal/mol
**
step 0
BOND = 348.4593* = 1458.89** ---> Similar to GROMACS results above
ANGLE = 1731.7886 = 7250.479** ---> Similar to GROMACS results above
DIHED = 2819.0630 = 11802.57** Similar to GROMACS results above (Proper Dih + Improper Dih)
 CMAP = 586.8384 = 3.32921e+03 ---> Different from GROMACS obtained with ParmEd (3.32921e+03)


* | mdin Single point
; see http://ringo.ams.stonybrook.edu/index.php/2012_AMBER_Tutorial_with_Biotin_and_Streptavidin
&cntrl
imin=0,
ntpr=1,
maxcyc=100,
ntmin=2,
ntb=0,
igb=0,
cut=999
/

on can say

1) The converted file generated by ParmEd is well recognized by GROMACS and in particular the CMAP term. IT is not the case from the file obtained bt ACEPYPE (since it does not support CMAP for AMBER yet)
2) CMAP energy in AMBER and GROMACS are different
2) the calculations of the energy for the other bonded and non bonded are the same (with rounding error)
2) WARNING !!! with the same inpcrd Parmed does not genertae a corrct GRO file. The box is to too small
   3.24900 3.12110 3.31620 ---> end gro file from ParmEd
   117.23400 182.43200 135.66400 ---> end gro file from ACEPYPE
3) Note that at this moment ParmEd seems not support mixing 1-4 LJ if the prmtop contains atom types from GRLYCAM06 and AMBER19SB ffields. The program crashes with the following error

Traceback (most recent call last):
  File "Parmed_Conversion_GROMACS_AMBER19SB.py", line 3, in <module>
    parm.save("ShuA_BOG_AMBER19SB_GLYCAM06_HMASS_From_Parmed.top") # .top saves to GROMACS topology file format
  File "/home/stephane_abel/amber20/lib/python3.9/site-packages/parmed/structure.py", line 1489, in save
    s = gromacs.GromacsTopologyFile.from_structure(self)
  File "/home/stephane_abel/amber20/lib/python3.9/site-packages/parmed/gromacs/gromacstop.py", line 1279, in from_structure
    raise GromacsError('Structure has mixed 1-4 scaling which is '
parmed.exceptions.GromacsError: Structure has mixed 1-4 scaling which is not supported by Gromacs

HTH and let me know if you have some comments

Best

Stéphane_____
De : ABEL Stephane
Envoyé : mercredi 17 novembre 2021 11:26
À : AMBER Mailing List
Objet : [PROVENANCE INTERNET] Re: [AMBER] Can Parmed convert AMBER files to GROMACS made with the AMBER19SB ?

Hi Jason

Thank you for taking the time to answer me. I

I have installed ParmEd with all the updates to modify my prmtop file for simulations with HMR and it works great. So I will try it for the conversion from AMBER to GROMACS

Best

Stéphane

----------------------------------------------------------
Stéphane Abel, Ph.D., HDR
CEA Centre de Saclay
DRF/JOLIOT/I2BC-S/SB2SM/LBMS
Bat 528, Office 138C
Gif-sur-Yvette, F-91191 FRANCE
Phone (portable) : +33 6 49 37 70 60
________________________________________
De : Jason Swails [jason.swails.gmail.com]
Envoyé : mercredi 17 novembre 2021 04:16
À : AMBER Mailing List
Objet : Re: [AMBER] Can Parmed convert AMBER files to GROMACS made with the AMBER19SB ?

Hi Stephane,

I believe that ParmEd will convert Amber19SB to GROMACS just fine -- I
recommend just trying it. It's possible that the use of different CMAPs
between the same atom type names will cause issues in ParmEd when
translating to GROMACS (for non-CMAP terms, ParmEd always supports the
option of writing parameters directly to each valence term, but I don't
remember how the CMAP terms work).

You will have to write a script to do the conversion:

```
import parmed as pmd
parm = pmd.load_file("ff19SB.parm7", "ff19SB.rst7")
parm.save("ff19SB.top") # .top saves to GROMACS topology file format
parm.save("ff19SB.gro")
```

And execute it as `python <file.py>`. Does that work for you?

I recommend running an energy evaluation with both GROMACS and Amber to
ensure the energies match, however.

HTH,
Jason

On Mon, Nov 15, 2021 at 8:12 AM ABEL Stephane <Stephane.ABEL.cea.fr> wrote:

> Dear all,
>
> from a previous mail of Jason, I can use parmed for my task. But I have a
> question about the conversion feature of Pamed.
>
> It is possible to use it for converting AMBER parameter files into GROMACS
> ones if use the AMBER19SB force field since this force field contains an
> additional CMAP term in contrast to other AMBER force fields for protein?
>
> Thanks for your reply.
>
> Best
>
> Stephane
>
> ________________________________________
> De : ABEL Stephane
> Envoyé : mardi 9 novembre 2021 18:35
> À : AMBER Mailing List
> Objet : [PROVENANCE INTERNET] Re: [AMBER] Hydrogen Mass Repartitioning
> without ParmEd
>
> Thanks again for your help Jason. I will try with miniconda
>
> As I have also seen that you posted a message in the ParmEd forum aboit
> the use of parm7 and rst7 (generated with CHARMM-GUI) files with acepype
>
> *https://githubmemory.com/repo/ParmEd/ParmEd/issues/1117
>
> Q: It is possible to use acepype with a parm7 and rst7 files generated by
> CHARMM-GUI with HMR
>
> A: These are the same file types. parm7 is a synonym for prmtop and rst7
> is a synonym for inpcrd. You can use ACPYPE as-is."
>
> AFAIK ut is not documented in the acepype website
>
> So Thank you a lot Jason (and Saleh) with your help I think I have
> resolved my problem.
>
> Best
>
> Stéphane
>
> ________________________________________
> De : Jason Swails [jason.swails.gmail.com]
> Envoyé : mardi 9 novembre 2021 17:46
> À : AMBER Mailing List
> Objet : Re: [AMBER] Hydrogen Mass Repartitioning without ParmEd
>
> On Tue, Nov 9, 2021 at 10:21 AM ABEL Stephane <Stephane.ABEL.cea.fr>
> wrote:
>
> > Hi Jason,
> >
> > >> tleap does not currently support writing a topology file with
> > repartitioned masses, and ParmEd is the only option for repartitioning
> > within Amber.
> > OK
> >
> > >> If you describe what you mean by "some reasons", maybe it's
> > straightforward
> > My aim is to prepare some Amber files for subsequent simulations with
> > GROMACS so I have installed a minimal version of AmberTools 21 on an HPC
> > machine without python support for using only tleap or Antechamber.
> >
> > ParmEd needs some python dependencies (e.g. setuptools, versioneer) and I
> > cannot download them to compile the code. That why I asked this question.
> > But as there is no other alternative to have the HMR functionality I
> will
> > try to install the full version of Ambertools on a desktop machine to
> have
> > ParmEd works
> >
>
> You can install a copy of Miniconda on the HPC most likely, then conda
> install parmed to make it available. The only runtime dependencies that
> ParmEd has is numpy (setuptools and versioneer are only needed to build).
>
> Note you can also use ParmEd to convert to GROMACS after repartitioning.
> This python script should do it:
>
> -------------------
>
> import parmed as pmd
>
> parm = pmd.load_file("prmtop", xyz="inpcrd")
> pmd.tools.HMassRepartition(parm, 3.024).execute()
> parm.save("gromacs.top", format="GROMACS")
> parm.save("gromacs.gro", format="GRO")
>
> -------------------
>
> HTH,
> Jason
>
> --
> Jason M. Swails
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>


--
Jason M. Swails
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Thu Nov 18 2021 - 04:30:03 PST
Custom Search