In my last post I described how to check convergence of metadynamics simulation. However, the 1D free energy surface showed in the previous post is not reweighted and corrected properly because I had both upper and lower wall in my metadynamics simulation. Hence we have to correct the free energy surface for the walls and we can do it using Plumed-2.4.1. It is important to mention that I ran my simulations with Plumed-1.3 and have to do some adjustments to get a 1ps sparsed COLVAR file whose header looks like following
#! FIELDS time d1 d2 d3 metad.rbias uwall.bias 1.0000 0.297346503 0.718053460 1.629577041 0.000000000 1.516799331 2.0000 0.312700063 0.676572204 1.619366407 0.500000000 1.017568588 3.0000 0.337963104 0.671517015 1.576508045 0.437842011 0.172387958 4.0000 0.304164201 0.674881995 1.545186877 0.981566668 0.022762353 5.0000 0.281478077 0.678171694 1.506772399 0.851512551 0.003159886 6.0000 0.343289763 0.649008214 1.564545155 1.156349301 0.086780988 7.0000 0.314589500 0.690125763 1.542209387 1.231884360 0.028832838 8.0000 0.297886729 0.662928343 1.551209688 1.714575529 0.034525435 9.0000 0.293216735 0.665717304 1.541974306 1.670215130 0.015825573
Worth mentioning that the original COLVAR output from Plumed 1.3 looks like following
#! FIELDS time cv1 cv2 cv3 vbias vwall 0.0000 0.331873775 0.679191828 1.648215294 0.000000000 2.416546345 0.1000 0.305706650 0.693865001 1.633054614 0.000000000 1.585586548 0.2000 0.289182395 0.736790895 1.625117898 0.000000000 1.509020209 0.3000 0.306852549 0.739897966 1.635408998 0.000000000 2.007535696 0.4000 0.289105147 0.741805196 1.605738878 0.000000000 0.980214059 0.5000 0.314303517 0.736100793 1.611924648 0.000000000 1.059433341 0.6000 0.316665620 0.693137527 1.657169104 0.000000000 3.068289280 0.7000 0.313074738 0.689189851 1.670769334 0.000000000 4.263953686 0.8000 0.281177700 0.716560900 1.678995728 0.000000000 5.230778694
I edited cv1, cv2, cv3 to d1, d2 and d3 and changed vbias to metad.rbias and vwall to uwall.bias just for naming convenience.
vwall in the original COLVAR defines the external wall potential as a combination of both upper and lower wall.
Note: You have to delete the values corresponds to 0.000 ps from the final COLVAR file.
The method is well described here.
See the following script which produces 1D reweighted free energy surface on CV1 which is d1 in the new COLVAR 1ps sparsed file. You will get two outputs fes.dat and fes_wall.dat and try to see how different the wall corrected (fes_wall.dat) and with wall free energy surfaces (fes.dat) look like. It will look something like following
Figure 1. 1D reweighted free energy surface on CV1.
I used the following script
# Read COLVAR file distance: READ FILE=COLVAR_sparse IGNORE_TIME VALUES=d1 metad: READ FILE=COLVAR_sparse IGNORE_TIME VALUES=metad.rbias uwall: READ FILE=COLVAR_sparse IGNORE_TIME VALUES=uwall.bias # Define weights weights1: REWEIGHT_METAD TEMP=300 ARG=metad.rbias weights2: REWEIGHT_BIAS TEMP=300 ARG=metad.rbias,uwall.bias # Calculate histograms HISTOGRAM ... ARG=distance GRID_MIN=0.0 GRID_MAX=1.0 GRID_BIN=100 BANDWIDTH=0.002 LOGWEIGHTS=weights1 LABEL=hh1 ... HISTOGRAM HISTOGRAM ... ARG=distance GRID_MIN=0.0 GRID_MAX=1.0 GRID_BIN=100 BANDWIDTH=0.002 LOGWEIGHTS=weights2 LABEL=hh2 ... HISTOGRAM # Print histograms to file ff1: CONVERT_TO_FES GRID=hh1 TEMP=300 DUMPGRID GRID=ff1 FILE=fes.dat DUMPGRID GRID=hh1 FILE=histo FMT=%24.16e ff2: CONVERT_TO_FES GRID=hh2 TEMP=300 DUMPGRID GRID=ff2 FILE=fes_wall.dat DUMPGRID GRID=hh2 FILE=histo_wall FMT=%24.16e
Ofcourse we can reweight the free energy surface on both CV1 and CV2 using the following script
# Read COLVAR file distance1: READ FILE=COLVAR_sparse IGNORE_TIME VALUES=d1 distance2: READ FILE=COLVAR_sparse IGNORE_TIME VALUES=d2 metad: READ FILE=COLVAR_sparse IGNORE_TIME VALUES=metad.rbias uwall: READ FILE=COLVAR_sparse IGNORE_TIME VALUES=uwall.bias # Define weights weights1: REWEIGHT_METAD TEMP=300 ARG=metad.rbias weights2: REWEIGHT_BIAS TEMP=300 ARG=metad.rbias,uwall.bias # Calculate histograms HISTOGRAM ... ARG=distance1,distance2 GRID_MIN=0.0,0.0 GRID_MAX=1.0,1.0 GRID_BIN=200,200 BANDWIDTH=0.05,0.05 LOGWEIGHTS=weights1 LABEL=hh1 ... HISTOGRAM HISTOGRAM ... ARG=distance1,distance2 GRID_MIN=0.0,0.0 GRID_MAX=1.0,1.0 GRID_BIN=200,200 BANDWIDTH=0.05,0.05 LOGWEIGHTS=weights2 LABEL=hh2 ... HISTOGRAM # Print histograms to file ff1: CONVERT_TO_FES GRID=hh1 TEMP=300 DUMPGRID GRID=ff1 FILE=fes.dat DUMPGRID GRID=hh1 FILE=histo FMT=%24.16e ff2: CONVERT_TO_FES GRID=hh2 TEMP=300 DUMPGRID GRID=ff2 FILE=fes_wall.dat DUMPGRID GRID=hh2 FILE=histo_wall FMT=%24.16e
Here I plotted fes.dat again which is two dimensional free energy surface reweighted on CV1 and CV2. You can plot fes.dat and fes_wall.dat in a 2D plot to see what was the effect of wall. Here I plotted fes_wall.dat which is wall corrected 2D FES.
Figure 2. 2D reweighted free energy surface on CV1 and CV2.
Note: Sorry for bad quality rendering (email me for better picture)
Reweighting algorithm described here
Pratyush Tiwary and Michele Parrinello. A time-independent free energy estimator for metadynamics. The Journal of Physical Chemistry B, 119(3):736–742, 2015. PMID: 25046020.
Manual of how to use reweighting metadynamics described here.
