Difference between revisions of "2017 DOCK tutorial 2 with PDB 3GPL NEW"

From Rizzo_Lab
Jump to: navigation, search
(→‎III. Generating Receptor Surface and Spheres)
(→‎3PGL)
 
(19 intermediate revisions by the same user not shown)
Line 15: Line 15:
  
 
    
 
    
−
In this tutorial we will use PDB code 3PGL, the deposited crystal structure of Scp1 in complex with rabeprazole.
+
In this tutorial we will use PDB code 3PGL, the deposited crystal structure of human small C-terminal domain phosphatase 1 in complex with rabeprazole.
  
−
[[Image:2017_3pgl_receptor_surface.png|thumb|center|450px| Scp1 in complex with rabeprazole. Waters and Magnesium ion in binding pocket are shown]]
+
[[Image:2017_3pgl_receptor_surface.png|thumb|center|450px| 3pgl in complex with rabeprazole. Waters and Magnesium ion in binding pocket are shown]]
  
 
===Organizing Directories===
 
===Organizing Directories===
Line 39: Line 39:
  
 
==II. Preparing the Receptor and Ligand==
 
==II. Preparing the Receptor and Ligand==
 +
 +
One issue that you will encounter when running Dock using the standard procedure on this system is that the ligand will move greatly. This is because normally all waters in the crystal structure are deleted, but they seem to be important for the positioning of the ligand. That causes the procedure in this tutorial to be slightly different from previous tutorials.
  
 
Download the PDB file (3PGL) and move file 3PGL.pdb to 00.files
 
Download the PDB file (3PGL) and move file 3PGL.pdb to 00.files
Line 47: Line 49:
 
raw_3pgl.pdb was opened with VI terminal editor
 
raw_3pgl.pdb was opened with VI terminal editor
  
−
    The header information, connect records, and waters were deleted. Change HETATM to ATOM, and delete all rows that have B in the fifth column
+
Delete the header information, connect records, and waters except for waters 2, 21, 32, 39, and 258 (these numbers are listed in the 6th column of the dpb file). Change HETATM to ATOM, and delete all rows that have B in the fifth column
 +
Change RZX A to LIG B to make the ligand a different chain
  
−
    Change RZX A to LIG B to make the ligand a different chain
+
Now open the system in chimera and use
 +
Tools-->Structure Editing-->Add Charge
 +
Select AMBER ff99SB for Standard Residues, and add a charge of +2 to Mg and 0 to everything else
  
 
'''Prepare ligand and receptor files for dock'''
 
'''Prepare ligand and receptor files for dock'''
  
−
''Create dockprep file''
+
All waters were deleted from the file except for the ones believed to interact with the ligand and the Mg atom. These are believed to be important for the proper binding of the ligand. To properly orient the waters, run the system through tleap, keeping everything fixed except the 5 water molecules (this process can be seen in the amber tutorials on this wiki). This new system will yield better docking results, and should be saved as 3pgl.dockprep.mol2.
−
 
 
−
Open raw_3pgl.pdb file in chimera, in Tools, find structure editing, then click AddH to add hydrogen in this system.
 
−
 
 
−
Also in structure editing, click Add charge, then, changing AMBER ff14SB to AMBER ff99SB and changing net charge to +2.
 
−
 
 
−
Save it as 3pgl.dockprep.mol2 in 01.dockprep.
 
−
 
 
−
''Create receptor file''
 
−
 
 
−
Open 3pgl.dockprep.mol2 in chimera, then click Select, Chain, B to choose ligand. Click Actions, Atom/Bonds, Delete and save the molecule as 3pgl.rec.mol2 in 01.dockprep.
 
  
−
''Create ligand file''
+
Open 3pgl.dockprep.mol2 in Chimera. Select the ligand and click Actions-->Atoms/Bonds-->delete. Save this file as 3pgl.rec.mol2.
  
−
Open 3pgl.dockprep.mol2 in chimera, using the same method delete the protein molecule. Save this as 3pgl.lig.mol2 in 01.dockprep.
+
With just the receptor open in Chimera, go to Select-->Chemistry-->element-->H and click Actions-->Atoms/Bonds-->delete. Save this file as 3pgl.rec.noH.mol2
  
−
Normally, a file containing the receptor with no hydrogens would be made, but we need to keep the hydrogens for this system because we are keeping waters as part of the receptor.
+
Open 3pgl.dockprep.mol2 in Chimera and select the receptor. Click Actions-->Atoms/Bonds-->delete to get rid of everything except the ligand. Save this as 3pgl.lig.mol2.
  
 
==III. Generating Receptor Surface and Spheres==  
 
==III. Generating Receptor Surface and Spheres==  
Line 121: Line 116:
 
The image that appears should resemble this:
 
The image that appears should resemble this:
  
−
[[Image:3pgl_all_spheres.png|thumb|center|1000px|4qmz receptor with all spheres]]
+
[[Image:3pgl_all_spheres.png|thumb|center|1000px|3PGL receptor with all spheres.]]
  
  
Line 138: Line 133:
 
