Question about how to set parameters for simulation annealing

GROMACS version: 2020.1
GROMACS modification: Yes/No
Dear Gromacs users,
I want to performing an Annealing Simulation with Molecular Dynamics (AS-MD) with no constrains of a small molecule in water.

The molecule is a nucleoside (approximately 35 atoms). My solute is center in a cubic box with at least 1nm from the edges. I will do the energy minimization and I am choosing the Andersen thermostat to avoid the flying ice cube problem. My simulation time will be 10ns and I will collect frames every 200fs. The simulation step is 1fs.

I want to explore the conformational space of this molecule and achieve the minimum. So After the minimization step I want to perform an AS-MD. But I am straggling to choose the annealing parameters for Gromacs.

Does anybody knows if there is a protocol or methodology to choose these parameters?
The parameters are:

  • Initial temperature
  • Final temperature
  • Annealing-npoints
  • Annealing-ntimes
  • Annealing-temp

Kindly,

Lazaro

Hi Lazaro,

You formulate two different goals within one simulation: exploring conformational space and achieving a minimum energy configuration.

I recommend you do this in a setup where you first explore conformational space in a plain MD simulation as you describe (though I’d recommend v-rescale instead of Anderson, and for almost all force-fields a 2fs time step and collecting frames a bit less frequently (every 2ps should be okay) and a dodecahedron box).

Then, if you are interested in an ensemble of minimum energy configurations, in a second step, you aextract frames from this trajectory with trjconv , using, e.g., -sep and start new simulations where you cool down your system. I reckon than cooling linearly from 300K to 0K in something like 500 ps should be fine.

These annealing simulations will give you final conformations and energies that you can read with gmx energy`. From there you can figure which one of all the final configuration from the annealing simulation had the least energy.

For more info about the annealing parameters have a look at
http://manual.gromacs.org/documentation/current/user-guide/mdp-options.html#simulated-annealing

At the bottom of this paragraph is a more thorough example (After Confused? …)

Dear Cblau

Thank you very much for the reply it is really appreciate it since I am new in Gromacs.

I just want to be sure that I understand: so you are proposing to do first like a nvt equilibration? But then how do I choose these frames? Do I choose the most stables? or the last one after the equilibration?

In the case of choosing the most stables how could I choose them? Using clustering for example or maybe just visualizing the PES (Potential Energy surface) for the system?

In the last step then it would be the AS?

I will follow your advise but do you know if there is any systematic way to choose the initial temperature and cooling parameters?

Kindly,

Lazaro

Hi Lazaro,

The question here is what you want to learn about the system?

The outline I gave you is a simple strategy to try to find a single minimum energy configuration.
What I had in mind, was that you’d run an equilibrium simulation to generate you an ensemble of different structures as a starting point for your cooling simulations.
Then you’d pick a frame every 0.5 ns or so and run a cooling simulation for all them. Then have a look at the energy of the very last frame of all these cooling simulations to identify the structure with the lowest overall energy.

But, as said, this a simple strategy to find the single minimum energy configuration. Most often this is not really what you’re after in MD, so it’d be useful to know what quantities you want to look at overall.

Dear Cblau

Thank you very much for the clarification. Now that I remember I have seen in some papers that people do that kind of procedure to explore the conformational space of small polypeptides. Then at the end they concatenate all the trajectories and do clustering or dihedral Principal Component Analysis to choose the most stable frame.

In my case I am looking for the most stable(s) conformers of my molecule because after I will reoptimized them with DFT, but I am looking for a sampling methodology that is not as expensive as DFT to explore as much as possible all the conformational space of the system.

In an initial text I arrange the parameters for the nvt for 100ps of running collecting every 2ps, but this mean 50 frames. If I run the annealing for 10ns after collecting every 2fs it would mean that it would be equivalent to 500ns simulation (5M frames). Even when my molecule is small and the 10ns runs really fast I think this is a lot.

Would you have any suggestions for the running times for the NVT and the annealing?
I noticed that my system equilibrates after 10ps.

For example could I do NVT at 10ps (is this too short?), collecting every 2ps (5 frames) and then run each at 10 ns (5*100000=500000 frames) and then concatenate all the trajectories and analyze the most stable frames?.

Kindly,

Lazaro