[gmx-users] The "correct" way to equilibrate a membrane / lipid bilayer w/water?

[Justin A. Lemkul]

Peter C. Lai wrote:
> Hi 
> 
> I haven't played with gromacs in years, so I'm pretty new at the ability
> to play with large systems like lipid bilayers.
> 
> What is the "correct" or lets say canonical method to equilibrate a 
> membrane patch of say 9x9 nm POPC with sufficient vertical waters to take
> care of pbc and long range interactions (> 1.5nm above/below the bilayer)
> using Tom's Charmm36 FF port.
> 
> I am seeing a lot of different methods out there to: construct and
> equilibrate such a system in prepration for g_embed or whatever.
> 
> I see from 
> http://www.mail-archive.com/gmx-users@gromacs.org/msg33812.html
> the ability to use genbox to replicate Tieleman's 128 patch 
> 
> but I also can build my bilayer inside VMD, which supports charmm36->pdb
> atom names (and things seem to work fine)
> 
> I have yet to see anybody's run parameters and constraints for running
> the equilibration, particularly using charmm36 (using vdwtype=switch).
> Should I constrain all heavy atoms in the lipid during the NVT run, or
> just the polar group (or even just P) and allow the tails to "melt" as seen 
> elsewhere? Should I only set Z axis constraints? How long should I run in 
> NVT before switching to NPT and for how long (10ns?) And should I change the
> constraints moving from NVT to NPT?

As a point of clarity, "constraints" and "restraints" serve very different 
functions in Gromacs:

http://www.gromacs.org/Documentation/Terminology/Constraints_and_Restraints

I presume you mean "restraints" in this context.  I see no point in restraining 
the lipids completely during equilibration.  It sort of defeats the purpose. 
Depending on how reasonable the starting configuration is, restraints on any 
group may not be necessary at all.  Poor starting configurations, particularly 
with inadequate solvation within the interfacial region, often require vertical 
position restraints on P atoms to prevent the membrane from separating. 
Restraining the entire headgroup probably restricts lipid motion too much.

There is no standard timeframe for equilibration.  Rotational relaxation occurs 
within a very short amount of time, translational within 10-20 ns for most 
lipids.  You should pay close attention to diffusion constants, membrane 
thickness, area per lipid, etc to judge the quality of equilibration.

> How should the thermostat be coupled? System or separate the water from the 
> lipid?
> 

This consideration is uniform among just about all systems and is a function of 
the thermostat rather than the system itself.  Couple water and lipids 
separately, just as you would couple protein and solvent separately for most 
thermostats.

COM motion should be removed separately for water and lipids as well, to prevent 
lateral sliding that is hidden through global COM motion removal.

> I experienced the same thing as Ng Hui Wen back in October with waters
> being extremely dehomogenized (and diffusing out of the box) after 100ps 
> NVT with all heavy POPC constrained at 1000 kJ/nm for all axes. I am sure
> that moving to NPT will restore the waters but still it was a bit scary to 
> see for the first run.
> 

That is probably an artifact of restraining all lipid atoms and not allowing for 
proper response between water and lipids.

> If a lot of these approaches haven't been published, maybe we can just
> stick them in the wiki to let everyone "pick their poison" when it comes
> to equilibrating their bilayers.
> 

I can think of at least a dozen plausible methods for equilibrating a lipid 
bilayer system, making it probably more confusing to simply list them and tell 
users to pick, especially if they don't know what they're doing.  As with any 
equilibration of any system, you have to judge (1) whether the desired 
thermodynamic ensemble is achieved and (2) whether the properties of the system 
of interest are established.  With these general rules in mind, the user is free 
to design their own procedure.

I wrote a membrane protein tutorial a long time ago that outlines many of these 
principles and discusses several more.

http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/membrane_protein/index.html

Even though your purposes may not involve a protein or the same force field, you 
may find some of the general information useful.

-Justin

> Thanks for any ideas/hints/suggestions
> 

-- 
========================================

Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

========================================