Peptide Folding with Temperature Replica Exchange Simulations

Objective and Overview

The objective of this tutorial is to learn how to use the MMTSB Tool Set in combination with CHARMM to carry out temperature replica exchange simulations.


In this tutorial we will fold the alpha-helical peptide (AAQAA)3 from an extended structure using implicit solvent and temperature replica exchange sampling to accelerate the crossing of kinetic barriers along the folding pathway.

1. Initial setup

We will begin with an extended structure for the sequence AAQAAAAQAAAAQAA that is available here.

Before starting the replica exchange simulation we need to briefly minimize the structure:
    
  minCHARMM.pl -par gb=bgb,blocked,param=19 aaqaa.extended.pdb > aaqaa.ext.min.pdb
Because we will use implicit solvent no further equilibration is necessary.

2. Replica exchange run

Replica exchange simulations involve multiple parallel molecular dynamics simulations that fluctuate in temperature according to the replica exchange criterion. A key input parameter for a replica exchange simulation is the number of replicas and range of temperatures to be covered by those replicas. The temperature difference between neighboring replicas needs to be tuned so that 20-80% of the attempted exchanges between neighbors are successful. Typically, the temperatures are spaced exponentially to achieve uniform acceptance ratios. For a small peptide in implicit solvent, eight replicas are sufficient to span a temperature range of 300 to 500K with exchange acceptance probabilities of about 40%. For larger systems and/or explicit solvent environments many more replicas are needed as a result of the larger number of degrees of freedom.

The second set of parameters unique to replica exchange simulations is how long constant temperature molecular dynamics simulations are run during each cycle before an exchange is attempted and how many of such cycles are run. Typically, exchanges are attempted every 0.1 to 2 ps, which corresponds to 200 to 1000 simulation steps with a typical time step of 2 fs, and it is often desirable to run nanoseconds or tens of nanoseconds per replica to achieve sufficient sampling.

Finally, replica exchange simulations are ideal for parallel computing platforms where they provide near-ideal scaling due to low communication overhead. Typically, each replica may be run in serial on a single core/CPU, but it is also possible to use multiple cores/CPUs for each replica if sufficiently large resources are available.

A replica exchange simulation for the extended AAQAA peptide is started with the following command:
    
  aarex.pl -temp 8:275:500 -n 1000 -mdpar gb=bgb,dynsteps=200,blocked,param=19 \
           -par archive,savedatafreq=100 -dir rex aaqaa.ext.min.pdb
In this command -temp specifies the number of replicas (8) and range of temperatures (275-500K), -n specifies the number of replica exchange cycles (1000), -mdpar provides the CHARMM simulation parameters (use of GB implicit solvent model, 200 time steps = 0.4 ps per replica exchange cycle, blocked termini, and use of the united-residue CHARMM19 force field), -par specifies additional replica exchange options (save replice exchange trajectories as archive files, save output data every 100 cycles), and the directory where all the output from the simulation is written is given with -dir. Finally, the starting conformation for each of the eight simulations is given as the last argument.

The previous command will work if you are on a computer with 8 or more cores/CPUs. If you have less than 8 CPUs on a single computer, you can use multiple nodes by providing a list of nodes with a file. As an example, you might have access to four computers (snoopy, goofy, clifford, scoobydoo) each with two CPUs. Then you would create a file called hostfile with the following contents:
    
  snoopy 2
  goofy 2
  clifford 2
  scoobydoo 2
The host file is then given to aarex.pl with the extra option
    
  -hosts hostfile
With that option the two replicas would be run on each of the four machines. This requires that you can login to each of the machines via ssh without being prompted for a password (by using ssh keys without a passphrase) and that you share a common directory where the data is written. It is also possible to run on machines that do not share common directories via the -mp option. See the aarex.pl documentation for more details on this option.

Once you have determined what the best option is to run in parallel in your environment, run the above command to complete the first 1000 cycles. It should take about 10 minutes to complete the simulations.

While we will learn more about analysis of replica exchange trajectories in a different tutorial, you may want to explore the output files in the rex subdirectories that were generated as a result of the replica exchange simulation. In this subdirectory you will find a number of data files beginning with rex... that contain various data related to the replica exchange as a whole.

Take a look at rexserver.ninx which contains information about the exchange acceptance fraction in the fifth column with:
    
  cat rexserver.ninx
Are the values between 20 and 80% as expected?

The aaN subdirectories contain the data from each of the replicas. The trajectory is stored in the prod.coor.archive file and the final structures at the end of each replica simulation are stored as final.pdb in each of the replica directories. Take a look at the final structures. Do any of those structures have helical sections?

If the peptide has not folded up, we should continue the simulation. Because all of the simulation parameters are stored, we can restart the replica exchange simulation and continue for another 19000 cycles to reach a total of 20000 cycles with the following command:
    
  aarex.pl -n 19000 -dir rex 
If you had to specify additional options for the parallel execution (such as the -hosts option) those options also have to be repeated when restarting the simulation, but the other options and the initial PDB file do not need to be repeated.

The longer replica exchange run will take more time, on the order of two hours. Once the simulation is completed inspect the final conformations again and continue with the replica exchange analysis tutorial.

Written by M. Feig