The selected spheres with the receptor surface should look similar to that as seen below:
 
The selected spheres with the receptor surface should look similar to that as seen below:
  
−
[[Image:Selectedspheres.png|thumb|center|1000px|4qmz receptor surface with the selected spheres within 10.0 Angstroms]]
+
[[Image:3pgl_selected_spheres.png|thumb|center|1000px|3PGL receptor surface with spheres within 8.0 Angstroms of known ligand]]
  
 
==IV. Generating Box and Grid==
 
==IV. Generating Box and Grid==
Line 168: Line 163:
 
Compute the energy grid
 
Compute the energy grid
  
−
create grid.in file using
+
create grid.in and grid out files using
 
     touch grid.in
 
     touch grid.in
 +
    touch grid.out
 +
   
 +
 +
 +
Compute the grid, and save the output in the file grid.out
 +
    grid -i grid.in > grid.out
  
 
When prpmped, answer the questions as follows:
 
When prpmped, answer the questions as follows:
Line 267: Line 268:
 
This structure can be visualized using chimera and loading the surface into chimera along with the minimized and unminimized ligand structures:
 
This structure can be visualized using chimera and loading the surface into chimera along with the minimized and unminimized ligand structures:
  
−
[[Image:image_min.png|thumb|center|1000px|4qmz with ligand before and after minimization]]
+
[[Image:01_minimized_lig.png|thumb|center|1000px|3pgl with ligand before and after minimization]]
  
 
===Footprint Calculation===
 
===Footprint Calculation===
Line 345: Line 346:
 
The output will look something like this:
 
The output will look something like this:
  
−
[[Image:figure_fps_minimized.png|thumb|center|1000px|Footprint comparison between unminimized and minimized ligand structure]]
+
[[Image:Fps_lig_minlig.png|thumb|center|1000px|Footprint comparison between unminimized and minimized ligand structure]]
  
 
===Rigid Docking===
 
===Rigid Docking===
Line 425: Line 426:
  
 
The best scoring result should look like
 
The best scoring result should look like
−
  Picture
+
  [[Image:02_rigid_bestscore.png|thumb|center|800px|Rigid docking results]]
  
 
===Flexible Docking===
 
===Flexible Docking===
Line 515: Line 516:
  
 
These can be viewed in the same way as before and the best scoring result should look like
 
These can be viewed in the same way as before and the best scoring result should look like
−
  Picture
+
  [[Image:03_flex_bestscore.png|thumb|center|800px|Flexible docking results]]
  
 
===Fixed Anchor Docking===
 
===Fixed Anchor Docking===
Line 603: Line 604:
 
  cluster_rmsd_threshold                                      2.0
 
  cluster_rmsd_threshold                                      2.0
 
  rank_ligands                                                no
 
  rank_ligands                                                no
 +
 +
The top scoring molecule is pictured below
 +
[[Image:04_fad_bestscore.png|thumb|center|800px|Fixed Anchor Docking results]]
  
 
==VI. Virtual Screening==
 
==VI. Virtual Screening==
 +
At this point we will begin docking potential inhibitors into the active site of 4qmz. The compound library we will use contains 25,000 molecules and was previously prepared for another docking experiment. The compound library contains molecules with rotatable bonds number between 7 and 15. The library used is titled small_ligand_library.mol2, this library can be found at /gpfs/projects/AMS536/virtual_screen_library on seawulf.
 +
 +
We began by copying all of our local directories outlined in the first section to a project directory on seawulf to allow use to expedite the docking process over the 25,000 compounds in our test library.
 +
 +
We will create a virtual-screen.in file that looks like the following:
 +
 +
  conformer_search_type                                        flex
 +
  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.0
 +
  pruning_conformer_score_scaling_factor                      1.0
 +
  use_clash_overlap                                            no
 +
  write_growth_tree                                            no
 +
  use_internal_energy                                          yes
 +
  internal_energy_rep_exp                                      12
 +
  internal_energy_cutoff                                      100.0
 +
  ligand_atom_file                                            small_ligand_library.mol2
 +
  limit_max_ligands                                            no
 +
  skip_molecule                                                no
 +
  read_mol_solvation                                          no
 +
  calculate_rmsd                                              yes
 +
  use_rmsd_reference_mol                                      yes
 +
  rmsd_reference_filename                                      ../01.dockprep/4qmz.lig.mol2
 +
  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
 +
  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/grid
 +
  multigrid_score_secondary                                    no
 +
  dock3.5_score_secondary                                      no
 +
  continuous_score_secondary                                  no
 +
  footprint_similarity_score_secondary                        no
 +
  pharmacophore_score_secondary                                no
 +
  descriptor_score_secondary                                  no
 +
  gbsa_zou_score_secondary                                    no
 +
  gbsa_hawkins_score_secondary                                no
 +
  SASA_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                                  500
 +
  simplex_grow_tors_premin_iterations                          0
 +
  simplex_random_seed                                          0
 +
  simplex_restraint_min                                        no
 +
  atom_model                                                  all
 +
  vdw_defn_file                                            /gpfs/projects/AMS536/zzz.programs/dock6/parameters/vdw_AMBER_parm99.defn
 +
  flex_defn_file                                              /gpfs/projects/AMS536/zzz.programs/dock6/parameters/flex.defn
 +
  flex_drive_file                                              /gpfs/projects/AMS536/zzz.programs/dock6/parameters/flex_drive.tbl
 +
  ligand_outfile_prefix                                        4qmz.virtualscreen
 +
  write_orientations                                          no
 +
  num_scored_conformers                                        1
 +
  write_conformations                                          no
 +
  cluster_conformations                                        yes
 +
  cluster_rmsd_threshold                                      2.0
 +
  rank_ligands                                                no
 +
 +
