Visualizing the path progression

Erika McCarthy¹, Şölen Ekesan¹, and Darrin M. York¹

¹Laboratory for Biomolecular Simulation Research, Institute for Quantitative Biomedicine and Department of Chemistry and Chemical Biology, Rutgers University, Piscataway, NJ 08854, USA

Learning objectives

Visualize the progression of the minimum free energy path on a free energy surface from a 2D reaction within MTR1

Relevant literature

Catalytic mechanism and pH dependence of a methyltransferase ribozyme (MTR1) from computational enzymology
Surface-Accelerated String Method for Locating Minimum Free Energy Paths
Extension of the Variational Free Energy Profile and Multistate Bennett Acceptance Ratio Methods for High-Dimensional Potential of Mean Force Profile Analysis

Make a path progression movie

Finally, let’s take a look at the convergence of the string over all of the iterations of string data provided in the outputs for the pathway you have chosen (linear, PT, or MT). The Example2d.py script has a feature for creating movies of the string and free energy surface progression. Each string directory contains a script called Make2dMovie.sh that runs the Example2d.py script in order to create a series of png files that are combined into a movie using ffmpeg. This exercise was performed on 50 string iterations to confirm convergence of the path. A run of this script has already been performed for you, and the output was written to the fullmovie directory:

[user@cluster] ls fullmovie
movie.sims.mp4  sims.00006.png  sims.00012.png  sims.00018.png  sims.00024.png  sims.00030.png  sims.00036.png  sims.00042.png  sims.00048.png
sims.00001.png  sims.00007.png  sims.00013.png  sims.00019.png  sims.00025.png  sims.00031.png  sims.00037.png  sims.00043.png  sims.00049.png
sims.00002.png  sims.00008.png  sims.00014.png  sims.00020.png  sims.00026.png  sims.00032.png  sims.00038.png  sims.00044.png  sims.00050.png
sims.00003.png  sims.00009.png  sims.00015.png  sims.00021.png  sims.00027.png  sims.00033.png  sims.00039.png  sims.00045.png  sims.00051.png
sims.00004.png  sims.00010.png  sims.00016.png  sims.00022.png  sims.00028.png  sims.00034.png  sims.00040.png  sims.00046.png  sims.00052.png
sims.00005.png  sims.00011.png  sims.00017.png  sims.00023.png  sims.00029.png  sims.00035.png  sims.00041.png  sims.00047.png  sims.00053.png

Transfer the movie.sims.mp4 to your local computer and open the file in google chrome or a supported media player:

[user@local] google-chrome fullmovie/movie.sims.mp4

Here is the movie depicting the string progression from the linear initial guess:

In the movie, you will see small circles along the path that represent the control points defining the path. The x’s indicate where the umbrella centers have been placed for simulations. You will notice that there is alternation between simulations (x’s) being placed explicitly on the path for refinement and off the path for exploration. This allows for rapid convergence as well as rigorous sampling on both sides of the path. Since we used the neqit=4 flag when running ndfes, string iterations 0-4 are eventually phased out once at least 5 previous iterations of sampling available to analyze. The reasoning is that the initial guess is very poor (the string evolves dramatically), so the short amount of sampling of these early iterations likely do not well-represent an equilibrium ensemble. Feel free to look at the movies generated for the other 2 paths to compare.

Visually, there is very little evolution of the minimum free energy path past iteration 10 or so. When the simulation data was generated, the MAKESIMS command in SASM.slurm outputs two convergence tests to analysis.CIT.log, where CIT is the current iteration. The first test is called the slope test that checks whether the slope of the reaction coordinates over at least 5 iterations is less than the threshold value. The RMSD test checks whether the RMSD of the overall path over at least 5 iterations in less the the threshold value. These are thresholds that have been established based on various test cases, but one should still visually inspect the path.

Check the convergence tests using the following commands:

[user@cluster] grep -r "Slopes appear to be converged" *log | head -n5
analysis.008.log:Slopes appear to be converged
analysis.009.log:Slopes appear to be converged
analysis.010.log:Slopes appear to be converged
analysis.011.log:Slopes appear to be converged
analysis.012.log:Slopes appear to be converged

[user@cluster] grep -r "Path is effectively CONVERGED" *log
analysis.010.log:Min slope in RMSD    1.687e-03 Path is effectively CONVERGED
analysis.011.log:Min slope in RMSD    6.127e-04 Path is effectively CONVERGED
analysis.012.log:Min slope in RMSD    2.864e-04 Path is effectively CONVERGED
analysis.013.log:Min slope in RMSD    2.337e-04 Path is effectively CONVERGED
analysis.014.log:Min slope in RMSD   -9.001e-05 Path is effectively CONVERGED

Based on this criteria, it appears that both tests were passed by iteration 10 for the linear initial guess. This is consistent with our visual inspection of the movie. Based on the string simulations performed in this tutorial, it is clear that the minimum free energy path for the MTR1 reaction involves a stepwise mechanism where a general acid protonates the bound ligand to initiate methyl transfer to the target adenine. This free energy barrier (~20 kcal/mol) is about 10 kcal/mol lower than the 1 dimensional case that we computed previously for just methyl transfer, and the 2 dimensional path where methyl transfer precedes proton transfer.