From: DAI, JIAN (jdai2_at_fsu.edu)
Date: Mon Jul 30 2012 - 20:54:15 CDT
Hi, James and Aron:
As can be seen from the definition for the CV below, the window is 2 angstrom by 2 angstrom and the binsize is 0.2.
I'm thinking of using smaller windows, but then I suspect that even with smaller windows, the sampling behavior will probably remain, i.e., with boundaries sampled better than the interior of the window.
Does this have to do with the initial configurations? In this case, we used targeted MD to direct the protein to change conformation from closed state to the open state and in this trajectory we pick initial configurations for each window.
Thanks.
Jian
Here goes my definition for the CV:
colvarsTrajFrequency 1000
colvarsRestartFrequency 1000
colvar {
name distance1
width 0.2
lowerboundary 7.5
upperboundary 9.5
lowerwallconstant 10.0
upperwallconstant 10.0
distance {
group1 {
atomNumbers {1870 1889 1903}
}
group2 {
atomNumbers {2790 2801}
}
}
}
colvar {
name distance2
width 0.2
lowerboundary 7.0
upperboundary 9.0
lowerwallconstant 10.0
upperwallconstant 10.0
distance {
group1 {
atomNumbers {152 169 186}
}
group2 {
atomNumbers {3240 3250}
}
}
}
abf {
name distance12abf
colvars distance1 distance2
fullSamples 500
}
As well as the config file:
amber on
parmfile ../../../01Prep/system.parm
ambercoor ../../../01Prep/system.crd
bincoordinates ../../../02Equil/md.1.coor
extendedSystem ../../../02Equil/md.1.xsc
set outputname md.1
firsttimestep 0
set temperature 310
temperature $temperature
exclude scaled1-4
1-4scaling 0.83333
cutoff 12.
switching off
timestep 1.0 ;# 2fs/step (only if needed to finish quickly)
rigidBonds water ;# needed for 2fs steps
nonbondedFreq 2
fullElectFrequency 4
stepspercycle 20
langevin on ;# do langevin dynamics
langevinDamping 1 ;# damping coefficient (gamma) of 5/ps
langevinTemp $temperature
langevinHydrogen off ;# don't couple langevin bath to hydrogens
wrapAll on
PME yes
PMEGridSpacing 1.0
useGroupPressure yes ;# needed for rigidBonds
useFlexibleCell no
useConstantArea no
langevinPiston on
langevinPistonTarget 1.01325 ;# in bar -> 1 atm
langevinPistonPeriod 100.
langevinPistonDecay 50.
langevinPistonTemp $temperature
outputName $outputname
restartfreq 1000 ;# 1000steps = every 1ps
dcdfreq 5000
xstFreq 1000
outputEnergies 1000
outputPressure 1000
colvars on
colvarsConfig ../colvars.tcl
run 2000000
________________________________________
From: Aron Broom [broomsday_at_gmail.com]
Sent: Monday, July 30, 2012 9:33 PM
To: JC Gumbart
Cc: DAI, JIAN; namd-l_at_ks.uiuc.edu
Subject: Re: namd-l: Uneven sampling when using ABF
Yeah, I knew that's what fullSamples meant, but wouldn't the initial estimate from filling all those samples still give you a good estimate of the PMF?
I had thought the accuracy only depended on how many samples you had for each window, presuming that you at least explore all the windows (i.e. reach fullSamples for each window). I had been thinking the only purpose for applying the biasing force was to make sure you explored past the initial local minimum in order to collect force samples over the whole reaction coordinate. Does the biasing force do something other than allow exploration?
Still, if there are some unbiased slow degrees of freedom that are causing the force estimate to converge very slowly, I think you can only get a good result from ABF by having fullSamples high enough to capture those unbiased degrees of freedom.
~Aron
On Mon, Jul 30, 2012 at 9:15 PM, JC Gumbart <gumbart_at_ks.uiuc.edu<mailto:gumbart_at_ks.uiuc.edu>> wrote:
I think you misunderstand what fullSamples means. The fullSamples parameter is how many steps per bin are required before the estimated biasing force is applied. While I agree that 500 is too small, you don't want it to be some absurdly large number either, because then you are not doing ABF, you're just doing pure Boltzmann sampling. If the initial estimate of the force is very bad, then it's possible that ABF can never correct for it, at least on a reasonable timescale. Try setting this parameter to, say, 10000.
Also, what is the bin size you're using? Are there any kind of restraints being applied to the system?
On Jul 30, 2012, at 6:10 PM, Aron Broom wrote:
for whatever reason, the link to your image kept timing out for me.
I have only very limited experience with ABF, but I think you are simply describing non-convergence, and the simplest solution is to sample for longer. That being said, if the force applied in certain areas is particularly far from the real average force over the whole ensemble, you could become trapped in a non-equilibrium condition for a VERY LONG time.
In my experience the solution to this is to increase fullSamples. You say it's currently set at 500, so that is 500 timesteps, or 500 fs. And you say the dimensions are 2 angstroms, which means that the system should fully explore all the relevant microstates within that 2 angstrom window, within 500 fs. In my experience, if you are looking at something like protein-ligand binding, you have a decorrelation time on the order of thousands of fs in explicit solvent, and I would expect if you are looking at protein folding it would be even higher. So my conclusion would be that your current "fullSamples" is only capturing one statistically relevant unique sample per window, and therefore, is highly unreliable. I would expect you'd need to increase that by at least several orders of magnitude to see good results.
Personally I stopped using ABF, because once you run into a problem like this, it seems like you more or less need to start over with a different fullSamples value. MetaDynamics may have similar problems, but if you use the default hill height and a deposit hills every 2000 fs or so, you can often get a rough idea, that can be refined later. For not wasting simulations, Umbrella sampling is possibly the best, and you can do it in 2D.
Keep in mind that the protein folding problem has a tremendous number of degrees of freedom, so there are many that ABF is not biasing against, meaning that it should take a lot of simulation time to get good data. From what you say it sounds like you currently get 162 (9x9x2) ns of simulation time. Maybe something to try just a sanity check would be to redo the ABF, but set it up such that at least half that time (81ns) is the time to fill all the "fullSamples" (in this case fullSamples would be 1,000,000 or 1 ns), and see how much different your PMF is from what you currently have with 500 as fullSamples.
~Aron
On Mon, Jul 30, 2012 at 1:38 PM, DAI, JIAN <jdai2_at_fsu.edu<mailto:jdai2_at_fsu.edu>> wrote:
Dear fellows:
We are trying to construct a free energy landscape of a protein using two dimensional ABF calculation. The landscape is divided into 9 equal sized windows, each with a dimension of 2 Angstrom by 2 Angstrom, where each dimension describes the center of mass between two groups of atoms.
Please see the figure by the link.
http://s663.photobucket.com/albums/uu358/djpittdj/?action=view¤t=count.png
The figure shows the number of counts with respect to those two order parameters, and they are not evenly distributed in each window, as I would expect. Instead, in four windows on the left bottom corner, the counts are heavily clustered in their respective boundaries, indicated by those red/yellow regions. Does anybody knows why and can possibly offer us a solution?
The fullSamples parameters is set to 500, and the landscape was obtained using abf_integrate after each window is run for 2 ns with a timestep of 1 fs.
Thanks a lot.
Jian
-- Aron Broom M.Sc PhD Student Department of Chemistry University of Waterloo -- Aron Broom M.Sc PhD Student Department of Chemistry University of Waterloo
This archive was generated by hypermail 2.1.6 : Tue Dec 31 2013 - 23:22:19 CST
