Difference between revisions of "2014 DOCK tutorial with HIV Protease"

From Rizzo_Lab
Jump to: navigation, search
(IV. Generating Box and Grid)
(VI. Virtual Screening)
Line 357: Line 357:
 
   num_scored_conformers                                        100
 
   num_scored_conformers                                        100
 
   rank_ligands                                                no
 
   rank_ligands                                                no
 +
 +
'''Analyze the Results'''
  
 
==VII. Running DOCK in Parallel on Seawulf==
 
==VII. Running DOCK in Parallel on Seawulf==

Revision as of 15:57, 3 March 2014

For additional Rizzo Lab tutorials see DOCK Tutorials. Use this link Wiki Formatting as a reference for editing the wiki. This tutorial was developed collaboratively by the AMS 536 class of 2013, using DOCK v6.6.

I. Introduction

DOCK

DOCK is a molecular docking program used in drug discovery. It was developed by Irwin D. Kuntz, Jr. and colleagues at UCSF (see UCSF DOCK). This program, given a protein binding site and a small molecule, tries to predict the correct binding mode of the small molecule in the binding site, and the associated binding energy. Small molecules with highly favorable binding energies could be new drug leads. This makes DOCK a valuable drug discovery tool. DOCK is typically used to screen massive libraries of millions of compounds against a protein to isolate potential drug leads. These leads are then further studied, and could eventually result in a new, marketable drug. DOCK works well as a screening procedure for generating leads, but is not currently as useful for optimization of those leads.

DOCK 6 uses an incremental construction algorithm called anchor and grow. It is described by a three-step process:

  1. Rigid portion of ligand (anchor) is docked by geometric methods.
  2. Non-rigid segments added in layers; energy minimized.
  3. The resulting configurations are 'pruned' and energy re-minimized, yielding the docked configurations.

HIV Protease

HIV protease (HIVPR) is a protein involved in viral maturation during the life cycle of HIV. HIVPR is an approximately 22 kDa homodimer with 99 residues per chain. Inhibition of this protein has been shown to be an effective form of treatment of HIV. Currently-available HIVPR inhibitors generally take the form of a symmetric cyclic urea compound. For more information, see Lam et al.

Organizing Directories

While performing docking, it is convenient to adopt a standard directory structure / naming scheme, so that files are easy to find / identify. For this tutorial, we will use something similar to the following:

~username/AMS536/dock-tutorial/00.files/
                              /01.dockprep/
                              /02.surface-spheres/
                              /03.box-grid/
                              /04.dock/
                              /05.mini-virtual-screen/
                              /06.virtual-screen/
                             

In addition, most of the important files that are derived from the original crystal structure will be given a prefix that is the same as the PDB code, '1HVR'. The following sections in this tutorial will adhere to this directory structure / naming scheme.

II. Preparing the Receptor and Ligand

Downloading the PDB Structure (1HVR)

