On 12/30/12 2:57 PM, Nash, Anthony wrote:
Dear Justin
Thanks for your reply.
I am studying a TM peptide, looking at how favourable the front and reverse
face are. Hence, different orientations but the same composition of amino
acids.
I then have a duplicate of this system (call it system B), but with a couple
of residue substitutions. I am comparing free energies of System A front
orientation with System A back orientation (and then B with B).
I have separated the peptides, and ran umbrella simulations every 1 angstrom
apart (I experimented with windows closer together, in some regions hence the
generalisation to half an angstrom). Papers I have looked through which
sample up to 2 nm apart have been running umbrella windows every 0.2 of an
angstrom. My PMF plateau is around 6.5-7.5 nm, which implies around 78
windows per system.
Without seeing the actual histograms and .mdp settings, I can't comment on this
approach directly. It seems like you should be able to use far fewer windows,
perhaps with a slightly less restrictive force constant to allow for slightly
broader distributions. That's a hand-waving guess, though, since I've seen none
of the data.
Taking one system as an example, if I normalise at 7.5 nm the difference at
their global minima is around 37.63 kJ mol−1 (8.99 kCal mol−1), but then when
I break down the PMF graph (0-1ns, 0-2ns, 0-3nm) the regions greater than 4.5
nm have not converged. Hence, normalising anywhere in this unconverged region
would be wildly inaccurate. If I normalise to zero at the point in which one
of the system pairs (A&A, or B&B) converges, the free energy different is 5
kJ mol−1.
I think it would be more productive to do your block analysis such that you
discard initial frames to see what effect the (presumably) unequilibrated part
of the trajectory is. Even with prior equilibration (i.e. with position
restraints), you still need to leave some frames behind as equilibration.
Analyzing 0-40 ns, 10-40 ns, 20-40 ns, etc will probably be more informative.
If you're always considering all frames in the WHAM analysis, that can throw off
the results considerably.
If by error estimate you mean from the inbuilt bootstrap analysis, then I am
very confident that the error along the entire reaction coordinate is very
small; which is confusing as I am not sure whether the slower converging
regions of the reaction coordinate (4.5 - 7.5) should shower greater error!?
That's not uncommon. At greater distance, you have two peptides floating around
rather independently, so the restraint potential may vary a bit more. At closer
COM distances, the interactions between the peptides keep the orientations a bit
more stable.
Yes it is heart breaking to hear that I may need to run the system for
longer. This however, isn't a massive problem, given the time of the year the
HPCs I am using are a bit empty. Would you suggest I need to run the 40 ns
for longer (continuation run), or begin brand new windows from scratch and
include those in the g_wham calculations? I guess the 40 ns one, given that
lipid reorientation needs to be taken into account.
There is no need to scrap the existing data. You may only need to extend the
simulations, but do carry out the analysis I suggested before to see if there is
undue influence from the initial frames.
I don't have numbers in front of me, but the system box size is big (I went
through your tutorial several times, and I have seen the PMF graph as a
result of a box too small). The longest i.e., the reaction coordinate, is
along the Y axes, Z is normal to the bilayer. When I calculate the PMF, I
have pull_dim = Y Y N.
Good to know. Lipid ordering effects can extend up to a few nm though, so
that's why I asked.
-Justin
--
========================================
Justin A. Lemkul, Ph.D.
Research Scientist
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
--
gmx-users mailing list [email protected]
http://lists.gromacs.org/mailman/listinfo/gmx-users
* Please search the archive at
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
* Please don't post (un)subscribe requests to the list. Use the
www interface or send it to [email protected].
* Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
