Speed up A to B transition sampling using PATH CV

389 views
Skip to first unread message

tLaw

unread,
Aug 10, 2023, 8:46:20 AM8/10/23
to PLUMED users
Dear Plumed users,

I am performing well-tempered metadynamics on a membrane transporter trying to sample the transition between the outward-facing (OF) and inward-facing (IF) configurations using PATH collective variable.
In particular, I am interested in calculating the energy barrier between the two states when the system passes through the "intermediate" state.

The transporter has an alternating gating mechanism ("rocker-switch") characterised by the rotation of two helix-bundle subunits (of2if.png)



of2if.png

I have generated a path of 10 points by linear interpolation of the OF and IF structures (available as X-ray PDBs) previously thermalized by standard MD and extracted the recommended lambda value (path.png).

path.png

By calculating PATH CVs on the unbaiased simulations I also got an estimate of s and z fluctuations to determine gaussian widths.

unbiased.png

The simulation was set up in the NVT ensemble using the following configuration. To note, I am using amber22 with ff19sb, OPC water and Langevin thermostat.

plumed.dat
----------

RESTART
p: PATHMSD REFERENCE=path.pdb LAMBDA=1001
METAD ...
    LABEL=r
    ARG=p.sss
    SIGMA=0.3
    HEIGHT=4.184*0.3
    BIASFACTOR=20
    TEMP=300.0
    PACE=500
    FILE=HILLS
    GRID_MIN=1
    GRID_MAX=10
    CALC_RCT
... METAD
uwall: UPPER_WALLS ARG=p.sss AT=9.8 KAPPA=2000 EXP=2 EPS=0.1 OFFSET=0.0
lwall: LOWER_WALLS ARG=p.sss AT=1.2 KAPPA=2000.0 EXP=2 EPS=0.1 OFFSET=0.0
PRINT ARG=p.sss,p.zzz,r.bias,r.rbias,r.rct,uwall.bias,lwall.bias STRIDE=500 FILE=colvar #every ps

plumed.rew.dat
--------------

p: READ FILE=colvar IGNORE_TIME VALUES=p.sss
m: READ FILE=colvar IGNORE_TIME VALUES=r.rbia
u: READ FILE=colvar IGNORE_TIME VALUES=uwall.
l: READ FILE=colvar IGNORE_TIME VALUES=lwall.
w: REWEIGHT_BIAS TEMP=300 ARG=m.rbias,u.bias,l.bias
h: HISTOGRAM ...
  ARG=p.sss
  GRID_MIN=1
  GRID_MAX=10
  GRID_BIN=50
  #BANDWIDTH=0.3
  KERNEL=DISCRETE
  LOGWEIGHTS=w
...
DUMPGRID GRID=h FILE=histo FMT=%24.16e 

I am controlling the state of the simulation using sum_hills and rct after skipping several picoseconds (hills.png, rct_rew.png)

hills.png

rct_rew.png

The problem I am facing is that the I cannot get the system to sample the transition enough times and I am concerned that my simulation is not converging even after long simulation times (should I even bother trying to perform block analysis?).
How can I improve my setup to better sample the transition between OF and IF states?

Thank you all for your attention :)

Best,
Tommaso

Davide Branduardi

unread,
Aug 11, 2023, 10:06:32 AM8/11/23
to plumed...@googlegroups.com
Ciao Tommaso

I think what you did is a great "first pass". 
You provided a completely fake interpolated initial path and still you got a transition. Well done you! :)
So now I would assume that during that simulation the system displayed the real path it wants to ride on so perhaps you can rechoose the frames from the set of snapshots you got from 200 to 400 ns and you may have more luck. 
What I fear is that even if you are able to cross the path with the interpolation once you end in state B and the system attempts to go back, the gradient is not right (because it comes from a fake interpolation) and gets stuck in a dead end which overlays in S/Z with your input path (especially because at Z~3Ang there may be  lots of accessible states with same S/Z): this messes up your free energy. 

HTH

D


