MM/GBSA

This text describes how to setup MM/GBSA simulations from a set of explicit simulations. It depends heavily on automated scripts that rather easily can be extended to new/special situations. For instance, there is extentions of these scripts to handle explicit water, but they are not discussed in this text.

Read here about program requisites. For the moment the calculations works best with Amber10, and it is known that Amber11 does not work at all because they changed the main script. This means that it only works on Platon at the moment.

Here is notes on setting up calculations in implicit solvent.

Samuel Genheden, 2011-2012

Setting up the MM/GBSA preparation folder

Start by setting up a "preparation" folder. This is a common respository folder for a set of ligands.

  • Create a folder called Mmgbsa-prep, at the same level as the "ligand-folders" in your directory tree.

  • Goto that folder, and copy all the contents from /away/bio/Setup_mmgbsa/, i.e. cp -p /away/bio/Setup_mmgbsa/* .

    Your folder now contains a number of different files that will be used to setup and extract information from the MM/GBSA calculations. The important files for your project are:

  • leapcom* - is a number of leapcom-files to setup prmtop-files for the complex, receptor, and ligand, for either mm_pbsa.in1 or mm_pbsa.in2.
  • ligands.dat - is a database file that contains information for leap for the various ligands in your project
  • mm_pbsa.in1 - input file for the enthalpy calculations
  • mm_pbsa.in2 - input file for the entropy calculations
  • prep_mmpbsa.pl - preparation script, used to setup the calculations.
  • leaprc.ff99SB - a modification of the Amber file that prevents leap to put terminals on truncated protein.

    Optionally, if you would like to use a different GB-model than OCBI:

  • Open all leapcom* files in an editor. In each of these files there will be a line set default pbradii mbondi2 that you need to change to set default pbradii bondi in order to get the correct radii for the GB-calculations.
  • Open mm_pbsa.in1. Find a line that reads IGB 2, and change it to IGB 7 Also find a line that reads MS 1 and change it to MS 0 to use the same SASA model as in the simulations. Note: this cannot be done with sed.

  • Open prep_mmgbsa.pl. At the top of the script there is a line that sets the $PREPPATH variable, viz. my $PREPPATH = "/temp1/paulius/Ferritin/Ex-wat/Prep-mmgbsa/";. Change the path to the absolute path of your preparation path.

  • Open center.bash and make sure that it will read the filenames of your mdcrd-files. Usually, one only needs to edit the line starting with trajin

    Editing the ligand database

  • Open ligands.dat.

    In this file there should be one line for each ligand.
    Each line consist of the ligand name (the residue name, with captial letters) followed by a comma (,), and a number of leap commands that loads force fields information into leap.

    Make sure that you have the correct relative or absolut path.

    The leapcom* files in the preparation folder contains a %%LOADSTR%% string that will be replaced by the information in the ligands.dat file. Insert %%BR%% between different commands to add a return in the leapcom-file.


    Setting up the MM/GBSA calculation

  • Goto a "ligand-folder" and into the folder where the simulation was done

  • Make sure that you have at least copied the *mdcrd files from Platon. It is best if you copy all of the files from Platon back to your local computer

  • Make a folder called Mmgbsa, and goto that folder.

  • Execute $PREPFOLDER/center.bash L02 ../prmtop
    and change $PREPFOLDER to the path of your Mmgbsa_prep-folder and L02 to the residue name of your ligand

    This will create a new set of mdcrd-files with the ligand centered in the box

  • Execute $PREPFOLDER/prep_mmpbsa.pl com L02 ../prmtop ../prmcrd 40 40 0
    and change $PREPFOLDER to the path of your Mmgbsa_prep-folder and L02 to the residue name of your ligand. 40 40 indicates the number of records in each independent trajectory, and the number of independent trajectories, respectively. The last flag (0) indicates that all of the independent trajectories have the same prmtop-file.

    The script has now created a sub-folder for each of the simulations and a number of other files.

    Type $PREPFOLDER/prep_mmpbsa.pl -h to get a help for the script. It can do a lot more than just this standard thing.

    If the script returns an error, it is in most cases because you have not edit the ligands.dat file properly, or that center.bash were unable to create mdcrd-files.

  • Execute bash get_prm1.bash and then bash get_prm2.bash and look out for error messages. Make sure that all the prmtop-files were created by executing ls -l *prmtop*, there should be five of them and none of them should have zero-size.
    The most common cause of error is that ligands.dat were not edited properly.

  • Make a Mmgbsa folder on Platon and copy the contents of the local folder. Make sure to add -r to recursively copy all the sub-folders as well.

  • Goto the Mmgbsa-folder and execute the following lines (copy all of them and paste the in the terminal-window).
    for X in {1..40}
    do
    cd R${X}_mmpbsa
    qsub q_r${X}_mmpbsa
    cd ..
    done
    and it will add 40 jobs to the queue. You could add this to a file as well.
    Note: this has to be done on Platon.

  • Check the results by executing
    grep -E "^40" R*/*all.out | grep -c 40
    it should return 240.


    Program requisite

    As described in Ulfs original notes, the nmode program needs to be edited. The modified nmode.f file is located in /away/bio/Setup_mmgbsa/Code/.
    In this folder, there is also a the mm_pbsa_loc.pl script that is also needed for the truncated entropy calculations, as well as mm_pbsa_calceneent.pm that contains some modification for radii of "unusual" atoms.

    The program also requires the maketrcrds program to create truncated coordinates. This program is located locally in /away/bio/Amber/MakeCrds/ and in /sw/pkg/bio/Amber/MakeCrds/ on Platon.

    MM/GBSA - from implicit simulations

    This text describes how to setup MM/GBSA simulations from a set of implicit simulations. It was written as instructions for a master student, and as such it not well-developed.

    Setting up the MM/GBSA preparation folder

  • Create a folder called Mmgbsa-prep, at the same level as the "ligand-folders" in your directory tree.

  • Goto that folder, and copy all the contents from /away/bio/Setup_mmgbsa/, i.e. cp -p /away/bio/Setup_mmgbsa/* .

    Your folder now contains a number of different files that will be used to setup and extract information from the MM/GBSA calculations. The important files for your project are:

  • leapcom* - is a number of leapcom-files to setup prmtop-files for the complex, receptor, and ligand, for either mm_pbsa.in1 or mm_pbsa.in2.
  • ligands.dat - is a database file that contains information for leap for the various ligands in your project
  • mm_pbsa.in1 - input file for the enthalpy calculations
  • mm_pbsa.in2 - input file for the entropy calculations
  • prep_mmpbsa_gb.pl - preparation script, used to setup the calculations.

  • Open all leapcom* files in an editor. In each of these files there will be a line set default pbradii mbondi2 that you need to change to set default pbradii bondi in order to get the correct radii for the GB-calculations. This can be done by executing: sed -i "s/set default pbradii mbondi2/set default pbradii bondi/" leapcom*

  • Open mm_pbsa.in1. Find a line that reads IGB 2, and change it to IGB 7 to use the same GB model as in the simulations. Also find a line that reads MS 1 and change it to MS 0 to use the same SASA model as in the simulations. Note: this cannot be done with sed.

  • Open prep_mmgbsa_gb.pl. At the top of the script there is a line that sets the $PREPPATH variable, viz. my $PREPPATH = "/temp1/paulius/Ferritin/Ex-wat/Prep-mmgbsa/";. Change the path to the absolute path of your preparation path.

    Editing the ligand database

  • Open ligands.dat.

    In this file there should be one line for each ligand.
    Each line consist of the ligand name (the residue name, with captial letters) followed by a comma (,), and a number of leap commands that loads force fields information into leap.

    For your project it will be enough to have a loadAmberPrep line for each of the ligands.
    Make sure that you have the correct relative or absolut path.

    Setting up the MM/GBSA calculation

  • Goto a "ligand-folder" and into the Sim_gb

  • Goto the Sim_gb-folder on Platon and execute gunzip *mdcrd5.gz to unzip all the mdcrd-files.

  • Make a folder called Mmgb, and goto that folder.

  • Execute $PREPFOLDER/prep_mmpbsa_gb.pl com L02 ../prmtop ../prmcrd 40 40 0
    and change $PREPFOLDER to the path of your Mmgbsa_prep-folder and L02 to the residue name of your ligand.
    Note: This has to be done locally in your mounted directory.

    The script has now created a sub-folder for each of the simulations and a number of other files.

  • Execute bash get_prm1.bash and then bash get_prm2.bash and look out for error messages. Make sure that all the prmtop-files were created by executing ls -l *prmtop*, there should be five of them and none of them should have zero-size.
    Note: This has to be done locally in your mounted directory.

  • Goto the Mmgb-folder and execute the following lines (copy all of them and paste the in the terminal-window).
    for X in {1..40}
    do
    cd R${X}_mmpbsa
    qsub q_r${X}_mmpbsa
    cd ..
    done
    and it will add 40 jobs to the queue. You could add this to a file as well.
    Note: this has to be done on Platon.

  • Check the results by executing
    grep -E "^40" R*/*all.out | grep -c 40
    it should return 240.