This script is submitted to the seawulf cluster using a script called screen.sh which contains the following:
 +
 +
  #! /bin/tcsh
 +
  #PBS -l nodes=4:ppn=28
 +
  #PBS -l walltime=24:00:00
 +
  #PBS -o jobstatus.out
 +
  #PBS -e joberror.err
 +
  #PBS -N large_vs
 +
  #PBS -V
 +
  #PBS -q long
 +
  cd $PBS_O_WORKDIR
 +
  mpirun -np 112 dock6.mpi -i virtual-screen.in -o virtualscreen1.out
 +
 +
 +
This screen.sh input specifies that we will use 112 processors to dock our virtual library. We have a wall clock for the job of 24 hours and we will output jobstatus.out, joberror.err and virtualscreen1.out. In addition to these output files there is an output for the work performed by each processor individually.
 +
 +
===Cartesian Minimization of the docked molecules===
 +
 +
Now we will do cartesian minimization of the docked molecules. Cartesian here is referring to the receptor being explicitly present during the minimization.
 +
 +
Write the minimization input script minimize-vs.in:
 +
 +
  conformer_search_type                                        rigid
 +
  use_internal_energy                                          yes
 +
  internal_energy_rep_exp                                      12
 +
  internal_energy_cutoff                                      100.0
 +
  ligand_atom_file                                            4qmz.virtualscreen_scored.mol2
 +
  limit_max_ligands                                            no
 +
  skip_molecule                                                no
 +
  read_mol_solvation                                          no
 +
  calculate_rmsd                                              no
 +
  use_database_filter                                          no
 +
  orient_ligand                                                no
 +
  bump_filter                                                  no
 +
  score_molecules                                              yes
 +
  contact_score_primary                                        no
 +
  contact_score_secondary                                      no
 +
  grid_score_primary                                          no
 +
  grid_score_secondary                                        no
 +
  multigrid_score_primary                                      no
 +
  multigrid_score_secondary                                    no
 +
  dock3.5_score_primary                                        no
 +
  dock3.5_score_secondary                                      no
 +
  continuous_score_primary                                    yes
 +
  continuous_score_secondary                                  no
 +
  cont_score_rec_filename                                      ../01.dockprep/4qmz.rec.mol2
 +
  cont_score_att_exp                                          6
 +
  cont_score_rep_exp                                          12
 +
  cont_score_rep_rad_scale                                    1
 +
  cont_score_use_dist_dep_dielectric                          yes
 +
  cont_score_dielectric                                        4.0
 +
  cont_score_vdw_scale                                        1
 +
  cont_score_es_scale                                          1
 +
  footprint_similarity_score_secondary                        no
 +
  pharmacophore_score_secondary                                no
 +
  descriptor_score_secondary                                  no
 +
  gbsa_zou_score_secondary                                    no
 +
  gbsa_hawkins_score_secondary                                no
 +
  SASA_score_secondary                                        no
 +
  amber_score_secondary                                        no
 +
  minimize_ligand                                              yes
 +
  simplex_max_iterations                                      1000
 +
  simplex_tors_premin_iterations                              0
 +
  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_random_seed                                          0
 +
  simplex_restraint_min                                        no
 +
  atom_model                                                  all
 +
  vdw_defn_file                                                /gpfs/projects/AMS536/zzz.programs/dock6/parameters/vdw_AMBER_parm99.defn
 +
  flex_defn_file                                              /gpfs/projects/AMS536/zzz.programs/dock6/parameters/flex.defn
 +
  flex_drive_file                                              /gpfs/projects/AMS536/zzz.programs/dock6/parameters/flex_drive.tbl
 +
  ligand_outfile_prefix                                        4qmz.virtualscreen1.minimized
 +
  write_orientations                                          no
 +
  num_scored_conformers                                        1
 +
  rank_ligands                                                no
 +
 +
Submit the minimization job using with the script "job-min.sh" following lines:
 +
 +
  #!/bin/tcsh
 +
  #PBS -l nodes=16:ppn=28
 +
  #PBS -l walltime=12:00:00
 +
  #PBS -o jobstatusmin.out
 +
  #PBS -e jobstatusmin.err
 +
  #PBS -N chunk0min
 +
  #PBS -V
 +
  #PBS -q medium
 +
  module load intel/mpi/64/2017/0.098
 +
  cd $PBS_O_WORKDIR
 +
  mpirun -np 448 /gpfs/projects/AMS536/zzz.programs/dock6/bin/dock6.mpi -i minimize-vs.in -o minimize-vs.out
 +
 +
 +
