- Contemporary messages sorted: [ by date ] [ by thread ] [ by subject ] [ by author ] [ by messages with attachments ]

From: Ross Walker <ross.rosswalker.co.uk>

Date: Thu, 28 Sep 2006 08:29:29 -0700

Dear Pietro

*> 1) after a careful reading of both AMBER9 manuals and Ross Walker's
*

*> Tutorial 7 from "Amber Workshop Tutorials" I understood that
*

*> NEB approach
*

*> can be performed using multiple frames from an approximate path as
*

*> starting guess for min. energy paths. In the case of our simulated
*

*> molecule, we obtained the desired conformational transition during a
*

*> standard MD simulation using implicit GB/SA solvation. Now,
*

*> we'd like to
*

*> obtain a NEB path using frames from this MD simulation. What
*

*> MD protocol
*

*> (in particular, temperature vs. time profile and frequency of frame
*

*> sampling from original MD run) is suggested in this case?
*

This is the million dollar question. The problem is that the simulated

annealing version of NEB in Amber 9 is new and because of this should be

considered still a prototype method. That is I don't think we have a full

understanding yet of the sampling issues involved or what is the best

protocol to use. I don't think you can do much damage by over sampling in

your simulated annealing but under sampling (running the simulated annealing

for too short a time) will cause problems. Where the balance between this

lies is very much dependent on your system size, degrees of freedom, number

of images etc etc. Suffice to say that the larger the system is and the more

images you have the longer you will need to run for. How long though would

only be a guess.

I would suggest you start your system at 0K and then do something along the

follow lines. Note the large values of gamma_ln (this is required for NEB)

to give the system viscosity. For a regular MD simulation these values are

orders of magnitude too large:

Step 1 -> Heat linearly from 0 to 300K dt=0.0005, skmin=10, skmax=10,

gamma_ln=100,

nstlim = 40,000 (20ps)

Step 2 -> Heat linearly to 1000K, dt=0.002, skmin=20, skmax=20,

gamma_ln=1000,

nstlim = 25,000 (50ps)

Step 3 -> Equilibrate at 1000K, dt=0.002, skmin=25, skmax=25, gamma_ln=1000,

nstlim = 50,000 (100ps)

Step 4 -> Cool 1000K to 300K linearly, dt=0.002, skmin=25, skmax=25,

gamma_ln=1000,

nstlim = 50,000 (100ps)

Step 5 -> Equilibrate at 300K, dt=0.002, skmin=25, skmax=25, gamma_ln=1000,

nstlim = 50,000 (100ps)

Step 6 -> Cool 300K -> 0K linearly, dt=0.002, skmin=25, skmax=25,

gamma_ln=1000,

nstlim = 50,000 (100ps)

Step 7 -> 0K, dt=0.002, skmin=25, skmax=25, gamma_ln=5000, nstlim = 50,000

(100ps),

vv=1, vfac=0.1, ntwx=25000

I have had success with this protocol. Although the time lengths are really

just a guess. You may want to increase the length of sampling in the middle

if you can afford it. You should also probably run 10 or so independent runs

with different random seeds to test how well your sampling is converged. You

can then see how many different reactions paths you get and compare them.

In terms of building your initial guess at the pathway there ar eno hard and

fast rules. You just have to provide the code with sufficient input

structures. And indeed this is very much a "guess" at the pathway so the

structures really just need to be reasonable representations of points along

the pathway. So I would just pick what you thing are representative of the

transition from your MD simulation and string them together.

You may also want to just try it with half the images piled up on one end

and half on the other end (I.e. no guess at the pathway). We have had some

success with this. You obviously need to equilibrate for longer but if you

are lucky it won't get trapped in a local minimum due to the high

temperature.

*> 2) looking at the source code of AMBER9 implementation of the quenched
*

*> velocity Verlet algorithm used for NEB minimization (see below), we
*

*> observed that, differently from what stated in the manual,
*

*> when the dot
*

*> product v dot f<0 , velocities are simply zeroed and not
*

*> scaled according
*

*> to eq. (6.34) v=x(v dot f)f and, consequently, the input
*

*> scaling factor
*

*> VFACT is never used (or, equivalently, the program behaves as if VFACT
*

*> value were hardwired to 0). I've found no explanation for
*

*> this discrepance
*

*> either in docs, or in man. upgrades, or in mailing list, or
*

*> in bugfixes.
*

*>
*

*> >>>> File: runmd.f
*

*> >>>> Subroutine: quench
*

*> > if (dotproduct>0.0d0) then
*

*> > v(1:3*natom) = dotproduct*f(1:3*natom)*force
*

*> > else
*

*> > v(1:3*natom) = 0.0d0
*

*> > end if
*

This is a typo in the manual. It should state:

if(v.f>=0) then v = (v.f)f

if(v.f<0) then v = 0