Go to PDB (Protein Data Bank) website (http://www.rcsb.org/pdb/home/home.do) enter the protein ID (1HVR), search for the PDB file and download it as a text form.

Preparing the ligand and receptor in Chimera

Put the 1HVR PDB file in 00.file/folder. If you are in the 00.files/directory, then tap the command:

 cp ~/Downloads/1HVR.pdb ./

When you are preparing you PDB files, you have to make some modifications on your original file. For example: we changed the atom name form "CSO" to "CYS" and deleted to lines "OD" and "HD". When you finish the modifications, save it as "1HVR.modified.pdb" in the 00.files/. And then we will create 4 files in 01.dockprep/ directory:

 1HVR.dockprep.mol2  
 1HVR.ligand.mol2  
 1HVR.receptor.mol2  
 1HVR.receptor.noH.pdb

Create the dockprep file

For the "1HVR.dockprep.mol2" file: open the 1HVR.modified.pdb in Chimera; delete the water molecules; delete the original hydrogen atoms; add the charge and add the hydrogen atoms manually. Or you can do all of the above by clicking Tools -> Structure Editing -> Dock Prep. Note when adding the charge to the ligand, you can choose AMBER ff99SB as the charge model and chose gasteiger as the charge method. In this 1HVR case, we set Net Charge to 0. (You may have to consider the chemistry of the ligand when assigning a charge state). Finally save the file as a mol2 format in Chimera.

Create the ligand file

To create the ligand file: Open 1HVR.dockprep.mol2 in Chimera and select the "1HVR" chain and delete it.

Create the receptor file

Select the residue LIG and delete it and save it as a mol2 file.

Create the receptor without hydrogen atoms

Open the "1HVR.receptor.mol2" file and delete all of the hydrogen atoms in Chimera and save it as a pdb file.

III. Generating Receptor Surface and Spheres

Generating the Receptor Surface

Check to make sure 02.surface-spheres directory exists under dock-tutorial. If not then make the following directory:

mkdir 02.surface-sphgen
cd 02.surface-sphgen

The following steps will be carried out to generate the receptor surface using Chimera:

Open Chimera by simply typing chimera into the terminal window

| Go File -> Open and choose the PDB file of the protein containing no hydrogens (1HVR.receptor.noH.pdb) from 01.dockprep

| Further, Actions -> Surface -> Show

| Go Tools -> Structure Editing -> Write DMS in order to obtain a dms file, which we will need to place spheres

| In the new window save the surface as 1HVR.receptor.dms

1HVR Receptor surface

Placing Spheres

We will be using SPHGEN to generate spheres: see the DOCK online owners manual for additional information:

<http://dock.compbio.ucsf.edu/DOCK_6/dock6_manual.htm>

The following steps will be used to place the spheres on the receptor surface:

1. Create a file called INSPH and fill it out as follows, then save it. This input file tells SPHGEN what to do, details of each line are below:

1HVR.receptor.dms
R
X
0.0
4.0
1.4
1HVR.receptor.sph

Input File Details:

 1HVR.receptor.dms - surface file from the previous step
 R - tells SPHGEN to place spheres either outside of the surface (R) or inside the surface (L)
 X - tells SPHGEN the subset of surface points to be used (X=all points)
 0.0 - prevents generation of large spheres with close surface contacts(defalut=0.0)
 4.0 - maximum sphere radius in angstroms (default=4.0)
 1.4 - minimum sphere radius in angstroms (default=radius of probe)
 1HVR.receptor.sph - clustered spheres file that we want to generate

2. Run the sphgen program from the terminal:

sphgen -i INSPH -o OUTSPH
-i tells sphgen where the input file INSPH is
INSPH tells sphgen what to do
-o tells sphgen what to call the oputput file
OUTSPH is the output file containing the sphere information

3. (optional) To look at the spheres generated, you need to put them into PDB format.

Run showsphere, by typing the follwoing into the terminal:

showsphere

You will be prompted with the following questions:

Enter name of sphere cluster file:
     1HVR.receptor.sph
Enter cluster number to process (<0 = all):
     -1
Generate surfaces as well as pdb files (<N>/Y)?
     N
Enter name for output file prefix:
     output_spheres
Process cluster 0 (contains ALL spheres) (<N>/Y)? 
     N

You can then open the receptor file in Chimera as well as the output_spheres.pdb file and should see many spheres placed all over the receptor surface.

1HVR Receptor surface with spheres

4. Run the sphere_selector program from in the terminal to select spheres of interest (note: not all of the spheres in the previous image are in the active site so we want to eliminate them):

sphere_selector 1HVR.receptor.sph ../01.dockprep/1HVR.ligand.mol2 8.0

This program goes to select the spheres within a user-defined radius (8.0 here) of a target molecule from a previously obtained file: 1HVR.receptor.sph. In turn, a new file selected_spheres.sph will be generated.

5. Run showsphere to visualize the spheres:

showsphere

When prompted on the command line, answering the questions as follows:

Enter name of sphere cluster file:
     selected_spheres.sph
Enter cluster number to process (<0 = all):
     -1
Generate surfaces as well as pdb files (<N>/Y)?
     N
Enter name for output file prefix:
     output_spheres_selected
Process cluster 0 (contains ALL spheres) (<N>/Y)? 
     N
  • Launch Chimera, choose File -> Open, choose 1HVR.receptor.noH.pdb
  • Go File -> Open, choosing output_spheres_selected.pdb
  • Go Select -> Residue -> SPH
  • Finally, Actions -> Atoms/Bonds -> sphere

The final image should look similar to the example below:

1HVR Receptor surface with spheres within 8A

IV. Generating Box and Grid

Mosavverul Arkin

1.) Make a new directory and name it: 03.box-grid/

      mkdir 03.box-grid