At this point it would be good to rank our docked ligands by score and extract the 100 best scoring (100 lowest scoring, most negative) ligands to investigate further.
  
 
==VIII. Frequently Encountered Problems==
 
==VIII. Frequently Encountered Problems==

Latest revision as of 16:59, 8 March 2017

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 a subsection of the AMS 536 class of 2017, using DOCK v6.8.

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.

3PGL

In this tutorial we will use PDB code 3PGL, the deposited crystal structure of human small C-terminal domain phosphatase 1 in complex with rabeprazole.

3pgl in complex with rabeprazole. Waters and Magnesium ion in binding pocket are shown

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-Spring2016/dock-tutorial/00.files/
                                         /01.dockprep/
                                         /02.surface-spheres/
                                         /03.box-grid/
                                         /04.dock/
                                         /05.large-virtual-screen/
                                         /06.virtual-screen/
                                         /07.footprint/
                                         /08.print_fps

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, '3pgl'. The following sections in this tutorial will adhere to this directory structure/naming scheme.

II. Preparing the Receptor and Ligand

One issue that you will encounter when running Dock using the standard procedure on this system is that the ligand will move greatly. This is because normally all waters in the crystal structure are deleted, but they seem to be important for the positioning of the ligand. That causes the procedure in this tutorial to be slightly different from previous tutorials.

Download the PDB file (3PGL) and move file 3PGL.pdb to 00.files

3PGL.pdb was copied to raw_3pgl.pdb with command 

cp 3PGL.pdb raw_3pgl.pdb

raw_3pgl.pdb was opened with VI terminal editor 

Delete the header information, connect records, and waters except for waters 2, 21, 32, 39, and 258 (these numbers are listed in the 6th column of the dpb file). Change HETATM to ATOM, and delete all rows that have B in the fifth column
Change RZX A to LIG B to make the ligand a different chain

Now open the system in chimera and use 

Tools-->Structure Editing-->Add Charge
Select AMBER ff99SB for Standard Residues, and add a charge of +2 to Mg and 0 to everything else

Prepare ligand and receptor files for dock

All waters were deleted from the file except for the ones believed to interact with the ligand and the Mg atom. These are believed to be important for the proper binding of the ligand. To properly orient the waters, run the system through tleap, keeping everything fixed except the 5 water molecules (this process can be seen in the amber tutorials on this wiki). This new system will yield better docking results, and should be saved as 3pgl.dockprep.mol2.

Open 3pgl.dockprep.mol2 in Chimera. Select the ligand and click Actions-->Atoms/Bonds-->delete. Save this file as 3pgl.rec.mol2.

With just the receptor open in Chimera, go to Select-->Chemistry-->element-->H and click Actions-->Atoms/Bonds-->delete. Save this file as 3pgl.rec.noH.mol2

Open 3pgl.dockprep.mol2 in Chimera and select the receptor. Click Actions-->Atoms/Bonds-->delete to get rid of everything except the ligand. Save this as 3pgl.lig.mol2.

III. Generating Receptor Surface and Spheres

Change directory to 02.surface-sphere Open 3pgl.rec.pdb in Chimera

To generate the molecular surface: 

  Action --> Surface --> Show


Save the .dms file 

  Tools --> Structure editing --> write DMS

.dms save to 3pgl.rec.dms



Create surface spheres

create input file INSPH 

    3pgl.rec.dms     
    R
    X
    0.0
    4.0
    1.4
    3pgl.rec.sph

line 1 designates input file line 2 designates the generated spheres will be outside the receptor surface line 3 designates that all points on the receptor will be used line 4 designates the maximum surface radius of the spheres line 5 designates the minimum surface radius of the spheres line 6 designates the output file name

run sphgen 

    sphgen -i INSPH -o OUTSPH

sphgen is the sphere generation program from dock -i desginates the input file: INSPH -o designates the output file

At this point it is beneficial to visualize the spheres that were created. This can be done with chimera: Open chimera from the terminal choose file --> open 3pgl.rec.mol2 choose file --> open 3pgl.rec.sph

The image that appears should resemble this:

3PGL receptor with all spheres.


Then we need to select the spheres pertinent to our docking experiment. Usually these speheres will be the closest N spheres to the native ligand molecule.

Run 

    sphere_selector 3pgl.rec.sph ../01.dockprep/3pgl.lig.mol2 8.0

This command will select all of the spheres within 8.0 angstroms of the ligand and output them to selected_spheres.sph

To visualize the spheres using Chimera as previously done: Launch Chimera, choose File -> Open, choose 3pgl.rec.pdb File -> Open, choose output_spheres_selected.pdb Select -> Residue -> SPH Actions -> Atoms/Bonds -> sphere The selected spheres with the receptor surface should look similar to that as seen below:

3PGL receptor surface with spheres within 8.0 Angstroms of known ligand

IV. Generating Box and Grid

enter directory 03.box-grid