This is what the code does and is what should happen. I will update the

manual errata page.

Note you may also want to read the following paper describing the NEB

implementation in Amber 9:

Mathews, D.H., Case, D.A. J. Mol. Biol. (2006) 357, 1683-1693

All the best

Ross

/\

\/

|\oss Walker

| HPC Consultant and Staff Scientist |

| San Diego Supercomputer Center |

| Tel: +1 858 822 0854 | EMail:- ross.rosswalker.co.uk |

| http://www.rosswalker.co.uk | PGP Key available on request |

Note: Electronic Mail is not secure, has no guarantee of delivery, may not

be read every day, and should not be used for urgent or sensitive issues.

-----------------------------------------------------------------------

The AMBER Mail Reflector

To post, send mail to amber.scripps.edu

To unsubscribe, send "unsubscribe amber" to majordomo.scripps.edu

Received on Sun Oct 01 2006 - 06:07:12 PDT

Date: Thu, 28 Sep 2006 08:29:29 -0700

Dear Pietro

This is the million dollar question. The problem is that the simulated

annealing version of NEB in Amber 9 is new and because of this should be

considered still a prototype method. That is I don't think we have a full

understanding yet of the sampling issues involved or what is the best

protocol to use. I don't think you can do much damage by over sampling in

your simulated annealing but under sampling (running the simulated annealing

for too short a time) will cause problems. Where the balance between this

lies is very much dependent on your system size, degrees of freedom, number

of images etc etc. Suffice to say that the larger the system is and the more

images you have the longer you will need to run for. How long though would

only be a guess.

I would suggest you start your system at 0K and then do something along the

follow lines. Note the large values of gamma_ln (this is required for NEB)

to give the system viscosity. For a regular MD simulation these values are

orders of magnitude too large:

Step 1 -> Heat linearly from 0 to 300K dt=0.0005, skmin=10, skmax=10,

gamma_ln=100,

nstlim = 40,000 (20ps)

Step 2 -> Heat linearly to 1000K, dt=0.002, skmin=20, skmax=20,

gamma_ln=1000,

nstlim = 25,000 (50ps)

Step 3 -> Equilibrate at 1000K, dt=0.002, skmin=25, skmax=25, gamma_ln=1000,

nstlim = 50,000 (100ps)

Step 4 -> Cool 1000K to 300K linearly, dt=0.002, skmin=25, skmax=25,

gamma_ln=1000,

nstlim = 50,000 (100ps)

Step 5 -> Equilibrate at 300K, dt=0.002, skmin=25, skmax=25, gamma_ln=1000,

nstlim = 50,000 (100ps)

Step 6 -> Cool 300K -> 0K linearly, dt=0.002, skmin=25, skmax=25,

gamma_ln=1000,

nstlim = 50,000 (100ps)

Step 7 -> 0K, dt=0.002, skmin=25, skmax=25, gamma_ln=5000, nstlim = 50,000

(100ps),

vv=1, vfac=0.1, ntwx=25000

I have had success with this protocol. Although the time lengths are really

just a guess. You may want to increase the length of sampling in the middle

if you can afford it. You should also probably run 10 or so independent runs

with different random seeds to test how well your sampling is converged. You

can then see how many different reactions paths you get and compare them.

In terms of building your initial guess at the pathway there ar eno hard and

fast rules. You just have to provide the code with sufficient input

structures. And indeed this is very much a "guess" at the pathway so the

structures really just need to be reasonable representations of points along

the pathway. So I would just pick what you thing are representative of the

transition from your MD simulation and string them together.

You may also want to just try it with half the images piled up on one end

and half on the other end (I.e. no guess at the pathway). We have had some

success with this. You obviously need to equilibrate for longer but if you

are lucky it won't get trapped in a local minimum due to the high

temperature.

This is a typo in the manual. It should state:

if(v.f>=0) then v = (v.f)f

if(v.f<0) then v = 0

This is what the code does and is what should happen. I will update the

manual errata page.

Note you may also want to read the following paper describing the NEB

implementation in Amber 9:

Mathews, D.H., Case, D.A. J. Mol. Biol. (2006) 357, 1683-1693

All the best

Ross

/\

\/

|\oss Walker

| HPC Consultant and Staff Scientist |

| San Diego Supercomputer Center |

| Tel: +1 858 822 0854 | EMail:- ross.rosswalker.co.uk |

| http://www.rosswalker.co.uk | PGP Key available on request |

Note: Electronic Mail is not secure, has no guarantee of delivery, may not

be read every day, and should not be used for urgent or sensitive issues.

-----------------------------------------------------------------------

The AMBER Mail Reflector

To post, send mail to amber.scripps.edu

To unsubscribe, send "unsubscribe amber" to majordomo.scripps.edu

Received on Sun Oct 01 2006 - 06:07:12 PDT

Custom Search