Locating waters in HSP90 binding cavity

In this tutorial, you will estimate the position of the waters in the binding cavity of N-terminal domain of HSP90 with a small fragment bound1. To achieve this we will use the JAWS2 technique.

To come back to the tutorial index, click here.

Prerequisites

Both the fragment and protein have been prepared for use with ProtoMS. For instance, they are correctly protonated. If using your own structures, you must make sure they are prepared correctly prior to simulation. You can read more about setting up structures here.


Simple setup

Setup
The simplest way to set up the simulation is by typing the following

python $PROTOMSHOME/protoms.py -s jaws1 -p protein.pdb -l fragment.pdb 

this sets up one JAWS stage 1 simulation, to locate the waters which can be placed 2 A around your ligand.

The simulations will run 5 m equilibration steps and 40 m production steps. Output will be printed to files every 100 k moves.

the -s jaws1 argument tells the script to setup a JAWS stage 1 simulation. The -p and -l arguments specifies the protein and the ligand, respectively.

You can read more about the files that the setup script creates further down on this page.

You can visualise the systems that will be simulated with for instance VMD:

vmd -m protein.pdb fragment.pdb water_clr.pdb jaws1_wat.pdb jaws1_box.pdb

Execution
To run the simulation you need to execute

$PROTOMSHOME/protoms3 run_jaws.cmd

The calculations will take approximately 10 h to complete.

Analysis
When the simulations are finished you need to analyse the simulation and conclude on a location for the waters in the binding site. We should start by separating the global pdb of your simulation into individual snapshots. Go to the directory called out_jaws and run:
python $PROTOMSHOME/tools/divide_pdb.py
and then, to visualize the positions of the waters at one given point in the simulation, for example, at the middle of your production run, type:
vmd snapshot_200.pdb
and try to see if there was any water in the binding site (the residue name of your jaws waters is "WA1"). Remember that whichever JAWS waters you may find, are the "on" — fully interacting — waters (with θ > 0.95, where θ corresponds to the scaling of their interactions). But looking at each snapshot individually is hard to get a whole picture of where your waters were distributed throughout the simulation. To get a more visual picture of this we produce and visualize a water density file. To generate such a file run:
python $PROTOMSHOME/tools/calc_density.py -f all.pdb -r wa1
and the file grid.dx will be produced. Such a file can be visualized in vmd. Try looking at it with:
vmd -m grid.dx ../protein.pdb ../fragment.pdb

You should be looking at something similar to this (after fiddling a bit with the graphical representation options in VMD):

A different visualization option are the clusters. These can be generated going to your out folder and using the following command

python $PROTOMSHOME/tools/calc_clusters.py -i all.pdb -m WA1 -a O00
and then visualized with
vmd -m clusters.pdb ../protein.pdb ../fragment.pdb

You should be looking at something similar to the image below, after tweaking the graphical representation options in VMD. It is probably interesting to display only the most populated clusters by showing only resid 1 to 5 , for example, since clusters are assigned residue ID in order of most populated clusters.

Exploring more options

Running JAWS on different areas of your simulation
By running protoms.py as above, protoms will try to locate the waters found in a box with edges at least 2A away from any of your ligand atoms. However, you may wish to find the location of waters in any area of your simulation (though normally it will include your binding site).

You can define any desired area (in box shape) of your simulation to run JAWS on, by specifying the box origin and dimentions with the flag --gcmcbox, such as it is done below

python $PROTOMSHOME/protoms.py -l fragment.pdb -p protein.pdb -s jaws1 -c big --gcmcbox 59 31 23 8 8 8

where the box will be a cube of 8A, with the "lower" corner in x = 59, y = 31, z = 23, and the "higher" corner in x = 67, y = 39, z = 31 (all values in A). The -c flag specifies the base-name of the generated .cmd files; it has been applied in this example to avoid overwriting the previously generated run_jaws.cmd file.

Running JAWS specifying the number of waters
When protoms is run such in the initial example (first section in the tutorial), or in the example above, the default is chosen for the amount of water molecules to include inside the box. This can be modified with the flag --gcmcwater. The amount of waters included as JAWS probes can have an effect in the results of your simulation. Feel free to play around with this parameter. It can be specified as shown below

python $PROTOMSHOME/protoms.py -l fragment.pdb -p protein.pdb -s jaws1 -c big --gcmcbox 59 31 23 8 8 8 --gcmcwater 10

Running longer simulations
There are two arguments that you can invoke to run a longer simulation by typing for instance

python $PROTOMSHOME/protoms.py -s jaws1 -p protein.pdb -l fragment.pdb --nequil 10E6 --nprod 50E6

you will run 10 m equiliration steps and 50 m production steps (instead of the 5 m and 40 m that is default)

Running independent repeats
Usually it is wise to run several independent repeats of your calculation to check for convergence. The argument that controls the number of independent repeats is called --repeats or just -r.

by typing for instance

python $PROTOMSHOME/protoms.py -s jaws1 -p protein.pdb -l fragment.pdb -r 5

you will create 5 input files for the JAWS calculation. Therefore, you also need to execute ProtoMS 5 times with the different input files. The output will be in 5 different folders, e.g. out1_jaws and out2_jaws.


Files created by the setup script

Files describing the fragment
You can read more about the setup of ligands here.
Files describing the protein
You can read more about the setup of proteins here.

Simulation specific files
You can read more about the input files in the ProtoMS manual. However, some sections of it is worth mentioning here:
jaws1 0