create input showbox.in 

    Y
    8.0
    ../02.surface-spheres/selected_spheres.sph
    1
    3pgl.box.pdb

This input designates to dock that: we want to create a box the box length should be 8.0 Angstroms use the selected spheres in the file desginated output the box to the file specified


to use this input type: 

    showbox < showbox.in

This box can be visualized in chimera using a similar approach as visualizing the spheres


Compute the energy grid

create grid.in and grid out files using 

    touch grid.in
    touch grid.out


Compute the grid, and save the output in the file grid.out 

    grid -i grid.in > grid.out

When prpmped, answer the questions as follows: 

  compute_grids             yes
  grid_spacing              0.4
  output_molecule           no
  contact_score             no
  energy_score              yes
  energy_cutoff_distance    9999
  atom_model                a
  attractive_exponent       6
  repulsive_exponent        9
  distance_dielectric       yes
  dielectric_factor         4
  bump_filter               yes
  bump_overlap              0.75
  receptor_file             ../01.dockprep/3pgl.rec.mol2
  box_file                  ../03.box-grid/3pgl.box.pdb
  vdw_definition_file       /opt/AMS536/dock6/parameters/vdw_AMBER_parm99.defn
  score_grid_prefix         grid


this script should output verbosely to the terminal and produce 2 output files, grid.bmp and grid.nrg, both binary files.

V. Docking a Single Molecule for Pose Reproduction

Minimization

create min.in file and run dock with this as the input 

 touch min.in
 dock6 -i min.in

When prompted, answer the questions as follows: 

  conformer_search_type                                        rigid
  use_internal_energy                                          yes
  internal_energy_rep_exp                                      12
  internal_energy_cutoff                                       100.0
  ligand_atom_file                                             ../01.dockprep/3pgl.lig.mol2
  limit_max_ligands                                            no
  skip_molecule                                                no
  read_mol_solvation                                           no
  calculate_rmsd                                               yes
  use_rmsd_reference_mol                                       yes
  rmsd_reference_filename                                      ../01.dockprep/3pgl.lig.mol2
  use_database_filter                                          no
  orient_ligand                                                no
  bump_filter                                                  no
  score_molecules                                              yes
  contact_score_primary                                        no
  contact_score_secondary                                      no
  grid_score_primary                                           no
  grid_score_secondary                                         no
  multigrid_score_primary                                      no
  multigrid_score_secondary                                    no
  dock3.5_score_primary                                        no
  dock3.5_score_secondary                                      no
  continuous_score_primary                                     yes
  continuous_score_secondary                                   no
  cont_score_rec_filename                                      ../01.dockprep/3pgl.rec.mol2
  cont_score_att_exp                                           6
  cont_score_rep_exp                                           12
  cont_score_rep_rad_scale                                     1
  cont_score_use_dist_dep_dielectric                           yes
  cont_score_dielectric                                        4.0
  cont_score_vdw_scale                                         1
  cont_score_es_scale                                          1
  footprint_similarity_score_secondary                         no
  pharmacophore_score_secondary                                no
  descriptor_score_secondary                                   no
  gbsa_zou_score_secondary                                     no
  gbsa_hawkins_score_secondary                                 no
  SASA_score_secondary                                         no
  amber_score_secondary                                        no
  minimize_ligand                                              yes
  simplex_max_iterations                                       1000
  simplex_tors_premin_iterations                               0
  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_random_seed                                          0
  simplex_restraint_min                                        no
  atom_model                                                   all
  vdw_defn_file                                                /opt/AMS536/dock6/parameters/vdw_AMBER_parm99.defn
  flex_defn_file                                               /opt/AMS536/dock6/parameters/flex.defn
  flex_drive_file                                              /opt/AMS536/dock6/parameters/flex_drive.tbl
  ligand_outfile_prefix                                        3pgl.min
  write_orientations                                           no
  num_scored_conformers                                        1
  rank_ligands                                                 no


this command will output 3pgl.lig.min_scored.mol2

This structure can be visualized using chimera and loading the surface into chimera along with the minimized and unminimized ligand structures:

3pgl with ligand before and after minimization

Footprint Calculation

At this point we are going to calculate the van der wals and electrostatic footprints for the unminimized and minimized ligand structures in relation to the receptor active site.

