Re: Umbrella sampling-NAMD

From: Aron Broom (broomsday_at_gmail.com)
Date: Fri Dec 16 2011 - 09:46:39 CST

Parisa,

If K=1 in NAMD, and that is following the convention of V = 1/2Kx^2
(according to the user guide), why would you scale K when going to WHAM?
According to the WHAM manual it also assumes the form V = 1/2Kx^2. Can you
elaborate on why you want to change the spring constant? If you are using
the NAMD colvars module to apply a harmonic restraint it doesn't matter
about CHARMM, that harmonic restraint isn't coming from your CHARMM
parameters.

In terms of number of bins when using WHAM, you should perhaps analyze all
of your data and see at what resolution you start to lose having overlap
between windows and sufficient data. There should be a statistically
meaningful amount of data in each bin in order for them to be linked
properly. My impression when using Alan Grossfield's WHAM code was that
the bin resolution actually had very little effect on the shape and
magnitude of my PMF, but perhaps I just had a lot of data and so didn't see
the issue. A simple choice would be to make the resolution the average of
the consecutive distances between your windows.

As I mentioned before though, you do need to manually trim your data for
each window, to remove equilibration time at the start and account for
decorrelation time for your system, otherwise your results will be
inaccurate and you will think that you have more data in each window than
you really have (usable data that is)

~Aron

On Fri, Dec 16, 2011 at 10:29 AM, Parisa Akhski <
Parisa.Akhshi_at_chem.queensu.ca> wrote:

