|
Difference between revisions of "Common applications of the MMTSB toolset"
Line 2: | Line 2: | ||
==Preparing protein simulations (from PDBfile to CHARMM trajectory):== | ==Preparing protein simulations (from PDBfile to CHARMM trajectory):== | ||
− | + | This tutorial illustrates the step-by-step instructions to perform explicit solvent molecular dynamics of a protein solvated in a cubic box of water at 0.15 mM ionic NaCl concentration using MMTSB toolset. The set of files and scripts provided below can be easily adapted to perform dynamics of any other protein of interest with different solvent box shapes and sizes. Here, we use RNAase peptide as an example. | |
+ | |||
+ | Step-by-step instructions follows: | ||
+ | |||
+ | 1. Download RNAse A crystal structure, 1RNU.pdb protein database (http://www.rcsb.org/pdb/home/home.do). | ||
+ | |||
+ | 2. Extract first 13 residues using convpdb.pl w/ charmm22 output, centered and with segnames: | ||
+ | |||
+ | convpdb.pl -center -sel 1:13 -segnames -out charmm22 -nohetero 1RNU.pdb > 1rnu_cpep.pdb | ||
+ | |||
+ | 3. Visualize peptide fragment with VMD to make sure its what you want. | ||
+ | |||
+ | 4. Minimize | ||
+ | |||
+ | minCHARMM.pl -par minsteps=0,sdsteps=100,sdstepsz=0.02 \ | ||
+ | -par trunc=switch,dielec=rdie,cutnb=12,cuton=8,cutoff=11 \ | ||
+ | -par param=22x,cmap\ | ||
+ | -par xtop=top_all27_prot_na.rtf \ | ||
+ | -par xpar=par_all27_prot_na.prm \ | ||
+ | -par blocked,nter=none,cter=ct3 \ | ||
+ | -cons heavy self 1:13_5 -log mmtsbmin.log \ | ||
+ | 1rnu_cpep.pdb > 1rnu_cpep_mmtsbmin.pdb | ||
+ | |||
+ | |||
+ | 5. Write protein structure file (PSF) for use in later steps involving: (1) viewing structures in VMD and (2) determining how | ||
+ | many ions to add for making the charge on the system neutral. | ||
+ | |||
+ | |||
+ | genPSF.pl -par param=22x,cmap\ | ||
+ | -par xtop=top_all27_prot_na.rtf \ | ||
+ | -par xpar=par_all27_prot_na.prm \ | ||
+ | -par blocked,nter=none,cter=ct3 \ | ||
+ | 1rnu_cpep_mmtsbmin.pdb > 1rnu_cpep_mmtsbmin.psf | ||
+ | |||
+ | In the following steps 6-8 of the tutorial, we will use the mmtsb tools to solvate our | ||
+ | minimized peptide. First we'll add solvent and ions to the system. Then we'll minimize the complete | ||
+ | system with harmonic restraints on the peptide. The system will then be ready for molecular dynamics runs. | ||
+ | |||
+ | 6. Add solvent box | ||
+ | |||
+ | First use convpdb.pl to solvate and add segment names - extract boxsize! | ||
+ | |||
+ | convpdb.pl -out charmm22 -solvate -cubic -cutoff 8 1rnu_cpep_mmtsbmin.pdb | \ | ||
+ | convpdb.pl -out charmm22 -segnames > 1rnu_cpep_solv.pdb | ||
+ | |||
+ | Find that boxsize is a: | ||
+ | |||
+ | In this case it is 38.743553 Å cube. Note down this number. | ||
+ | |||
+ | 7. Add ions at 0.15 mM concentration to the cubic box of water generated above | ||
+ | |||
+ | Counterions are added to solvated system by specifying the number of positive (SOD) | ||
+ | and/or negative ions (CLA). Determine how many and which type of ions we need to | ||
+ | neutralize this system using following CHARMM script. | ||
+ | |||
+ | $CHARMMEXEC pdbfile=1rnu_cpep_mmtsbmin psffile=1rnu_cpep_mmtsbmin boxsize=38.74355 <how_many_ions_to_add.inp> how_many_ions_to_add.log & | ||
+ | |||
+ | |||
+ | Now grep information about the number of ions from the log file "how_many_ions_to_add.log" | ||
+ | |||
+ | Number of positive ions (SOD) to add: | ||
+ | |||
+ | grep "NPOS ->" how_many_ions_to_add.log | tail -1 | ||
+ | |||
+ | Number of negative ions (CLA) to add: | ||
+ | |||
+ | grep "NNEG ->" how_many_ions_to_add.log | tail -1 | ||
+ | |||
+ | Now finally add ions to the system with waterbox: | ||
+ | |||
+ | convpdb.pl -ions SOD:5=CLA:6 1rnu_cpep_solv.pdb > 1rnu_cpep_solvions.pdb | ||
+ | |||
+ | |||
+ | Generate complete PSF of the system including waterbox and ions. | ||
+ | |||
+ | genPSF.pl -par param=22x,cmap\ | ||
+ | -par xtop=top_all27_prot_na.rtf \ | ||
+ | -par xpar=par_all27_prot_na.prm \ | ||
+ | -par blocked,nter=none,cter=ct3 \ | ||
+ | 1rnu_cpep_solvmin.pdb >1rnu_cpep_solvmin.psf | ||
+ | |||
+ | |||
+ | 8. Now minimize with restraints and SHAKE constraints using PME electrostatics | ||
+ | |||
+ | minCHARMM.pl -par minsteps=0,sdsteps=500,sdstepsz=0.02 \ | ||
+ | -par trunc=switch,cutnb=12,cuton=8,cutoff=11 \ | ||
+ | -par param=22x,cmap\ | ||
+ | -par xtop=top_all27_prot_na.rtf \ | ||
+ | -par xpar=par_all27_prot_na.prm \ | ||
+ | -par blocked,nter=none,cter=ct3 \ | ||
+ | -par shake,boxsize=38.743553,nblisttype=bycb \ | ||
+ | -cons heavy self 1:13_5 \ | ||
+ | -cmd mmtsbminsolvate.inp -log mmtsbminsolvate.log \ | ||
+ | 1rnu_cpep_solvions.pdb > 1rnu_cpep_solvions_min.pdb | ||
+ | |||
+ | 9. Running dynamics on this peptide in water using the mmtsb tool mdCHARMM.pl. | ||
+ | |||
+ | In this section we will take the configuration we just generated and run 2 ps of | ||
+ | molecular dynamics on this system. | ||
+ | |||
+ | mdCHARMM.pl -par dynsteps=1000,dynens=NPT,dynitemp=298,dyneqfrq=1000 \ | ||
+ | -par dynnose=1,dynoutfrq=10,dynpress=1,echeck=20000 \ | ||
+ | -par trunc=switch,cutnb=12,cuton=8,cutoff=11 \ | ||
+ | -par param=22x,cmap\ | ||
+ | -par xtop=top_all27_prot_na.rtf \ | ||
+ | -par xpar=par_all27_prot_na.prm \ | ||
+ | -par blocked,nter=none,cter=ct3 \ | ||
+ | -par shake,boxsize=38.743553,nblisttype=bycb \ | ||
+ | -cmd mmtsbdynsolvate.inp -log mmtsbdynsolvate.log \ | ||
+ | -enerout 1rnu_cpep_d0.ene -trajout 1rnu_cpep_d0.dcd \ | ||
+ | -restout 1rnu_cpep_d0.res -final 1rnu_cpep_d0.pdb \ | ||
+ | 1rnu_cpep_solvions_min.pdb | ||
+ | |||
+ | 10. Restart dynamics and run perform production dynamics using NVT ensemble: | ||
+ | |||
+ | mdCHARMM.pl -par dynsteps=1000,dyntemp=298,\ | ||
+ | -restart 1rnu_cpep_d0.res -final - 1rnu_cpep_solvions_min.pdb \ | ||
+ | -par trunc=switch,cutnb=12,cuton=8,cutoff=11 \ | ||
+ | -par param=22x,cmap\ | ||
+ | -par xtop=top_all27_prot_na.rtf \ | ||
+ | -par xpar=par_all27_prot_na.prm \ | ||
+ | -par blocked,nter=none,cter=ct3 \ | ||
+ | -par shake,boxsize=38.743553,nblisttype=bycb \ | ||
+ | -cmd mmtsbdynsolvate1.inp -log mmtsbdynsolvate1.log \ | ||
+ | -enerout 1rnu_cpep_d1.ene -trajout 1rnu_cpep_d1.dcd \ | ||
+ | -restout 1rnu_cpep_d1.res | ||
+ | |||
+ | 11. Visualization of dynamics trajectory. | ||
+ | |||
+ | After running the dynamics, visualize the trajectory with vmd. To look at molecular dynamics trajectories with vmd, you need to read in a pdb and a dcd file. The following command will | ||
+ | accomplish this. | ||
+ | |||
+ | vmd 1rnu_cpep_solvions_min.pdb 1rnu_cpep_d1.dcd | ||
+ | |||
+ | |||
+ | 12. Analyze trajectory | ||
+ | |||
+ | |||
+ | Process dcd trajectory to extract root mean square deviation of the peptide CA atoms with respect to the starting conformation as a function of time. | ||
+ | |||
+ | processDCD.pl -psf 1rnu_cpep_solvmin.psf 1rnu_cpep_d1.dcd -rms CA 1rnu_cpep_solvmin.pdb >rms.out | ||
+ | |||
+ | |||
+ | gnuplot | ||
+ | pl ‘rms.out’ u 2:3 w lp | ||
==Preparing protein:DNA simulations (from PDBfile to CHARMM trajectory):== | ==Preparing protein:DNA simulations (from PDBfile to CHARMM trajectory):== |
Revision as of 00:53, 11 July 2009
Page currently under construction!!!
Contents
Preparing protein simulations (from PDBfile to CHARMM trajectory):
This tutorial illustrates the step-by-step instructions to perform explicit solvent molecular dynamics of a protein solvated in a cubic box of water at 0.15 mM ionic NaCl concentration using MMTSB toolset. The set of files and scripts provided below can be easily adapted to perform dynamics of any other protein of interest with different solvent box shapes and sizes. Here, we use RNAase peptide as an example.
Step-by-step instructions follows:
1. Download RNAse A crystal structure, 1RNU.pdb protein database (http://www.rcsb.org/pdb/home/home.do).
2. Extract first 13 residues using convpdb.pl w/ charmm22 output, centered and with segnames:
convpdb.pl -center -sel 1:13 -segnames -out charmm22 -nohetero 1RNU.pdb > 1rnu_cpep.pdb
3. Visualize peptide fragment with VMD to make sure its what you want.
4. Minimize
minCHARMM.pl -par minsteps=0,sdsteps=100,sdstepsz=0.02 \ -par trunc=switch,dielec=rdie,cutnb=12,cuton=8,cutoff=11 \ -par param=22x,cmap\ -par xtop=top_all27_prot_na.rtf \ -par xpar=par_all27_prot_na.prm \ -par blocked,nter=none,cter=ct3 \ -cons heavy self 1:13_5 -log mmtsbmin.log \ 1rnu_cpep.pdb > 1rnu_cpep_mmtsbmin.pdb
5. Write protein structure file (PSF) for use in later steps involving: (1) viewing structures in VMD and (2) determining how
many ions to add for making the charge on the system neutral.
genPSF.pl -par param=22x,cmap\
-par xtop=top_all27_prot_na.rtf \ -par xpar=par_all27_prot_na.prm \ -par blocked,nter=none,cter=ct3 \
1rnu_cpep_mmtsbmin.pdb > 1rnu_cpep_mmtsbmin.psf
In the following steps 6-8 of the tutorial, we will use the mmtsb tools to solvate our minimized peptide. First we'll add solvent and ions to the system. Then we'll minimize the complete system with harmonic restraints on the peptide. The system will then be ready for molecular dynamics runs.
6. Add solvent box
First use convpdb.pl to solvate and add segment names - extract boxsize!
convpdb.pl -out charmm22 -solvate -cubic -cutoff 8 1rnu_cpep_mmtsbmin.pdb | \ convpdb.pl -out charmm22 -segnames > 1rnu_cpep_solv.pdb
Find that boxsize is a:
In this case it is 38.743553 Å cube. Note down this number.
7. Add ions at 0.15 mM concentration to the cubic box of water generated above
Counterions are added to solvated system by specifying the number of positive (SOD) and/or negative ions (CLA). Determine how many and which type of ions we need to neutralize this system using following CHARMM script.
$CHARMMEXEC pdbfile=1rnu_cpep_mmtsbmin psffile=1rnu_cpep_mmtsbmin boxsize=38.74355 <how_many_ions_to_add.inp> how_many_ions_to_add.log &
Now grep information about the number of ions from the log file "how_many_ions_to_add.log"
Number of positive ions (SOD) to add:
grep "NPOS ->" how_many_ions_to_add.log | tail -1
Number of negative ions (CLA) to add:
grep "NNEG ->" how_many_ions_to_add.log | tail -1
Now finally add ions to the system with waterbox:
convpdb.pl -ions SOD:5=CLA:6 1rnu_cpep_solv.pdb > 1rnu_cpep_solvions.pdb
Generate complete PSF of the system including waterbox and ions.
genPSF.pl -par param=22x,cmap\
-par xtop=top_all27_prot_na.rtf \ -par xpar=par_all27_prot_na.prm \ -par blocked,nter=none,cter=ct3 \
1rnu_cpep_solvmin.pdb >1rnu_cpep_solvmin.psf
8. Now minimize with restraints and SHAKE constraints using PME electrostatics
minCHARMM.pl -par minsteps=0,sdsteps=500,sdstepsz=0.02 \ -par trunc=switch,cutnb=12,cuton=8,cutoff=11 \ -par param=22x,cmap\ -par xtop=top_all27_prot_na.rtf \ -par xpar=par_all27_prot_na.prm \ -par blocked,nter=none,cter=ct3 \ -par shake,boxsize=38.743553,nblisttype=bycb \ -cons heavy self 1:13_5 \ -cmd mmtsbminsolvate.inp -log mmtsbminsolvate.log \ 1rnu_cpep_solvions.pdb > 1rnu_cpep_solvions_min.pdb
9. Running dynamics on this peptide in water using the mmtsb tool mdCHARMM.pl.
In this section we will take the configuration we just generated and run 2 ps of molecular dynamics on this system.
mdCHARMM.pl -par dynsteps=1000,dynens=NPT,dynitemp=298,dyneqfrq=1000 \
-par dynnose=1,dynoutfrq=10,dynpress=1,echeck=20000 \ -par trunc=switch,cutnb=12,cuton=8,cutoff=11 \ -par param=22x,cmap\ -par xtop=top_all27_prot_na.rtf \ -par xpar=par_all27_prot_na.prm \ -par blocked,nter=none,cter=ct3 \ -par shake,boxsize=38.743553,nblisttype=bycb \ -cmd mmtsbdynsolvate.inp -log mmtsbdynsolvate.log \ -enerout 1rnu_cpep_d0.ene -trajout 1rnu_cpep_d0.dcd \ -restout 1rnu_cpep_d0.res -final 1rnu_cpep_d0.pdb \ 1rnu_cpep_solvions_min.pdb
10. Restart dynamics and run perform production dynamics using NVT ensemble:
mdCHARMM.pl -par dynsteps=1000,dyntemp=298,\ -restart 1rnu_cpep_d0.res -final - 1rnu_cpep_solvions_min.pdb \ -par trunc=switch,cutnb=12,cuton=8,cutoff=11 \ -par param=22x,cmap\ -par xtop=top_all27_prot_na.rtf \ -par xpar=par_all27_prot_na.prm \ -par blocked,nter=none,cter=ct3 \ -par shake,boxsize=38.743553,nblisttype=bycb \ -cmd mmtsbdynsolvate1.inp -log mmtsbdynsolvate1.log \ -enerout 1rnu_cpep_d1.ene -trajout 1rnu_cpep_d1.dcd \ -restout 1rnu_cpep_d1.res
11. Visualization of dynamics trajectory.
After running the dynamics, visualize the trajectory with vmd. To look at molecular dynamics trajectories with vmd, you need to read in a pdb and a dcd file. The following command will accomplish this.
vmd 1rnu_cpep_solvions_min.pdb 1rnu_cpep_d1.dcd
12. Analyze trajectory
Process dcd trajectory to extract root mean square deviation of the peptide CA atoms with respect to the starting conformation as a function of time.
processDCD.pl -psf 1rnu_cpep_solvmin.psf 1rnu_cpep_d1.dcd -rms CA 1rnu_cpep_solvmin.pdb >rms.out
gnuplot
pl ‘rms.out’ u 2:3 w lp
Preparing protein:DNA simulations (from PDBfile to CHARMM trajectory):
replace text
Preparing protein simulations for replica-exchange:
replace text
Visualizing the electrostatic surface potential of a macromolecule:
First, we select the protein from a PDB file (1enh.pdb), add hydrogen atoms, and center the molecule at the origin. Then, we calculate the electrostatic potential on a grid and generate the molecular surface of the protein. Finally, we view the molecule in VMD, projecting the electrostatic potential onto the molecular surface.
convpdb.pl -nsel protein 1ENH.pdb | complete.pl | convpdb.pl -center > 1enh.center.pdb
pbCHARMM.pl -emap phi.dx 1enh.center.pdb
pbCHARMM.pl -dx -grid epsx grid.dx 1enh.center.pdb
vmd 1enh.center.pdb phi.dx grid.dx
In vmd, select:
Graphics/Representations/Drawing Method [Surf]
Graphics/Representations/Coloring Method [Volume]
For a more thorough description, see the MMTSB Tool Set - Continuum electrostatics calculations tutorial.
For more examples, download and follow the tutorials prepared for past MMTSB workshops.