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

From: Jason Swails <jason.swails.gmail.com>

Date: Thu, 4 Jul 2013 10:20:00 -0400

On Wed, Jul 3, 2013 at 8:47 PM, Biao Ma <jackyma1981.gmail.com> wrote:

*> Hello, amber users,
*

*>
*

*> Now, I am reading a paper that is "Insights into structural properties of
*

*> denatured human prion 121-230 at melting temperature studied by replica
*

*> exchange molecular dynamics". Base on the describe of paper, I write the
*

*> process and the input file of amber simulation, but I cannot sure that is
*

*> correct.****
*

*>
*

*> Who can help me check it?****
*

*>
*

I'll offer some help and comments, but on the whole there is too much here

to offer an in-depth evaluation of your procedure in a reasonable amount of

time. In general, the best way to validate a specific scheme is to try it

on some systems, especially ones where you already know what the answer

should be.

*> **
*

*> Thank a lot.
*

*> Biao Ma
*

*> **
*

*> pumed url of the paper: http://www.ncbi.nlm.nih.gov/pubmed/22339436****
*

*>
*

*> ** **
*

*>
*

*> Blow is a part of method in the paper:****
*

*>
*

*> =====start
*

*> In this work, the temperature setting used for sheep PrP 125-230 such that
*

*> MD simulations were exchanged at 320.0, 322.0,324.0, 326.0, 328.1, 330.2,
*

*> 332.3, 334.4, 336.5, 338.6, 340.8, 343.0, 345.2, 347.4, 349.6, 351.8,
*

*> 354.0, 356.2, 358.5, 360.8, 363.1, 365.4, 367.7, and 370.0 K was used. All
*

*>
*

These temperatures appear to be evenly spaced. I would suggest using a

geometric temperature series to optimize the performance of replica

