[gmx-users] FEP

David Mobley dmobley at gmail.com
Thu Mar 23 19:57:55 CET 2006


 Hi all.
> I'm new to gromacs, and I would be glad to get some advice from more
> experienced usres.
> I would like to estimate the difference in free energy between two
> systems.
>  First system (A) consists of  DPPC membrane, TM-protein and some Zn ions
> interacting with the protein.  Second (B) has the same elements but Zn
> ions
> are far away from the protein (in the solvent).  Is it possible to
> estimate
> free energy difference between those two systems in gromacs?
> I thought about using dummy atoms and slow-growth methods for free energy
> calculations.

  During the simulation, Zn ions from system (A) would be perturbated into
> dummy atoms
>  (charge=0, mass = 0) and at the same time equal number of dummy atoms in
> the solvent
> would be perturbated to Zn atoms?  Is it possible to achieve?

I am rather nervous about doing FEP/TI for disappearing charged molecules,
as it is not at all clear to me that it is possible do this correctly with
current methods. Perhaps someone else may be able to comment more, but at
least with long range electrostatics (PME), systems are required to be
neutral. It IS possible to apply PME to a "charged" system, but I think this
is typically done by spreading a uniform neutralizing charge throughout the
box. It isn't clear that this is really correct for FEP/TI, I don't think,
since in that case you will have a charged molecule overlapping with some of
the uniform neutralizing charge, which could give strange energy artifacts.

Perhaps things are better if you don't use PME electrostatics, but there may
be boundary effects.

I need to personally dig in to the literature on this a bit more, but I
would suggest doing so yourself unless you're very sure that doing FEP/TI
will give you something meaningful on this system.

Also, if you *do* go that route, I would recommend against slow growth, as
it tends to be hysteretic. It's probably better to run a bunch of separate
simulations at different lambda values.

And even further, it is probably *not* a good idea to turn off the charges
and LJ interactions at the same time. It works better to turn off the
electrostatics first (without using soft core, as this generally is fairly
smooth with the charge), and then turn off the LJ interactions using soft
core (recommend sc-alpha=0.5, as sc-alpha=1.5 is usually too large for LJ).

However, if I were doing it, I would probably try to do it by computing the
PMF for pulling the zinc ions (or one zinc ion, which is easier) away from
the protein. You can compute the free energy difference from the PMF. The
unfortunate thing about this is that it involves sampling a lot of stuff you
don't care about, but at least it is clear the electrostatics is correct.
Alternatively, dig into the literature and try and figure out if there is a
correct way to do the electrostatics in this case. (And please let me know
what you find out, if you do!)

Oh, and one more thing: If DO compute the free energy by annhilating the
zinc ions, you probably should use restraints to keep them in a particular
region as you turn off the interactions, otherwise they will have to sample
the entire simulation box in order for you to get converged results.

Good luck.

> For the beginning  I tried to estimate the difference in free energy
> between two systems.
> 1) first system consisted of two Cl- and one Zn2+ ions + solvent
> 2) second only solvent molecules
> During the 100 ps MD I perturbated  chare and mass of Cl and Zn
> ions to zero.  The energy values that I've obtained are wrong 960200,46.
>  I'm new to gromacs so I would appreciate your help a lot.
> Here are my input files (I'm using ffgmx )
> --------------md.mdp-------------
> pp                 =  /lib/cpp
> define              =  -DFLEX_SPC
> constraints         =  none
> integrator          =  md
> dt                  =  0.001
> nsteps              =  100000
> nstxout             =  100000
> nstvout             =  100000
> nstfout             =  100000
> nstlog              =  100
> nstenergy           =  100
> nstxtcout           =  100
> nstlist             =  10
> ns_type             =  grid
> rlist               =  1.0
> free_energy = yes
> init_lambda = 0
> delta_lambda = 0.00001
> sc-alpha = 1.5
> sc-sigma = 0.3
> ....
> ...
> ....
> ---------------------------------------
> --------ZN.itp-----------------------
> [ moleculetype ]
> ; molname       nrexcl
> ZnB              1
> [ atoms ]
> ; id    at type res nr  residu name     at name  cg nr
> charge              typeB chargeB massB
> 1       Zn      1       Zn              Zn       1      2     65.37000
> Zn       0  0
> ---------------------------------------
> -----------Cl.itp-----------------
> [ moleculetype ]
> ; molname       nrexcl
> ClB              1
> [ atoms ]
> ; id    at type res nr  residu name     at name  cg nr  charge typeB
> chargeB massB
> 1       Cl      1       Cl              Cl       1      -1        35.45300
> Cl  0   0
> ----------------------------
> All the best.
> :)
> _______________________________________________
> gmx-users mailing list    gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20060323/e00abe7a/attachment.html>

More information about the gromacs.org_gmx-users mailing list