create the input file: footprint.in 

  conformer_search_type                                        rigid
  use_internal_energy                                          no
  ligand_atom_file                                             ../3pgl.lig.min_scored.mol2
  limit_max_ligands                                            no
  skip_molecule                                                no
  read_mol_solvation                                           no
  calculate_rmsd                                               no
  use_database_filter                                          no
  orient_ligand                                                no
  bump_filter                                                  no
  score_molecules                                              yes
  contact_score_primary                                        no
  contact_score_secondary                                      no
  grid_score_primary                                           no
  grid_score_secondary                                         no
  multigrid_score_primary                                      no
  multigrid_score_secondary                                    no
  dock3.5_score_primary                                        no
  dock3.5_score_secondary                                      no
  continuous_score_primary                                     no
  continuous_score_secondary                                   no
  footprint_similarity_score_primary                           yes
  footprint_similarity_score_secondary                         no
  fps_use_footprint_reference_mol2                             yes
  fps_footprint_reference_mol2_filename                        ../01.dockprep/3pgl.lig.mol2
  fps_foot_compare_type                                        Euclidean
  fps_normalize_foot                                           no
  fps_foot_comp_all_residue                                    yes
  fps_receptor_filename                                        ../01.dockprep/3pgl.rec.mol2
  fps_vdw_att_exp                                              6
  fps_vdw_rep_exp                                              12
  fps_vdw_rep_rad_scale                                        1
  fps_use_distance_dependent_dielectric                        yes
  fps_dielectric                                               4.0
  fps_vdw_fp_scale                                             1
  fps_es_fp_scale                                              1
  fps_hb_fp_scale                                              0
  pharmacophore_score_secondary                                no
  descriptor_score_secondary                                   no
  gbsa_zou_score_secondary                                     no
  gbsa_hawkins_score_secondary                                 no
  SASA_score_secondary                                         no
  amber_score_secondary                                        no
  minimize_ligand                                              no
  atom_model                                                   all
  vdw_defn_file                                                /opt/AMS536/dock6/parameters/vdw_AMBER_parm99.defn
  flex_defn_file                                               /opt/AMS536/dock6/parameters/flex.defn
  flex_drive_file                                              /opt/AMS536/dock6/parameters/flex_drive.tbl
  ligand_outfile_prefix                                        fps.min.output
  write_footprints                                             yes
  write_hbonds                                                 yes
  write_orientations                                           no
  num_scored_conformers                                        1
  rank_ligands                                                 no


this footprint calculation can be run in the same way as the minimization: 

  dock6 -i footprint.in

This will output three files: fps.min.output_scored.mol2 

                             fps.min.output_scored_footprint_scored.txt
                             fps.min.output_hbond_scored.txt

With an in house script, plot_footprint_single_magnitude.py the two sets of footprints will be plotted for comparison

this script will be run using: 

  python plot_footprint_single_magnitude.py fps.min.output_scored_footprint_scored.txt 50

The output will look something like this:

Footprint comparison between unminimized and minimized ligand structure

Rigid Docking

To run rigid docking, create rigid.in file and run dock with it as an input, use commands 

touch rigid.in
dock6 -i rigid.in

And answer the prompts as follows 

  conformer_search_type                                        rigid
  use_internal_energy                                          yes
  internal_energy_rep_exp                                      12
  internal_energy_cutoff                                       100.0
  ligand_atom_file                                             ../01.dockprep/3pgl.lig.mol2
  limit_max_ligands                                            no
  skip_molecule                                                no
  read_mol_solvation                                           no
  calculate_rmsd                                               yes
  use_rmsd_reference_mol                                       yes
  rmsd_reference_filename                                      ../01.dockprep/3pgl.lig.mol2
  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
  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/grid
  multigrid_score_secondary                                    no
  dock3.5_score_secondary                                      no
  continuous_score_secondary                                   no
  footprint_similarity_score_secondary                         no
  pharmacophore_score_secondary                                no
  descriptor_score_secondary                                   no
  gbsa_zou_score_secondary                                     no
  gbsa_hawkins_score_secondary                                 no
  SASA_score_secondary                                         no
  amber_score_secondary                                        no
  minimize_ligand                                              yes
  simplex_max_iterations                                       1000
  simplex_tors_premin_iterations                               0
  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_random_seed                                          0
  simplex_restraint_min                                        yes
  simplex_coefficient_restraint                                10.0
  atom_model                                                   all
  vdw_defn_file                                                /opt/AMS536/dock6/parameters/vdw_AMBER_parm99.defn
  flex_defn_file                                               /opt/AMS536/dock6/parameters/flex.defn
  flex_drive_file                                              /opt/AMS536/dock6/parameters/flex_drive.tbl
  ligand_outfile_prefix                                        3pgl_rigid_restraint
  write_orientations                                           no
  num_scored_conformers                                        40
  write_conformations                                          no
  cluster_conformations                                        yes
  cluster_rmsd_threshold                                       2.0
  rank_ligands                                                 no

The output file is 3pgl.rigid_restraint_scored.mol2, and is saved in the directory from where you ran dock.

To vies the results, you can open the files of the receptor and ligand in CHIMERA from the 01.dockprep directory. The open your rigid docking results using viewdock: 

Tools --> Surface/Binding Analysis --> ViewDock

Select the output file and set the file type to Dock 4,5,or 6. This should allow you to view the docked conformers ranked by grid score

The best scoring result should look like

Rigid docking results

Flexible Docking

To Run flexible docking, create flex.in and run dock6 with command 

touch flex.in
dock6 -i flex.in

When prompted, answer the questions as follows 