The jaws1 line indicates the initial value of θ for the JAWS waters.

originx 59.774
originy 31.007
originz 20.931
x 13.02
y 8.188
z 8.402

These lines include the dimensions of the JAWS box, to which the JAWS waters are constrained. The originx, originy, originz, x, y, z parameters indicate the origin and length of the JAWS box in each spatial coordinate.


Setting up the simulation with individual tools

In this advanced section we will go through the setup of the simulations step-by-step using individual setup scripts rather than protoms.py. This is relevant in cases where you are interested in specifying certain parameters which are pre-determined by protoms.py.
Setting up the ligand
We will start setting up the ligand from the initial fragment.pdb.

  • First we need to make sure that the very first line of the fragment.pdb contains a directive, telling ProtoMS the name of the solute. The line should read HEADER XDK, where xdk is the residue name in the fragment, and can be added by typing

    sed -i "1iHEADER XDK" fragment.pdb

  • Thereafter, we will create the force field for the fragment using AmberTools. The force field will be GAFF with AM1-BCC charge. Type

    python $PROTOMSHOME/tools/ambertools.py -f fragment.pdb -n XDK

    and this will execute the AmberTools programs antechamber and parmchk, creating the files fragment.prepi and fragment.frcmod, respectively.

  • These files are in Amber format and in order to use them in ProtoMS we need to reformat them into a ProtoMS template file. This file will also contain a z-matrix that describes how the ligand will be sampled during the simulation. To do this, you can type

    python $PROTOMSHOME/tools/build_template.py -p fragment.prepi -f fragment.frcmod -o fragment.tem -n XDK

    this will creates the files fragment.tem containing the ProtoMS template file and fragment.zmat. It is a good idea to check this files to see if the script has defined the molecule properly.

  • Now we have setup the fragment, but there is another "ligand" that we have not mentioned: the water molecules used for the JAWS calculation. These are considered by ProtoMS as "grand canonical solutes", and will be generated later, but it is important to mention here that they will, as well, need a template. However, there is no need to generate this template, since it contains some particularities and it can be found in $PROTOMSHOME/data/gcmc_tip4p.tem.

    Setting up the protein
  • First, we will change the name of atom and residue names in the protein PDB-file such that they agree with the ProtoMS naming convention. Type

    python $PROTOMSHOME/tools/convertatomnames.py -p protein.pdb -o protein_pms.pdb -s amber -c $PROTOMSHOME/data/atomnamesmap.dat

    The converted structure will be in protein_pms.pdb. This execution assumes that the Amber naming convention is used in protein.pdb.

  • We have crystallographic water in the protein structure and we need to convert them to the water model we will be using in the simulation (TIP4P). This can be done by

    python $PROTOMSHOME/tools/convertwater.py -p protein_pms.pdb -o protein_pms_t4p.pdb

    creating protein_pms_t4p.pdb.

  • Next step is to truncate the protein, creating a scoop. This will enable us to solvate the protein is a droplet and thereby reducing the number of simulated particles. The command is

    python $PROTOMSHOME/tools/scoop.py -p protein_pms_t4p.pdb -l fragment.pdb -o protein_scoop.pdb

    The protein scoop is centred on the fragment molecule and all residue further than 20 A are cut-away. The scoop is written to protein_scoop.pdb

    Setting up the solvation and jaws-waters
  • First we will create a droplet of TIP4P water molecules around both ligand and protein. To this, type:

    python $PROTOMSHOME/tools/solvate.py -b $PROTOMSHOME/data/wbox_tip4p.pdb -s fragment.pdb -pr protein_scoop.pdb -o water.pdb -g droplet

    this will create a droplet with 30 A radius centred on the fragment molecule. The droplet is written to water.pdb

  • The solvate.py script adds the crystallographic waters from the scoop to the droplet. Therefore, we need to remove them from the scoop PDB-file, to avoid duplicates.

    sed -i -e "/T4P/d" -e "/TER/d" protein_scoop.pdb

  • Next we need to create the set of waters that will form our JAWS waters. This can be done by typing:
    python $PROTOMSHOME/tools/distribute_waters.py -b 59 31 23 8 8 8 -m 10 -o jaws1_wat.pdb 

    with this command, the file jaws1_water.pdb will be generated. It will contain 10 water molecules (flag -m). The box dimensions to which the waters will be constrained are specified with the flag -b, and this information on the box limits is stored in the first line of the jaws1_water.pdb file.

  • The last step in this section consists of clearing the water molecules from the solvation droplet which overlap with the JAWS-box. Do this by typing:
    python $PROTOMSHOME/tools/clear_gcmcbox.py -b jaws1_wat.pdb -s water.pdb -o water_clr.pdb

    Feel free to check the difference between water.pdb and water_clr.pdb, there may be a few less waters in the later, but there might be no difference at all. Now we have all the files need to run the simulation. As you noticed, there is some difference in the files produced with this step-by-step procedure and those created with protoms.py.

    Making ProtoMS input files
    To make the input files for ProtoMS type

    python $PROTOMSHOME/tools/generate_input.py -s jaws1 -p protein_scoop.pdb -l fragment.pdb -t fragment.tem -pw water_clr.pdb --gcmcwater jaws1_wat.pdb -o run

    creating run_jaws.cmd.

    References

    1. Murray et.al. J. Med. Chem., 2010, 53(16), pp 5942-5955

    2. Michel et.al. J. Phys. Chem. B, 2009, 113(40), pp 13337-13346