2.) Make a new file in this directory and name it showbox.in

     vim showbox.in

3.) This will automatically open the file showbox.in. Edit the file showbox.in as follows:

   Y                                               #Yes, generate a box
   8.0                                             #Size of the box in Angstroms
   ../02.surfaces-spheres/selected_spheres.sph     #Sphere.sph file
   1                                               #Cluster number
   1HVR.box.pdb                                    #Name of the output file


4.) Save the file by the inputting the command:

    :wq

5.) Run the command:

    showbox < showbox.in
Receptor surface with generated box

6.) Make a new file in the same directory (03.box-grid/) and name it: grid.in

   vim grid.in


This will automatically open the file grid.in. Edit this file with the following parameters:

    compute_grids                  yes                   #compute scoring grids (yes)                 
    grid_spacing                   0.4                   #what is the distance between grid points along each axis (in Angstroms).                 
    output_molecule                no                    #write up coordinates of the receptor into a new file                 
    contact_score                  no                    #compute contact grid? default is no
    energy_score                   yes                   #compute energy score? yes - we are using this method to compute force fields on probes
    energy_cutoff_distance         9999                  #the max distance between atoms for the energy contribution to be computed
    atom_model                     a                     #atom_model u means united atom model where atoms are attached to hydrogens, and a stands for all-atom model, where hydrogens on carbons are treated separately
    attractive_exponent            6                     #attractive component stands for exponent of the attractive LJ term in VDW potential
    repulsive_exponent             9                     #repulsive component stands for exponent in the repulsive LJ term in VDW potential
    distance_dielectric            yes
    dielectric_factor              4
    bump_filter                    yes
    bump_overlap                   0.75
    receptor_file                  ../01.dockprep/1HVR.receptor.mol2
    box_file                       1HVR.box.pdb
    vdw_definition_file            ../00.files/vdw_AMBER_parm99.defn
    score_grid_prefix              1HVR.grid






distance dielectric stands for the dielectric constant to be linearly dependent on distance distance dielectric factor is the coefficient of the dielectric bump filter flag determines if we want to screen orientation for clashes before scoring and minimization bump_overlap stands for the fraction of allowed overlap where 1 corresponds to no allowed overlap and 0 corresponds to full overlap being permitted. our receptor file the box file we generated in the Box section VDW parameters file Prefix for the grid file name. All the extensions will be generated automatically.

V. Docking a Single Molecule for Pose Reproduction

Jess Junjie Kai Rigid Ligand DOCKing

Rigid Ligand Docking is usually used to check the validity of docking program and user-defined parameters. It also help users to identify any mistakes in

Rigid docking result;RMSD = 0.9848; green: experiment; purple: predicted

Flexible Ligand DOCKing

In the case of flexible ligand DOCKing, the ligand is allowed to be flexible. This type of docking allows the ligand to structurally rearrange in response to the receptor. Initially the largest rigid substructure of the ligand is identified; then the flexible layers are identified. The orientations are ranked according to their score, and are grouped by root mean squared deviation (RMSD). To run the docking calculation, use the program dock6 that is distributed with DOCK in the /bin directory. You need to generate an input file by answering questions interactively or manually via text file. The program will generate an output file summarizing the parameters used in the run, any warning or error messages, and summary information about the best scoring pose. It will also produce a structure (.mol2) file, containing the geometric coordinates of the best pose as well as a summary of the interaction energy of that pose.

Best scored ligand result
Worst scored ligand result

VI. Virtual Screening

Virtual Screening Introduction

Virtual screening is a method used to predict most favorable ligand binding to a target protein within a ligand database. It also allows for comparison of both qualitative (e.g. position in binding site) and quantitative (e.g. grid scores, internal energy) data pertaining to the each screened ligand with the originally docked molecule.

To perform virtual screening, we use HIVPR.ligands.005.mol2, a mol2 file which contains 5 small molecules to be the virtual library. This is a reasonable computational cost for a quick search, so we can conduct it on own computer . After which, we may able to conduct virtual screening within a larger database HIVPR.ligands.100.mol2. Since the computational cost of virtual screening is much higher, it is better to run it on Seawulf.

Running Virtual Screen