conformer_search_type                                        flex
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.0
pruning_conformer_score_scaling_factor                       1.0
use_clash_overlap                                            no
write_growth_tree                                            no
use_internal_energy                                          yes
internal_energy_rep_exp                                      12
internal_energy_cutoff                                       100.0
ligand_atom_file                                             ../01.dockprep/3pgl.lig.mol2
limit_max_ligands                                            no
skip_molecule                                                no
read_mol_solvation                                           no
calculate_rmsd                                               yes
use_rmsd_reference_mol                                       yes
rmsd_reference_filename                                      ../3pgl.min_scored.mol2
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
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/grid
multigrid_score_secondary                                    no
dock3.5_score_secondary                                      no
continuous_score_secondary                                   no
footprint_similarity_score_secondary                         no
pharmacophore_score_secondary                                no
descriptor_score_secondary                                   no
gbsa_zou_score_secondary                                     no
gbsa_hawkins_score_secondary                                 no
SASA_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                                  500
simplex_grow_tors_premin_iterations                          0
simplex_random_seed                                          0
simplex_restraint_min                                        yes
simplex_coefficient_restraint                                10.0
atom_model                                                   all
vdw_defn_file                                                /opt/AMS536/dock6/parameters/vdw_AMBER_parm99.defn 
flex_defn_file                                               /opt/AMS536/dock6/parameters/flex.defn
flex_drive_file                                              /opt/AMS536/dock6/parameters/flex_drive.tbl
ligand_outfile_prefix                                        3pgl_flex
write_orientations                                           no
num_scored_conformers                                        40
write_conformations                                          no
cluster_conformations                                        yes
cluster_rmsd_threshold                                       2.0
rank_ligands                                                 no

The output file is 3pgl.flex_scored.mol2, and is saved in the directory from where you ran dock.

These can be viewed in the same way as before and the best scoring result should look like

Flexible docking results

Fixed Anchor Docking

Along the same line as above, use 

touch fad.in
dock6 -i fad.in

When prompted, answer the questions as follows 

Fix This

conformer_search_type                                        flex
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.0
pruning_conformer_score_scaling_factor                       1.0
use_clash_overlap                                            no
write_growth_tree                                            no
use_internal_energy                                          yes
internal_energy_rep_exp                                      12
internal_energy_cutoff                                       100.0
ligand_atom_file                                             ../01.dockprep/3pgl.lig.mol2
limit_max_ligands                                            no
skip_molecule                                                no
read_mol_solvation                                           no
calculate_rmsd                                               yes
use_rmsd_reference_mol                                       yes
rmsd_reference_filename                                      ../3pgl.min_scored.mol2
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
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/grid
multigrid_score_secondary                                    no
dock3.5_score_secondary                                      no
continuous_score_secondary                                   no
footprint_similarity_score_secondary                         no
pharmacophore_score_secondary                                no
descriptor_score_secondary                                   no
gbsa_zou_score_secondary                                     no
gbsa_hawkins_score_secondary                                 no
SASA_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                                  500
simplex_grow_tors_premin_iterations                          0
simplex_random_seed                                          0
simplex_restraint_min                                        yes
simplex_coefficient_restraint                                10.0
atom_model                                                   all
vdw_defn_file                                                /opt/AMS536/dock6/parameters/vdw_AMBER_parm99.defn 
flex_defn_file                                               /opt/AMS536/dock6/parameters/flex.defn
flex_drive_file                                              /opt/AMS536/dock6/parameters/flex_drive.tbl
ligand_outfile_prefix                                        3pgl_flex
write_orientations                                           no
num_scored_conformers                                        40
write_conformations                                          no
cluster_conformations                                        yes
cluster_rmsd_threshold                                       2.0
rank_ligands                                                 no

The top scoring molecule is pictured below

Fixed Anchor Docking results

VI. Virtual Screening

At this point we will begin docking potential inhibitors into the active site of 4qmz. The compound library we will use contains 25,000 molecules and was previously prepared for another docking experiment. The compound library contains molecules with rotatable bonds number between 7 and 15. The library used is titled small_ligand_library.mol2, this library can be found at /gpfs/projects/AMS536/virtual_screen_library on seawulf.

We began by copying all of our local directories outlined in the first section to a project directory on seawulf to allow use to expedite the docking process over the 25,000 compounds in our test library.

We will create a virtual-screen.in file that looks like the following: 

  conformer_search_type                                        flex
  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.0
  pruning_conformer_score_scaling_factor                       1.0
  use_clash_overlap                                            no
  write_growth_tree                                            no
  use_internal_energy                                          yes
  internal_energy_rep_exp                                      12
  internal_energy_cutoff                                       100.0
  ligand_atom_file                                             small_ligand_library.mol2
  limit_max_ligands                                            no
  skip_molecule                                                no
  read_mol_solvation                                           no
  calculate_rmsd                                               yes
  use_rmsd_reference_mol                                       yes
  rmsd_reference_filename                                      ../01.dockprep/4qmz.lig.mol2
  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
  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/grid
  multigrid_score_secondary                                    no
  dock3.5_score_secondary                                      no
  continuous_score_secondary                                   no
  footprint_similarity_score_secondary                         no
  pharmacophore_score_secondary                                no
  descriptor_score_secondary                                   no
  gbsa_zou_score_secondary                                     no
  gbsa_hawkins_score_secondary                                 no
  SASA_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                                  500
  simplex_grow_tors_premin_iterations                          0
  simplex_random_seed                                          0
  simplex_restraint_min                                        no
  atom_model                                                   all
  vdw_defn_file                                            /gpfs/projects/AMS536/zzz.programs/dock6/parameters/vdw_AMBER_parm99.defn
  flex_defn_file                                               /gpfs/projects/AMS536/zzz.programs/dock6/parameters/flex.defn
  flex_drive_file                                              /gpfs/projects/AMS536/zzz.programs/dock6/parameters/flex_drive.tbl
  ligand_outfile_prefix                                        4qmz.virtualscreen
  write_orientations                                           no
  num_scored_conformers                                        1
  write_conformations                                          no
  cluster_conformations                                        yes
  cluster_rmsd_threshold                                       2.0
  rank_ligands                                                 no