--
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/E2A1515D-4034-4BCC-AC04-E108F97F0528%40gmail.com.

tLaw

unread,
Aug 11, 2023, 11:19:09 AM8/11/23
to plumed...@googlegroups.com
Ciao Davide,

Thank you for your prompt and kind reply! I have been struggling with this for a few months but I decided to reach out to the experts after countless attempts :)

So, I have reweighed the simulation on s and z using this input file:

p: READ FILE=colvar IGNORE_TIME VALUES=p.sss,p.zzz
m: READ FILE=colvar IGNORE_TIME VALUES=r.rbias
u: READ FILE=colvar IGNORE_TIME VALUES=uwall.bias 
l: READ FILE=colvar IGNORE_TIME VALUES=lwall.bias
w: REWEIGHT_BIAS TEMP=300 ARG=m.rbias,u.bias,l.bias
h: HISTOGRAM ...
  ARG=p.sss,p.zzz
  GRID_MIN=1,0
  GRID_MAX=10,0.05
  GRID_BIN=100,100
  #BANDWIDTH=0.3
  KERNEL=DISCRETE
  LOGWEIGHTS=w
...
f: CONVERT_TO_FES GRID=h TEMP=300
DUMPGRID GRID=f FILE=fes2d FMT=%24.16e 

And extrapolated the frames lying on the minimum-energy-path (found using Dijkstra algorithm).

mep_fes2d.png

Although the FES is still rough, my guess would be to extract 10 frames from this collection, make them equidistant using path_tools and re-run the simulation, does this sound right?

Thanks again,
Tommaso

tLaw

unread,
Aug 11, 2023, 11:34:09 AM8/11/23
to PLUMED users
Adding to my last reply, I realised I am running into the same problem of ill-definedness as there may be multiple structures at a given (s, z) when z is large. I could potentially calculate an RMSD matrix on all frames and get the medoid around each point of the path I could sample so far. I could also be over-thinking this as, anyway, I cannot trust my FES yet. So probably extracting 10 frames equidistant on S from the 200-400ns interval would do just fine at this stage.

Davide Branduardi

unread,
Aug 12, 2023, 2:43:40 AM8/12/23
to plumed...@googlegroups.com
Ciao Tommaso

that sounds sensible but a bit an overkill as the fes and bias you got is not exactly right. anyway a better job is better than a worse one so it all sounds good!
D

--
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.

Davide Branduardi

unread,
Aug 12, 2023, 2:46:56 AM8/12/23
to plumed...@googlegroups.com

another trick to avoid the degeneration at high z is running with a wall in z so you limit the sampling. the drawback there is that 1) if the path is fake there may not be any physical path at low z or may be high energy - a gentle steered md rather than full metà may ensure that the transition to occur in a reasonable computing time 2) you may need multiple iterations for the path to properly relax

good luck with it. it looks like an exciting project and you are on the right path (pun intended! :-)

D

--
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.

Lv Wade

unread,
Oct 10, 2023, 5:34:37 AM10/10/23
to PLUMED users
Dear Davide, 
    I'm working on similar problem and read your suggestions in above replies. I'm a bit confused  on this sentence: "a gentle steered md rather than full metà may ensure that the transition to occur in a reasonable computing time ". Does this means that one have to regenerate the path reference nodes before running metadynamics on PCV?
Best,
Wade

Davide Branduardi

unread,
Oct 10, 2023, 1:21:40 PM10/10/23
to plumed...@googlegroups.com
Yes. Without reference conformations you can't provide an input for the path cvs. 
But there are plenty of ways to do that, depending on the problem. 
Sometimes even a geometrically interpolated path may work but my experience is that it is in general a bad approximation which requires some refinement anyways. 
So my approach is generally have a way to "stimulate" the transition I am interested in so to have a path which is "as realistic as possible" from which I can calculate the free energy. 
HTH

D

Lv Wade

unread,
Oct 11, 2023, 10:20:25 PM10/11/23
to PLUMED users
Many thanks Davide!
Reply all
Reply to author
Forward
0 new messages