GROMACS version: 2021
GROMACS modification: Yes/No
Hi everybody,

I performed an umbrella sampling using as CV the rotation angle of protein helix. But I didn’t use pulling md to produce configurations. I used enforced rotation to generate the configurations. And my question is: if can I use WHAM analysis tool with the output files generated from enforced rotation as input files?

Note: To perform umbrella sampling I based on gromacs tutorial of umbrella sampling and performed at each configuration the enforced rotation with 200kj/mol of applied force and 0 degrees/ps of angle rate. So the helix will not rotate at the simulations.

No, I didn’t receive any answer, but I tried use Alan Grossfield’s WHAM software but I failed because the files with the rotation angle have only zeros. I don’t know if this were possible to do. Now, I’ m trying with an alternative CV with metadynamics with plumed. Here, I left the link. https://groups.google.com/g/plumed-users/c/JoV6hZrFTwk/m/IVSEh0K1BwAJ

Hi, i tried enforced rotation with rot-fit-method0 = rmsd ; (rmsd, norm, or potential).
It does not give rotation angle as 0 in output files but fluctuating about 0 (method=potential gives all zeros). So, I’m currently trying umbrella sampling with rot-fit-method0 = rmsd and then will use Alan Grossfield’s WHAM to get PMF. I’ll let you know if it works.

Hi, I had been using rmsd fit method, but in somes simulations at rotation angle of 100º or upper, I had only zeros in rotation angles and in column fit-angle values with variation of 180º. My idea was if the variation of fit-angle were around zero, I could sum the inital rotation angle (10º, … , 90, … , 180º ) to create the histograms and then calculate WHAM with Alan Grossfield softwhare, but I couldn’t. But after you wrote me these, I tried to use potential fit method and I had only zeros in rotation angle column and values around zero in the fit-angle column. And now, I could sum the initial rotation angle and got the histograms. Thanks to your help I could do this when I’ ve given up with to do umbrella sampling using enforced rotation and I was abandoning one of my PhD thesis objetives. I hope this could help you with your work.

Hi, sorry for I didn’t answer earlier. I had been encreasing the sampling so I couldn’t perform WHAM before that. I performed WHAM and the PMF curve makes sense, but the magnitude values are slower than I have expected. Maybe I used the wrong force constant units. Were you able to perfom WHAM? if yes what constant units did you use? Or are you having the same problem? Let me know if you can.

One thing you might want to consider is that the angle you get from enforced rotations logs is not exactly the same angle that would produce the corresponding bias within a harmonic potential, due to how the module is implemented (at least for some variants).

We were facing this issue and opted for mBAR, which should give a much more correct estimate in this case.

Hi, thank you for the advice. I have been reading the umbrella-sampling-fes example in pymbar github project and I have a question, how did you provide springs constants units using pymbar? Could you help me with that?

There are no spring constants, as I’m talking about doing mBAR, and not WHAM. You can calculate the value of the bias for every frame in every Hamiltonian (= every US window), and you can get the optimized weight of each frame. Then you can rebuild your probability distribution with the frames’ weights, and Boltzmann-invert to obtain the PMF.

Hello, thank you for your reply. So, as I understand I have to use the potential energy from each US window to calculate the reduced potential and use it to perform mBAR. My only question is: at the energy module when I calculate the potential energy, does it include the bias potential energy?

I don’t remember the technical details now, but one of the extra output files should report the energy of the external potential.

If you’re in doubt about it, there’s also the generic option of doing two reruns, one with a .tpr including the original enforced rotation settings, and one without this module.

I verified the potential energy and the COM-pulling-energy and the last has the same value of the energy reported on the output file. But I am bit confused with the pymbar documentation and the umbrella sampling example at github. In the last to initialize FES object they use a matrix KxKxN where K is the numbre of ‘windows’ and N the numbrer of uncorrelated snapshots from each ‘window’, but in the documentation the matrix to be used is KxN. Do you remember how you did the calculation?

The matrix of reduced potentials that you have to provide was either KxN or NxK (not sure which), with K the number of windows and N the total number of samples, aggregated from all windows.

I recall there was a separate vector that specified which sample came from which thermodynamic state, but that might be obsolete for some applications.