This script is submitted to the seawulf cluster using a script called screen.sh which contains the following: 

  #! /bin/tcsh
  #PBS -l nodes=4:ppn=28
  #PBS -l walltime=24:00:00
  #PBS -o jobstatus.out
  #PBS -e joberror.err
  #PBS -N large_vs
  #PBS -V
  #PBS -q long
  cd $PBS_O_WORKDIR
  mpirun -np 112 dock6.mpi -i virtual-screen.in -o virtualscreen1.out


This screen.sh input specifies that we will use 112 processors to dock our virtual library. We have a wall clock for the job of 24 hours and we will output jobstatus.out, joberror.err and virtualscreen1.out. In addition to these output files there is an output for the work performed by each processor individually.

Cartesian Minimization of the docked molecules

Now we will do cartesian minimization of the docked molecules. Cartesian here is referring to the receptor being explicitly present during the minimization.

Write the minimization input script minimize-vs.in: 

  conformer_search_type                                        rigid
  use_internal_energy                                          yes
  internal_energy_rep_exp                                      12
  internal_energy_cutoff                                       100.0
  ligand_atom_file                                             4qmz.virtualscreen_scored.mol2
  limit_max_ligands                                            no
  skip_molecule                                                no
  read_mol_solvation                                           no
  calculate_rmsd                                               no
  use_database_filter                                          no
  orient_ligand                                                no
  bump_filter                                                  no
  score_molecules                                              yes
  contact_score_primary                                        no
  contact_score_secondary                                      no
  grid_score_primary                                           no
  grid_score_secondary                                         no
  multigrid_score_primary                                      no
  multigrid_score_secondary                                    no
  dock3.5_score_primary                                        no
  dock3.5_score_secondary                                      no
  continuous_score_primary                                     yes
  continuous_score_secondary                                   no
  cont_score_rec_filename                                      ../01.dockprep/4qmz.rec.mol2
  cont_score_att_exp                                           6
  cont_score_rep_exp                                           12
  cont_score_rep_rad_scale                                     1
  cont_score_use_dist_dep_dielectric                           yes
  cont_score_dielectric                                        4.0
  cont_score_vdw_scale                                         1
  cont_score_es_scale                                          1
  footprint_similarity_score_secondary                         no
  pharmacophore_score_secondary                                no
  descriptor_score_secondary                                   no
  gbsa_zou_score_secondary                                     no
  gbsa_hawkins_score_secondary                                 no
  SASA_score_secondary                                         no
  amber_score_secondary                                        no
  minimize_ligand                                              yes
  simplex_max_iterations                                       1000
  simplex_tors_premin_iterations                               0
  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_random_seed                                          0
  simplex_restraint_min                                        no
  atom_model                                                   all
  vdw_defn_file                                                /gpfs/projects/AMS536/zzz.programs/dock6/parameters/vdw_AMBER_parm99.defn
  flex_defn_file                                               /gpfs/projects/AMS536/zzz.programs/dock6/parameters/flex.defn
  flex_drive_file                                              /gpfs/projects/AMS536/zzz.programs/dock6/parameters/flex_drive.tbl
  ligand_outfile_prefix                                        4qmz.virtualscreen1.minimized
  write_orientations                                           no
  num_scored_conformers                                        1
  rank_ligands                                                 no

Submit the minimization job using with the script "job-min.sh" following lines: 

  #!/bin/tcsh
  #PBS -l nodes=16:ppn=28
  #PBS -l walltime=12:00:00
  #PBS -o jobstatusmin.out
  #PBS -e jobstatusmin.err
  #PBS -N chunk0min
  #PBS -V
  #PBS -q medium
  module load intel/mpi/64/2017/0.098 
  cd $PBS_O_WORKDIR
  mpirun -np 448 /gpfs/projects/AMS536/zzz.programs/dock6/bin/dock6.mpi -i minimize-vs.in -o minimize-vs.out


At this point it would be good to rank our docked ligands by score and extract the 100 best scoring (100 lowest scoring, most negative) ligands to investigate further.

VIII. Frequently Encountered Problems

Retrieved from "https://ringo.ams.stonybrook.edu/index.php?title=2017_DOCK_tutorial_2_with_PDB_3GPL_NEW&oldid=10277"