> **
>
> Hi Chris,
>
> Thanks very much for your e-mail.
>
> I have already seen your previous discussion in the mailing list. K=100
> dose not produce the expected results.
>
> Here are the issues I am running to:
>
> 1. According to the user manual:
>
> forceConstant < (harmonic) Harmonic force constant (kcal/mol) >
> Acceptable Values: positive decimal
> Default Value: 1.0
> Description: This defines the constant for the harmonic potential. To
> ensure consistency
> for multidimensional restraints, it is scaled by the square of the
> specific width for each colvar
> involved (which is 1 by default), so that all colvars are internally
> rendered dimensionless and
> of comparable size. For instance, a force constant of 10 kcal/mol acting
> on a dihedral colvar
> with a width of 2 degrees will amount effectively to 2.5 kcal/mol×degree-2.
>
> In my previous e-mails I was thinking that NAMD ALWAYS scales K. But, now
> I see that the reason to scale K is to fix the dimensions in
> multidimensional harmonics. Now, my understanding is that for 1D (applying
> harmonic on one colvar) NAMD does not scale K since it already has correct
> dimensions. For 2 and 3 colvars, it scales accordingly. Is this correct? If
> so, I am using K=1 and I only have one colvar (1D) so NAMD applies K=1.
> This tells me that I have to use K=2 in WHAM which gives me very small
> energy barriers and very far from my expectations.
>
> In response to Jerome's question, going from K=200 to K=100 and K=50
> produces PMF with half and one-fourth height, respectively. But, for
> smaller values of K (K=2) the dependence is not linear.
>
> 2. Number of bins in WHAM: I am still not clear how to choose the "number
> of bins" in "WHAM"? My results are highly dependent on this parameter
> especially when the force constant is high (K=200). How should I decide on
> the "number of bins" in "WHAM"?
>
> Cheers,
>
> Parisa
>
>
>
> -----Original Message-----
> From: Chris Harrison [mailto:charris5_at_gmail.com <charris5_at_gmail.com>]
> Sent: Thu 12/15/2011 3:07 PM
> To: Parisa Akhski
> Cc: Aron Broom; namd-l_at_ks.uiuc.edu
> Subject: Re: namd-l: Umbrella sampling-NAMD
>
> Parisa,
>
> Please refer to Jerome's response as it is absolutely correct, and has
> been previously discussed on the namd-l list. I explicitly remember
> answering the same questions less than 2 months ago.
>
> So as Jerome asked, how does the profile look when WHAM is used at K=100?
>
>
> Best,
> Chris
>
>
> --
> Chris Harrison, Ph.D.
> Theoretical and Computational Biophysics Group
> NIH Resource for Macromolecular Modeling and Bioinformatics
> Beckman Institute for Advanced Science and Technology
> University of Illinois, 405 N. Mathews Ave., Urbana, IL 61801
>
> char_at_ks.uiuc.edu Voice: 773-570-0329
> http://www.ks.uiuc.edu/~char Fax: 217-244-6078
>
>
> Parisa Akhski <Parisa.Akhshi_at_chem.queensu.ca> writes:
> > Date: Thu, 15 Dec 2011 13:51:57 -0500
> > From: Parisa Akhski <Parisa.Akhshi_at_chem.queensu.ca>
> > To: Aron Broom <broomsday_at_gmail.com>
> > CC: namd-l_at_ks.uiuc.edu
> > Subject: RE: namd-l: Umbrella sampling-NAMD
> >
> > Thanks Aron for your comment,
> >
> > Yes, I have checked the convergence of ABF results and they are OK. For
> umbrella sampling the only criteria I know to check is the overlap between
> windows which is OK. Is there anything else I could check for US to make
> sure it has converged?
> >
> > Many thanks,
> >
> > Parisa
> >
> >
> >
> >
> >
> >
> > -----Original Message-----
> > From: Aron Broom [mailto:broomsday_at_gmail.com <broomsday_at_gmail.com>]
> > Sent: Thu 12/15/2011 12:11 PM
> > To: Parisa Akhski
> > Cc: namd-l_at_ks.uiuc.edu
> > Subject: Re: namd-l: Umbrella sampling-NAMD
> >
> > Parisa,
> >
> > If you had a very simple system like two molecules in the gas phase then
> I
> > would be worried by the disagreement between Umbrella Sampling and ABF,
> but
> > given the complexity of your system, are you sure that the ABF
> calculations
> > have given an accurate PMF? Have you checked that doubling the sampling
> > parameter and/or doubling the simulation time has no effect on your ABF
> > PMF? I would intuitively suspect that the difference is because one of
> > your analysis methods has not really converged properly.
> >
> > ~Aron
> >
> > On Thu, Dec 15, 2011 at 10:24 AM, Parisa Akhski <
> > Parisa.Akhshi_at_chem.queensu.ca> wrote:
> >
> > > **
> > >
> > > Dear NAMD experts,
> > >
> > > I am using NAMD to run Umbrella Sampling on an ion channel. I am using
> a
> > > "harmonic force constant" to keep the ion in the center of bin in each
> > > window and I know that this value is scaled by the square of the
> specific
> > > width as mentioned in the manual. So, if I define "forceConstant=1"
> and use
> > > width=0.1, NAMD actually applied K=100 in the calculations.
> > >
> > > Then, I am processing the trajectory files by using WHAM to produce
> PMF. I
> > > know that WHAM assumes that the biasing potential is V=1/2K(X-X0)^2
> whereas
> > > CHARMM just defines it as V=K(X-X0)^2. This means if I am using K=100
> in
> > > NAMD, I should define K=200 in WHAM.
> > >
> > > Since I already have the PMF profile from ABF calculations for the same
> > > system, I know how the profile should look like and what are the energy
> > > barriers in my system. Comparing theses results with the ones from
> Umbrella
> > > Sampling and WHAM, the shape of the profile makes sense but not the
> amount
> > > of the barriers. However, if I change the K value from 200 in WHAM to
> 50,
> > > everything makes sense. I double checked this with WHAM experts to make
> > > sure I am using the correct K (200) in WHAM.
> > >
> > > Could you please let me know whether you have any suggestions,... to
> solve
> > > this issue?
> > >
> > > Many thanks for your time,
> > >
> > > Parisa
> > >
> >
>
>
>
>
>

This archive was generated by hypermail 2.1.6 : Mon Dec 31 2012 - 23:21:04 CST