Before runing virtual screen, a nre dictionary 05.mini-virtual-screen is created and copy ligand files HIVPR.ligands.005.mol2 and HIVPR.ligands.100.mol2 from wjallen. Since most parameters of dock.in and new virtual screen inputare the same, we can do few modifications based on dock.in but rename it as mini-virtual-screen.in.

  mkdir 05.mini-virtual-screen
  cd 05.mini-virtual-screen
  cp ~wjallen/AMS536/multi-mol2/files/HIVPR.ligands.005.mol2 ./
  cp ../04.dock/dock.in/ ./
  mv dock.in mini-virtual-screen.in

Here is the content of mini-virtual-screen.in

  ligand_atom_file                                             HIVPR.ligands.005.mol2
  limit_max_ligands                                            no
  skip_molecule                                                no
  read_mol_solvation                                           no
  calculate_rmsd                                               no
  use_database_filter                                          no
  orient_ligand                                                yes
  automated_matching                                           yes
  receptor_site_file                                           ../02.surface-spheres/selected_spheres.sph
  max_orientations                                             1000
  critical_points                                              no
  chemical_matching                                            no
  use_ligand_spheres                                           no
  use_internal_energy                                          yes
  internal_energy_rep_exp                                      12
  flexible_ligand                                              yes
  user_specified_anchor                                        no
  limit_max_anchors                                            no
  min_anchor_size                                              5
  pruning_use_clustering                                       yes
  pruning_max_orients                                          1000
  pruning_clustering_cutoff                                    100
  pruning_conformer_score_cutoff                               100
  use_clash_overlap                                            no
  write_growth_tree                                            no
  bump_filter                                                  no
  score_molecules                                              yes
  contact_score_primary                                        no
  contact_score_secondary                                      no
  grid_score_primary                                           yes
  grid_score_secondary                                         no
  grid_score_rep_rad_scale                                     1
  grid_score_vdw_scale                                         1
  grid_score_es_scale                                          1
  grid_score_grid_prefix                                       ../03.box-grid/1HVR.grid
  multigrid_score_secondary                                    no
  dock3.5_score_secondary                                      no
  continuous_score_secondary                                   no
  descriptor_score_secondary                                   no
  gbsa_zou_score_secondary                                     no
  gbsa_hawkins_score_secondary                                 no
  SASA_descriptor_score_secondary                              no
  amber_score_secondary                                        no
  minimize_ligand                                              yes
  minimize_anchor                                              yes
  minimize_flexible_growth                                     yes
  use_advanced_simplex_parameters                              no
  simplex_max_cycles                                           1
  simplex_score_converge                                       0.1
  simplex_cycle_converge                                       1.0
  simplex_trans_step                                           1.0
  simplex_rot_step                                             0.1
  simplex_tors_step                                            10.0
  simplex_anchor_max_iterations                                500
  simplex_grow_max_iterations                                  1000
  simplex_grow_tors_premin_iterations                          0
  simplex_random_seed                                          0
  simplex_restraint_min                                        no
  atom_model                                                   all
  vdw_defn_file                                                ../00.files/vdw_AMBER_parm99.defn
  flex_defn_file                                               ../00.files/flex.defn
  flex_drive_file                                              ../00.files/flex_drive.tbl
  ligand_outfile_prefix                                        1HVR.output
  write_orientations                                           no
  num_scored_conformers                                        100
  rank_ligands                                                 no

Analyze the Results

VII. Running DOCK in Parallel on Seawulf

Fengfei Lu

VIII. Frequently Encountered Problems

Tianao

Mike

The most common problem experienced in UNIX is having the wrong file path. ALWAYS CHECK TO MAKE SURE THAT YOUR FILE PATH IS CORRECT This can be facilitated with the Tab button which will auto complete what you have started to type.

Ivan

Make sure you are opening the correct files for each visualization step (e.g. opening the receptor file with the ligand still there may make analyzing your docking results a bit more challenging).

Use the auto-fill (tab) function but double check to make sure the correct file was auto-filled (e.g. if you type 1H then hit tab you might get 1HVR_selected_spheres but the file you want actually is 1HVR_selected_spheres.sph).

Junjie

Kai

For users with no UNIX experience and limited programming experience - it's a good idea to run through a couple tutorials before you start using UNIX. Even if you type in the commands without fully understanding their implications, become familiar with the basic commands and shortcuts. It will help a lot when you actually start using UNIX!

Arkin

Yan

Yao

Lu

Fengfei

Mosavverul

Joe

Write some text here..

 command or input file