exchange. I have found this server [http://folding.bmc.uu.se/remd/] to be

remarkably accurate.

*> of the replicas were equilibrated for 20 ns without exchanging temperatures
*

*> and then extended for 65 ns of REMD simulation. The generalized Born model
*

*> used in this study modified the calculation of Born radii and improved the
*

*> accuracy in the solvent polarization for macromolecules. The combinational
*

*> use of the all-atom point-charge force-field (also known as ff03) and the
*

*> generalized Born model led to successful folding of several proteins. The
*

*> AMBER 11 simulation package 26 was used in both REMD simulation and data
*

*> analysis. The melted huPrP 121-230 was computed starting from an extended
*

*> huPrP. To generate the initial extended structure, a heating method was
*

*> applied to a known NMR structure (PDBcode 1hjn,15Figure 1A), enabling it to
*

*> unfold at 600 K for 40 ns of MD simulation to result in an extended
*

*> conformation (Figure 1B) as described previously.****
*

*>
*

*> During this simulation, the disulfide covalent bond between residues 179
*

*> and 214 was preserved. In total, 24 replicas with duration of 65 ns and
*

*> with an integration time step of 2 fs were computed based on the extended
*

*> huPrP with different random number seeds to generate the initial
*

*> conditions. A 16 Å force-shifted non-bonded cutoff and generalized Born
*

*> solvent models with salt concentration of 0.2 M were applied.****
*

*>
*

*> **=====end
*

*> **
*

*>
*

*> **
*

*>
*

*> **
*

*>
*

*> Simulation Procedure:****
*

*>
*

*> 1. system building****
*

*>
*

*> 2. system minimization****
*

*>
*

*> 3. heating system****
*

*>
*

*> 4. generate the extend conformation****
*

*>
*

*> 5. local minimization after heating system****
*

*>
*

*> 6. equilibrate the every replica****
*

*>
*

*> 7. REMD simulation****
*

*>
*

*> ** **
*

*>
*

*> ** **
*

*>
*

*> *1.leap.inp for system building*
*

*>
*

*> ** **
*

*>
*

*> source pdb: 1ag2****
*

*>
*

*> use the the ff03 (Duan et al.) force field****
*

*>
*

*> ** **
*

*>
*

*> leap.inp
*

*>
*

*> source leaprc.ff03.r1****
*

*>
*

*> set default PBradii mbondi2****
*

*>
*

*> # load pdb file****
*

*>
*

*> 1ag2 = loadPdb input.pdb****
*

*>
*

*> savePdb 1ag2 1ag2.pdb****
*

*>
*

*> bond 1ag2.179.SG <http://1ag2.179.sg/> 1ag2.214.SG <http://1ag2.214.sg/
*

*> >****
*

*>
*

*> # save 1ag2 to prmtop and inpcrd files****
*

*>
*

*> saveAmberParm 1ag2 1ag2.prmtop 1ag2.inpcrd****
*

*>
*

*> # finish****
*

*>
*

*> quit****
*

*>
*

*> ** **
*

*>
*

*> 2. system minimization****
*

*>
*

*> minimisation for heated system****
*

*>
*

*> &cntrl****
*

*>
*

*> imin=1, maxcyc=1000, ncyc=500,****
*

*>
*

*> igb=5, ntb=0,****
*

*>
*

*> cut = 16, rgbmax = 16, saltcon = 0.2,****
*

*>
*

This cutoff is too small IMO. I've seen people use a cutoff of 30 Å when

employing a salt concentration like the one you have, but I would suggest

using an infinite cutoff with a GB simulation.

*>
*

*> ntpr=100****
*

*>
*

*> /****
*

*>
*

*> ~****
*

*>
*

*> ~****
*

*>
*

*> ** **
*

*>
*

*> *3. heat the system*
*

*>
*

*> heating system from 0 K to 600K during 100 ps.****
*

*>
*

*> &cntrl****
*

*>
*

*> nstlim = 50000, dt = 0.002,****
*

*>
*

At 600 K, a 2 fs time step may be too short.

*> ntt = 1, tautp = 1.0,****
*

*>
*

I would suggest using Langevin dynamics to maintain your temperature,

since it is more theoretically rigorous than the Berendsen thermostat.

ntt = 3, gamma_ln=1.0, ig=-1,

*> tempi = 0, temp0 = 600, ntc =2, ntf = 2,****
*

*>
*

*> ntpr =100, ntwx = 100,****
*

*>
*

*> ntb = 0, igb = 5, ****
*

*>
*

*> cut = 16, rgbmax = 16, saltcon = 0.2,****
*

*>
*

*> nmropt = 1,****
*

*>
*

*> /****
*

*>
*

*> &wt****
*

*>
*

*> type = 'TEMP0', istep1 = 0, istep2 =50000, value1 = 0, value2=600,****
*

*>
*

*> /****
*

*>
*

*> &wt type = 'END'****
*

*>
*

*> /****
*

*>
*

*> ** **
*

*>
*

*> *4**. **generate the full unfolded conformation*
*

*>
*

*> ** **
*

*>
*

*> To generate the initial extend structure, a heating method was used to a
*

*> known NMR structure (PDB code:1ag2), enabling it to unfold at 600 K for 40
*

*> ns of LD simulation to result in an extended conformation.****
*

*>
*

*> here, I do not sure whether I do the restart MD(set irest = 1, ntx = 5, or
*

*> irest = 0, ntx = 1), Simulation results is difference, which setting should
*

*> I use?
*

*>
*

Since the last simulation ended at 600 K, you should just restart. The

results will actually be different between any two runs (generally, unless

you run 2 identical sets of inputs on the same model GPU).

*>
*

*> 40nsld.inp****
*

*>
*

*> enabling the heated NMR structure to unfold at 600 K for 40ns of LD
*

*> simulation****
*

*>
*

*> &cntrl****
*

*>
*

*> irest = 1, ntx = 5,****
*

*>
*

*> nstlim = 20000000, dt = 0.002,****
*

*>
*

*> ntt = 3, gamma_ln = 1.0,****
*

*>
*

Make sure you set an explicit value for ig (the random number seed). It

should be different for every simulation or you will begin to see

synchronization artifacts. http://pubs.acs.org/doi/abs/10.1021/ct800573m

Setting ig=-1 will pull the random seed from the wall clock, ensuring that

they are different for each run.

*>
*

*> tempi = 600,temp0 = 600,****
*

*>
*

*> ntb = 0, igb = 5,****
*

*>
*

*> ntpr = 500, ntwx = 1000, ntwr = 2000000,****
*

*>
*

*> ntc = 2, ntf = 2,****
*

*>
*

*> cut = 16, rgbmax = 16, saltcon = 0.2,****
*

*>
*

*> /
*

*>
*

*>
*

*> 600kmd.inp****
*

*>
*

*> enabling the heated NMR structure to unfold at 600 K for 40ns of LD
*

*> simulation****
*

*>
*

*> &cntrl****
*

*>
*

*> irest = 0, ntx = 1,****
*

*>
*

*> nstlim = 20000000, dt = 0.002,****
*

*>
*

*> ntt = 3, gamma_ln = 1.0,****
*

*>
*

*> tempi = 600,temp0 = 600,****
*

*>
*

*> ntb = 0, igb = 2,****
*

*>
*

*> ntpr = 500, ntwx = 1000, ntwr = 2000000,****
*

*>
*

*> ntc = 2, ntf = 2,****
*

*>
*

*> cut = 16, rgbmax = 16, saltcon = 0.2,****
*

*>
*

*> /
*

*>
*

*> 40nsld.out irest = 1, ntx = 5,****
*

*>
*

*> 600kmd.out irest = 0, ntx = 1****
*

*>
*

*> [image: 40nsld.etot.png]****
*

*>
*

*> [image: 600kmd.etot.png]
*

*>
*

*>
*

*> ** **
*

*>
*

*> 5. local minimization after heating system****
*

*>
*

*> use the mdin file which is same to step 3 above。****
*

*>
*

*> ** **
*

*>
*

*> 6. equilibrate the every replica****
*

*>
*

*> ** **
*

*>
*

*> equilibrate.mdin****
*

*>
*

*> ** **
*

*>
*

*> equilibration 20 ns, every 10ps save output. <- Is this must to do
*

*> equilibrate for 20 ns?**
*

*>
*

*> equilibration****
*

*>
*

*> &cntrl****
*

*>
*

*> irest=0, ntx=1,****
*

*>
*

*> nstlim=10000000, dt=0.002,****
*

*>
*

*> irest=0, ntt=3, gamma_ln=1.0,****
*

*>
*

*> temp0=XXXXX, ig=RANDOM_NUMBER,****
*

*>
*

You can just use ig=-1.

*>
*

*> ntc=2, ntf=2, nscm=1000,****
*

*>
*

*> ntb=0, igb=5,****
*

*>
*

*> cut = 16, rgbmax = 16, saltcon = 0.2,****
*

*>
*

*> ntpr=5000, ntwx=5000, ntwr=10000000,****
*

*>
*

*> nmropt=1,****
*

*>
*

*> /****
*

*>
*

*> &wt TYPE='END'****
*

*>
*

*> /****
*

*>
*

*> DISANG=system_chir.dat****
*

*>
*

*> ** **
*

*>
*

*> 7. REMD simulation****
*

*>
*

*> remd 65ns exchange every 2ps <- Here I have been confused, how often to
*

*> exchange is more appropriate ?****
*

*>
*

For T-REMD, you don't lose any efficiency by attempting exchanges more

often. I would attempt exchanges at least every 100 steps, although we've

shown that convergence gets better for even more rapid exchange attempts.

*>
*

*> ** **
*

*>
*

*> remd.mdin****
*

*>
*

*> remd 65ns exchange every 2ps****
*

*>
*

*> &cntrl****
*

*>
*

*> irest=0, ntx=1,****
*

*>
*

*> nstlim=1000, dt=0.002,****
*

*>
*

*> irest=0, ntt=3, gamma_ln=1.0,****
*

*>
*

*> temp0=XXXXX, ig=RANDOM_NUMBER,****
*

*>
*

*> ntc=2, ntf=2, nscm=1000,****
*

*>
*

*> ntb=0, igb=5,****
*

*>
*

*> cut = 16, rgbmax = 16, saltcon = 0.2,****
*

*>
*

*> ntpr=100, ntwx=1000, ntwr=100000,****
*

*>
*

*> nmropt=1,****
*

*>
*

*> numexchg=32500,****
*

*>
*

*> /
*

*>
*

HTH,

Jason

Date: Thu, 4 Jul 2013 10:20:00 -0400

On Wed, Jul 3, 2013 at 8:47 PM, Biao Ma <jackyma1981.gmail.com> wrote:

I'll offer some help and comments, but on the whole there is too much here

to offer an in-depth evaluation of your procedure in a reasonable amount of

time. In general, the best way to validate a specific scheme is to try it

on some systems, especially ones where you already know what the answer

should be.

These temperatures appear to be evenly spaced. I would suggest using a

geometric temperature series to optimize the performance of replica

exchange. I have found this server [http://folding.bmc.uu.se/remd/] to be

remarkably accurate.

This cutoff is too small IMO. I've seen people use a cutoff of 30 Å when

employing a salt concentration like the one you have, but I would suggest

using an infinite cutoff with a GB simulation.

At 600 K, a 2 fs time step may be too short.

I would suggest using Langevin dynamics to maintain your temperature,

since it is more theoretically rigorous than the Berendsen thermostat.

ntt = 3, gamma_ln=1.0, ig=-1,

Since the last simulation ended at 600 K, you should just restart. The

results will actually be different between any two runs (generally, unless

you run 2 identical sets of inputs on the same model GPU).

Make sure you set an explicit value for ig (the random number seed). It

should be different for every simulation or you will begin to see

synchronization artifacts. http://pubs.acs.org/doi/abs/10.1021/ct800573m

Setting ig=-1 will pull the random seed from the wall clock, ensuring that

they are different for each run.

You can just use ig=-1.

For T-REMD, you don't lose any efficiency by attempting exchanges more

often. I would attempt exchanges at least every 100 steps, although we've

shown that convergence gets better for even more rapid exchange attempts.

HTH,

Jason

-- Jason M. Swails Quantum Theory Project, University of Florida Ph.D. Candidate 352-392-4032 _______________________________________________ AMBER mailing list AMBER.ambermd.org http://lists.ambermd.org/mailman/listinfo/amberReceived on Thu Jul 04 2013 - 07:30:02 PDT

Custom Search