Department of Chemistry, Faculty of Science, University of Kurdistan |
Molecular dynamics simulations using AMBER After setting up the protein and parametrizing hetero groups (cf. AMBER-setup and AMBER-Parm pages) we can start to set up MD simulations. 1. Protonate and solvate the structureRun tleap as the following tleap -s -f tleap.in where the tleap.in file is as: source leaprc.ff14SB source leaprc.gaff source leaprc.water.tip3p addPath /lunarc/nobackup/projects/teobio/Mehdi/glxii/Parms loadAmberPrep zn.in loadAmberPrep xoh.in loadAmberPrep o2x.in loadAmberPrep slg.in loadAmberPrep gsh.in loadAmberPrep gbp.in loadAmberParams zn.dat loadAmberParams xoh.dat loadAmberParams o2x.dat loadAmberParams slg.dat loadAmberParams gsh.dat loadAmberParams gbp.dat x=loadpdb 01-tided_up solvateOct x TIP3PBOX 10 saveamberparm x prmtop prmcrd savepdb x pdbout quit
In this step, note the following points: · Add "TER" between the chains and hetero groups coordinates in the final PDB file (01-tided_up). · If Cys-Cys crosslinks exist in the PDB structure, add extra bonds in the tleap.in file. (e.g., bond x.22.SG x.65.SG; note that the atom names must be capital; also note that tleap renumbers the residues, especially if there are several subunits). · You need to change and sort the atom names of the hetero groups (column 3) as they appear in the .in files. For example, we edited the lines containing GBP as follows. HETATM 4083 C1 GBP B 464 27.140 47.038 60.878 1.00 43.76 C HETATM 4084 O1 GBP B 464 27.407 47.659 61.940 1.00 45.46 O HETATM 4085 O2 GBP B 464 26.979 47.636 59.785 1.00 45.08 O HETATM 4086 N1 GBP B 464 28.149 44.845 60.180 1.00 42.86 N HETATM 4087 C2 GBP B 464 26.999 45.486 60.865 1.00 42.41 C HETATM 4088 C3 GBP B 464 26.858 44.817 62.259 1.00 39.40 C HETATM 4089 C4 GBP B 464 23.765 40.788 63.655 1.00 24.26 C HETATM 4090 S1 GBP B 464 24.903 40.411 62.262 1.00 27.99 S HETATM 4091 C5 GBP B 464 25.645 45.211 63.151 1.00 34.52 C HETATM 4092 C6 GBP B 464 24.354 44.388 63.029 1.00 32.05 C HETATM 4093 O3 GBP B 464 23.381 44.874 62.439 1.00 32.66 O HETATM 4094 N2 GBP B 464 24.355 43.170 63.578 1.00 26.25 N HETATM 4095 C7 GBP B 464 23.215 42.228 63.508 1.00 24.01 C HETATM 4096 C8 GBP B 464 22.119 42.465 64.618 1.00 21.59 C HETATM 4097 O4 GBP B 464 21.071 41.810 64.582 1.00 20.46 O HETATM 4098 N3 GBP B 464 22.477 43.285 65.617 1.00 20.07 N HETATM 4099 C9 GBP B 464 21.677 43.493 66.841 1.00 21.16 C HETATM 4100 C10 GBP B 464 21.340 44.944 67.176 1.00 21.22 C HETATM 4101 O5 GBP B 464 21.720 45.837 66.388 1.00 25.56 O HETATM 4102 O6 GBP B 464 20.693 45.180 68.211 1.00 18.65 O HETATM 4103 N4 GBP B 464 24.505 40.113 59.521 1.00 35.22 N HETATM 4104 O7 GBP B 464 25.829 40.195 59.532 1.00 37.31 O HETATM 4105 C11 GBP B 464 23.951 40.271 60.735 1.00 31.76 C HETATM 4106 O8 GBP B 464 22.741 40.219 60.938 1.00 33.64 O HETATM 4107 C12 GBP B 464 23.856 39.902 58.227 1.00 35.58 C HETATM 4108 C13 GBP B 464 22.626 40.498 57.912 1.00 36.70 C HETATM 4109 C14 GBP B 464 24.500 39.089 57.282 1.00 36.16 C HETATM 4110 C15 GBP B 464 22.032 40.300 56.668 1.00 36.01 C HETATM 4111 C16 GBP B 464 23.909 38.894 56.039 1.00 36.35 C HETATM 4112 C17 GBP B 464 22.689 39.488 55.748 1.00 36.44 C HETATM 4113 Br1 GBP B 464 21.893 39.204 54.066 1.00 35.21 BR
· You may want to add counter ions into your system to neutralize it or add an ionic strength to it. If so, add the following to your tleap.in file after solvateOct. addIonsRand x Na+ 96 Cl- 110 4 And also, add the following at the beginning of the tleap.in file to read the ion parameters. loadAmberParams frcmod.ionsjc_tip3p However, since it replaces any water molecule with a counter ion, the ions may end up inside the protein. This can be avoided by running tleap both without and with counter ions and saving pdbout files of both calculations. Then run changepdb, command fci (fix counter ions) on the pdbout file with counter ions (pdbout-ion), and read in the pdbout file without counter ions. The program overwrites the prmcrd file with counterions that are more than 4 Å from the protein (and also from each other).
changepdb <<EOF pdbout-ion fci prmcrd-ion pdbout prmcrd q EOF \mv prmtop-ion prmtop
2. Treating metal centersIf there are metal centers in your protein, treat them as described in section 3 of the AMBER-Parm page)
3. Run MD simulationsRun MD simulations with the following commands and the sander.in files listed below. The following uses the GPU-accelerated pmemd code of AMBER for the simulations. However, you can instead use the sander module of AmberTools if you have not purchased the AMBER program (replace pmemd.cuda with mpirun -bind-to core -np 18 sander.MPI). Where "18" is the number of CPU cores that I use in my simulations (you have to replace it with the number of your CPU cores), note that simulations with the sander module on CPUs are much slower than those with the pmemd.cuda on GPUs.
pmemd.cuda -O -i sander.in1 -o sander.out1 -r mdrest1 -p prmtop -c prmcrd -ref prmcrd pmemd.cuda -O -i sander.in2 -o sander.out2 -r mdrest2 -p prmtop -c mdrest1 -ref prmcrd pmemd.cuda -O -i sander.in3 -o sander.out3 -r mdrest3 -p prmtop -c mdrest2 -ref prmcrd pmemd.cuda -O -i sander.in4 -o sander.out4 -r mdrest4 -p prmtop -c mdrest3 pmemd.cuda -O -i sander.in5 -o sander.out5 -r mdrest5 -p prmtop -c mdrest4 -x mdcrd5
The sander.in1 file Restrained minimization &cntrl irest=0,ntx=1, imin=1,maxcyc=1000,drms=0.0001,ntmin=2, ntc=2,ntf=2, cut=8.0, ntpr=100,ntwx=0,ntwv=0,ntwe=0, ipol=0,igb=0, ntr=1,restraint_wt=100, restraintmask="!:WA= & !@H=" nmropt=1 &end
&wt type='END' &end DISANG=rst
The sander.in2 file Restrained NVT equilibration &cntrl irest=0,ntx=1,ig=-1, nstlim=10000,dt=0.002, temp0=300.0,ntt=3,gamma_ln=2.0, ntc=2,ntf=2, cut=8.0, ntpr=500,ntwx=0,ntwv=0,ntwe=0, ntb=1, ipol=0,igb=0, ntr=1,restraint_wt=50, restraintmask="!:WA= & !@H=" nmropt=1 &end
&wt type='END' &end DISANG=rst
The sander.in3 file Restrained NPT equilibration &cntrl irest=1,ntx=5, nstlim=10000,dt=0.002, temp0=300.0,ntt=3,gamma_ln=2.0, ntc=2,ntf=2, cut=8.0, ntpr=500,ntwx=0,ntwv=0,ntwe=0, ntb=2,ntp=1,pres0=1.0,taup=1.0, ipol=0,igb=0, ntr=1,restraint_wt=50, restraintmask="!:WA= & !@H=" nmropt=1 &end
&wt type='END' &end DISANG=rst
The sander.in4 file Equilibration, 1 ns NPT &cntrl irest=1,ntx=5, nstlim=500000,dt=0.002, temp0=300.0,ntt=3,gamma_ln=2.0, ntc=2,ntf=2, cut=8.0, ntpr=500,ntwx=0,ntwv=0,ntwe=0, ntb=2,ntp=1,pres0=1.0,taup=1.0, ipol=0,igb=0, ntr=0 nmropt=1 &end
&wt type='END' &end DISANG=rst
The sander.in5 file Production 100 ns, sampling frequency 10 ps &cntrl irest=1,ntx=5, nstlim=50000000,dt=0.002, temp0=300.0,ntt=3,gamma_ln=2.0, ntc=2,ntf=2, cut=8.0, ntpr=5000,ntwx=5000,ntwv=0,ntwe=0, ntb=2,ntp=1,pres0=1.0,taup=1.0, ipol=0,igb=0, ntr=0 nmropt=1 &end
&wt type='END' &end DISANG=rst
The lines that are shown in boldface in the sander.in files are only necessary if you have a metal center in your protein (cf. section 2). |
Mehdi Irani Teaching duties Methods |