Diverging behaviour: sum_hills vs reweighting with WALKERS_MPI

[Alexander Zlobin]

Hello, Plumed community!

I recently got troubled with a strange issue regarding well-tempered metadynamics with multiple walkers. In brief, fes obtained with sum_hills and REWEIGHT_METAD are absolutely different and I do not know which one to trust.

Now in detail:

I performed a run with 14 walkers, with this setup, using Plumed 2.6.1:

METAD ...
  LABEL=M
  ARG=CV
  PACE=100
  HEIGHT=2.0
  SIGMA=0.02
  GRID_MIN=-1
  GRID_MAX=1
  GRID_SPACING=0.004
  BIASFACTOR=40

  WALKERS_MPI

  CALC_RCT
  RCT_USTRIDE=50
 
... METAD

I obtained two profiles. One with classic sum_hills:

And I have reason to believe this since the process I am biasing is a reaction with a stable tetrahedral intermediate state.

The second profile I obtained with plumed driver giving it these instructions:

CV: READ FILE=CV IGNORE_TIME VALUES=CV  (CV is a file with all individual CVs from walkers concatenated)
RW: READ FILE=CV IGNORE_TIME VALUES=M.rbias
weights: REWEIGHT_METAD TEMP=300 ARG=RW.rbias

HISTOGRAM ...
  ARG=CV
  GRID_MIN=-1
  GRID_MAX=1
  GRID_SPACING=0.004
  BANDWIDTH=0.02
  LOGWEIGHTS=weights
  LABEL=histo
... HISTOGRAM

ff: CONVERT_TO_FES GRID=histo TEMP=300
DUMPGRID GRID=ff FILE=fes.dat

This is what I got:

As you can see it looks completely different with no INT state whatsoever and inverted dG0.

So my question is as follows: what am I doing wrong and what way to treat such simulation in a proper way? I feel that CALC_RCT works somewhat erroneous with multiple walkers since I have such curious jumps in rbias in CV files (last column; biased variable is the one before it):

 2.480000 0.107724 0.153952 0.160004 0.105840 0.255623 0.134849 -0.128711 44.255887
 2.490000 0.101515 0.146502 0.162732 0.103308 0.251178 0.143521 -0.122094 42.075668
 2.500000 0.116463 0.123083 0.155535 0.107078 0.252551 0.137729 -0.156660 55.215274
 2.510000 0.101646 0.139324 0.140886 0.114134 0.255720 0.134210 -0.110585 -56.412790
 2.520000 0.125303 0.115797 0.131888 0.111162 0.249621 0.143026 -0.136827 -46.547447
 2.530000 0.097534 0.147590 0.142016 0.101316 0.243125 0.139240 -0.094528 -61.328360

[Michele Invernizzi]

Hi Alexander,

Standard metadynamics reweighing can be performed only once you reach an adiabatic regime (i.e. the bias changes slowly). The jumps you describe are quite normal at the beginning of the trajectory when the simulation is out of equilibrium and the c(t) estimate is still very rough. You should discard this initial transient, and use for reweighting only the adiabatic part of the trajectory.

However, looking at your CV file I am afraid you are still far from convergence. It seems that none of the 14 walkers has done a full transition, going back and forth between the two basins. If you do not have any transition between the basins, then neither the sum_hills nor the reweighting free energy estimate can be trusted.

You should run your simulations longer, and if you still don't see transitions you should consider trying to improve your CV or adding a second one.

Regards,

Michele

--
You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/5bb9028c-2693-462d-b60e-7f7287d80c10%40googlegroups.com.

[Alexander Zlobin]

Hi Michelle,

thank you for your answer! I was aware that my simulation is far from convergence, however, I had a prior intuition that both ways of constructing FES would yield similar (yet both far from "real") profiles. Well, as it turns out, this was a wrong intuition.

I took some time to run another simulation, this time without well tempering at all (and for another substrate, yet the reaction itself is absolutely similar). I observed full recrossing in half the walkers (going back and forth between all three basins) and partial in remaining ones. This is what I got:

Again, sum_hills looks really good in terms of interpretability and idea drawing, and barrier heights are reasonable. Yet reweighting just yields a straight line. I suppose it might be due to this strange lines in the output:

 14.830000 0.253625 0.101506 0.270833 0.107773 0.253797 0.133992 -0.434984 -1348.504252
 14.840000 0.254861 0.098580 0.263720 0.110672 0.260209 0.143319 -0.426219 -1347.627908
 14.850000 0.275446 0.096871 0.270998 0.112210 0.280903 0.135577 -0.482690 -1359.299084
 14.860000 0.283788 0.092843 0.313674 0.107372 0.297964 0.132401 -0.562811 -inf
 14.870000 0.298803 0.091307 0.313440 0.107292 0.301869 0.138387 -0.577125 -inf
 14.880000 0.315248 0.093265 0.296801 0.108408 0.299627 0.135538 -0.574465 -inf

The last column is rbias, and it spontaneously went infinite pretty early on, leaving almost half of simulation time useless. I then performed reweighting only on the frames before rbias turn infinite, yet the result is still nonsense:

Now my plan is to use the sum_hills profile to generate sets of conformations for each basin and then to train neural network CV construction as proposed by Bonati et. al, 2020 to use it in well-tempered profile "polishing". However I still not sure whether everything is okay with this part of Plumed code, the whole situation looks quite shady.

I also ran some tests on a dummy model (simple ch3cl + cl- in vacuum) for which I can be totally sure about convergence and sampling. I got the same situation with rbias going to infinity for non well-tempered simulation, effectively thrashing any reweighting (the term which affects it is rct - it turns inf while bias has reasonable values. It is not the exclusive problem of walkers thought, single run yields exactly the same). Profiles for well-tempered run though were very similar either being constructed with sum_hills or by reweighting. Is my assumption correct that one should not perform reweighting without tempering?

[Michele Invernizzi]

Dear Alexander,

yes, I know, at first the reweighting procedure for metadynamics can be quite confusing (I have been through this as well) but I assure you there is nothing shady in what the PLUMED code does. I suggest you to look at one of the nice tutorials on the plumed website, such as this one.

One key idea is that you can use that reweighting procedure only when the simulation is approaching convergence and the bias is changing only adiabatically. Unfortunately it is not always easy to understand when this is the case. Generally speaking the higher the BIASFACTOR the longer the initial transient, and for non-tempered METAD you might simply never reach this regime (and never converge). This is why the rct in your simulation goes to inf. (In general I would use non-tempered METAD only as an exploration tool, or for quickly testing a new CV, because it does not converge even if you run it for very long)

Being more gentle with the biasing can be a better strategy, especially if your CV is suboptimal ( https://pubs.acs.org/doi/pdf/10.1021/acs.jpclett.0c00497 ).

Here you can see in purple the reweighted FES I obtained with your CV file, by discarding the firs half of the simulation, and using only the second half. In green instead is the sum_hills estimate, that evolves with time. They qualitatively agree, but you can't say much more than that, still too far form convergence.

Hope this helps.

Regards,

Michele

--

You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/e26c08bf-520b-47de-bf26-574b8f726b91%40googlegroups.com.

[Pratyush Tiwary]

The bias, the force experienced by the biased variable, and the histogram of the biased variable (adequately weighted) are 3 estimators of the free energy. They are not entirely independent - indeed the histogram sampled till a point affects the bias, and the bias itself affects the weight used in the histogram. That said, they are not entirely correlated either. For instance, in an unbiased simulation (or as the bias factor and hill height go to 0) - the bias based estimator would be trivially wrong while the histogram would be more reliable.

Thus, agreement between sumhills (which uses bias) and reweighting (which uses histograms and bias) is a necessary thought not sufficient test for ascertaining convergence at any point of time in the simulation. It should work in non-tempered metadynamics as well though the profiles would be choppier. I have not yet thought carefully about applicability to multiple walker, bias exchange or parallel bias schemes.

--

Pratyush Tiwary

Assistant Professor, Department of Chemistry and Biochemistry & Institute for Physical Science and Technology 

University of Maryland, College Park, MD 20742

sites.google.com/site/pratyushtiwary/

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/CAGBhpifnywB-RghuHxXfDK0_9uVd3%2BirJdhY8UFMuom8tRythw%40mail.gmail.com.