2023.05.08.6GJ8_dimer - tbalius/teb_docking_test_sets GitHub Wiki
Dimer docking tutorial using KRAS in complex with BI-2852 (PDB ID 6GJ8)
Here we describe docking to the 6GJ8 dimer.
For docking to the monomer, see 2023.02.24.6GJ8.
Organizing file directories
Make a directory to house the files for this tutorial:
>> mkdir 6GJ8_dimer_tutorial
Inside the created 6GJ8_tutorial directory, make subdirectories to keep files organized in:
>> mkdir 001_structure 002_surface_spheres 003_gridbox 004_dock
Dock Prep
Download the 6GJ8 PDB file from the RCSB Protein Data Bank (PDB). Save the file in the 001_structure directory.
Prepare the dimer
Open 6GJ8.pdb in Chimera.
To show the dimer, select Tools -> Higher-Order Structure -> Unit Cell in the menu bar and select "Make Copies" in the "Unit Cell" dialog box that opens. Delete the extra KRAS models by first selecting the two structures we wish to keep, select models #0 and #2 by typing select #0 #2 in the command line. The selected models should be highlighted in green. Next, invert the selection by typing select invert in the command line. Delete the extra models by selecting Actions -> Atoms/Bonds -> delete in the menu bar or by typing delete selected in the command line.
Save the structure as 6GJ8_dimer.pdb and close the session.
Run the following command in the directory containing the created 6GJ8_dimer.pdb file.
>> sed -e 's/HETATM/ATOM /g' 6GJ8_dimer.pdb | grep ^ATOM | grep -v HOH > 6GJ8_dimer_mod.pdb
The first command will change the labels of the nonstandard residues (i.e. nonprotein atoms such as the ligand and cofactor atoms) from HETATM to ATOM. The output is then piped to two grep commands that returns all the lines beginning with ATOM, excluding the lines specifying the coordinates for the water atoms. The output is saved as PDB file called 6GJ8_dimer_mod.pdb.
Open the created 6GJ8_dimer_mod.pdb to edit:
>> vim 6GJ8_dimer_mod.pdb
The fifth column (22nd character in each line) of the PDB file contains the chain identifier. Both proteins, and their corresponding ligand and cofactors, are currently labelled as the same chain (chain A). We want to label the two structures as separate chains. In the Vim editor, search for where the second protein begins by typing /GLY A 0 to search for glycine 0. There should be two residues labeled as glycine 0, one for each protein. The first glycine is the first residue for the first protein, and will be at the beginning of the PDB file. The second glycine is the first residue for the second protein, and will be in the middle of the PDB file (line 1455). Cycle through the search results until you reach the second glycine. Enter visual block mode by pressing "Ctrl+v", the bottom left corner of the terminal should display the mode as "VISUAL BLOCK". Once in visual block mode, move the cursor over to the fifth column and jump to the bottom of the file by pressing the "]" key twice. By pressing the right arrow key, change the highlighted text so that only the fifth column is highlighted. Once the fifth column for the second protein is highlighted in visual block mode, type "rB" and hit enter to replace the highlighted characters with "B", changing the chain identifier for the second protein from "A" to "B".
The excerpt below shows the 6GJ8_dimer_mod.pdb file after the changes.
...
ATOM 1451 N1 GCP A 203 16.147 29.785 1.285 1.00 7.75 N
ATOM 1452 C2 GCP A 203 15.233 30.796 1.298 1.00 10.10 C
ATOM 1453 N2 GCP A 203 15.232 31.705 2.315 1.00 9.32 N
ATOM 1454 N3 GCP A 203 14.233 30.873 0.368 1.00 8.29 N
ATOM 1455 C4 GCP A 203 14.154 29.962 -0.647 1.00 8.45 C
ATOM 1 N GLY B 0 4.779 9.509 -27.962 1.00 25.57 N
ATOM 2 CA GLY B 0 5.485 9.817 -26.724 1.00 24.56 C
ATOM 3 C GLY B 0 6.592 10.842 -26.880 1.00 26.10 C
ATOM 4 O GLY B 0 6.549 11.675 -27.794 1.00 25.83 O
ATOM 5 N MET B 1 7.586 10.791 -25.968 1.00 20.33 N
...
Save the changes to the 6GJ8_dimer_mod.pdb file and exit the Vim editor by typing :wq and hitting enter.
Prepare the receptors and cofactors
This process is virtually the same as the process for preparing the monomer, the only differences being the steps required for preparing the monomer must be simultaneously applied to both of the receptors and cofactors for the dimer.
Open 6GJ8_dimer_mod.pdb in Chimera.
Isolate the receptors and cofactors
Select the ligands by typing select #0: F0K in the command line and delete the ligands.
Change GCP to GTP
Select carbon C3B of both GCP cofactors. To select the second carbon without deselecting the first, hold "Ctrl+Shift" while clicking the second carbon after the first carbon is already selected. Change the carbons to oxygens using Build Structure. Open Build Structure by selecting Tools -> Structure Editing -> Build Structure in the menu bar. Change the mode to "Modify Structure", and change the element to oxygen, change the number of bonds to 2, set the atom names to O3B, and set the modified residue names to GTP before applying the changes.
Add hydrogens and charges
Chimera will misplace a hydrogen on both of the GTPs when adding hydrogens, the misplaced hydrogens must be corrected before charges can be added. Open Dock Prep by first selecting Surface/Binding Analysis -> Dock Prep in the menu bar. Deselect the "Add charges" and "Write Mol2 file" options so that Dock Prep stops after adding hydrogens and charges are not added. Make sure that the delete solvent, delete alternate side chain locations, and mutate incomplete side chains steps are selected in addition to adding hydrogens so that all the necessary Dock Prep steps are run. In the "Add Hydrogens for Dock Prep" dialog box that opens, select "Unspecified (determined by method)" and click OK. Select the misplaced hydrogen (H3) on both GTPs and delete them. Select the nitrogen that should be protonated (N1) on both GTPs and open Build Structure. Switch to "Modify Structure", change the element to nitrogen, change the number of bonds to 3, and apply. The newly added hydrogens are named incorrectly, select the hydrogens (H1) and reopen Build Structure. Switch to "Modify Structure", change the element to hydrogen, set the number of bonds to 1, set the atom name to HN1, and apply. After renaming, the bond length between the hydrogen and nitrogen for both GTPs should change.
Now, finish running Dock Prep by adding charges and saving the structure. Select Surface/Binding Analysis -> Dock Prep in the menu bar and select the "Add charges" and "Write Mol2 file" steps before rerunning Dock Prep. The magnesium ion should have a net charge of +2. Save the structure as kras_cof_dimer.mol2 and close the session.
Prepare the ligands
Reopen 6GJ8_dimer_mod.pdb in Chimera.
Isolate the ligands
Select the two ligands by typing select #0: F0K in the command line. Invert the selection by typing select invert in the command line. Delete the selected receptors, cofactors, and water molecules, leaving only the ligands.
Run Dock Prep
Add hydrogens and partial charges by first selecting Tools -> Surface/Binding Analysis -> Dock Prep in the menu bar and clicking OK in the "Dock Prep" dialog box that opens. Select "Unspecified (determined by method)" in the "Add Hydrogens for Dock Prep" dialog box that opens and click OK to add hydrogens. Leave the default selections and click OK in the "Assign Charges for Dock Prep" and the "Specify Net Charges" dialog boxes to assign partial charges. The ligands should have a net charge of +1.
Save both of the prepared ligands as bi2852_dimer.mol2. Delete one of the ligands and save the remaining ligand as bi2852.mol2 and close the session.
Surface spheres
Generating DMS
Open kras_cof_dimer.mol2.
Run surface calculation by selecting Actions -> Surface -> Show in the menu bar or by typing surf #0 in the command line.
Save the surface as a DMS file by selecting Tools -> Structure Editing -> Write DMS in the menu bar. Save as kras_cof_dimer.dms.
Creating surface spheres
Change directories to the 002_surface_spheres directory:
>> cd ../002_surface_spheres
Create a file called INSPH by running the command and format the file as follows:
../001_structure/kras_cof_dimer.dms
R
X
0.0
4.0
1.4
6GJ8_dimer.sph
This file will be the input file for the sphere generation program, sphgen. The first line of the file specifies the input DMS file. In the second line, R indicates that spheres should be generated outside the surface. In the third line, X indicates that all surface points should be used. The fourth line specifies the acceptable steric clash in angstroms, which in this case is set to 0.0 or no steric clash allowed. The fifth and sixth lines specify the maximum and minimum sphere radii in angstroms, respectively. The last line specifies the output sphere file name.
Run sphgen to generate spheres based on the parameters given in the INSPH file:
>> ${dock6_path}/bin/sphgen -i INSPH -o OUTSPH
- ${dock6_path} is a variable (placeholder) for the DOCK 6 installation path (e.g. /home/baliuste/zzz.github/dock6)
The surface spheres generated by sphgen can be visualized in Chimera. In a new Chimera session open 6GJ8.pdb and show the surface. In a Chimera session with kras_cof_dimer.mol2 open and the surface shown, open the output sphere file (6GJ8_dimer.sph) created by sphgen to visualize the surface spheres.
Selecting relevant spheres
To select spheres based on the distance from the ligand, run the command:
>> ${dock6_path}/bin/sphere_selector 6GJ8_dimer.sph ../001_structure/bi2852_dimer.mol2 3.0
The first two arguments passed to sphere_selector specify the surface sphere and ligand files, respectively. The last argument, 3.0, specifies that only spheres within 3.0 angstroms of the ligand should be selected.
The selected_spheres.sph file generated by sphere_selector can be visualized in Chimera.
Generating box/grid
Change directories to the 003_gridbox directory.
Box generation
Create an input file for box generation called showbox.in and format the file as follows:
Y
10.0
../002_surface_spheres/selected_spheres.sph
1
6GJ8_dimer.box.pdb
In the first line of the file, Y (Yes) indicates that showbox should create a box. The second line specifies the distance, in angstroms, between the spheres and the box borders. The third line specifies the input sphere file. The fourth line indicates the cluster number. The last line specifies the output box file name.
Generate a box based on the parameters given in showbox.in by running showbox as follows:
>> ${dock6_path}/bin/showbox < showbox.in
The 6GJ8.box.pdb file generated by showbox can be visualized in Chimera.
Grid generation
Create the input file for grid generation called grid.in:
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
allow_non_integral_charges no
bump_filter yes
bump_overlap 0.75
receptor_file ../001_structure/kras_cof_dimer.mol2
box_file 6GJ8_dimer.box.pdb
vdw_definition_file ${dock6_path}/parameters/vdw_AMBER_parm99.defn
score_grid_prefix grid
Generate grids by running:
>> ${dock6_path}/bin/grid -i grid.in -o grid.out
Docking and auxiliary calculations
Change directories to the 004_dock directory.
Single point calculation
This calculation performs an energy evaluation on the crystallographic ligand conformation without any movement.
Create the input file for the single point calculation called single_point.in:
conformer_search_type rigid
use_internal_energy yes
internal_energy_rep_exp 12
internal_energy_cutoff 100.0
ligand_atom_file ../001_structure/bi2852.mol2
limit_max_ligands no
skip_molecule no
read_mol_solvation no
calculate_rmsd yes
use_rmsd_reference_mol yes
rmsd_reference_filename ../001_structure/bi2852.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 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 ../003_gridbox/grid
gist_score_secondary no
multigrid_score_secondary no
dock3.5_score_secondary no
continuous_score_secondary no
footprint_similarity_score_secondary no
pharmacophore_score_secondary no
hbond_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 ${dock6_path}/parameters/vdw_AMBER_parm99.defn
flex_defn_file ${dock6_path}/parameters/flex.defn
flex_drive_file ${dock6_path}/parameters/flex_drive.tbl
ligand_outfile_prefix single_point.out
write_orientations no
num_scored_conformers 1
rank_ligands no
To run the single point calculation:
>> ${dock6_path}/bin/dock6 -i single_point.in -o single_point.out
Energy minimization
This calculation performs a minimization run on the crystallographic ligand conformation, allowing the ligand to subtly adjust its conformation (internal degrees of freedom) and orientation (3 translational and 3 orientational degrees of freedom).
Create the input file for energy minimization called min.in:
conformer_search_type rigid
use_internal_energy yes
internal_energy_rep_exp 12
internal_energy_cutoff 100.0
ligand_atom_file ../001_structure/bi2852.mol2
limit_max_ligands no
skip_molecule no
read_mol_solvation no
calculate_rmsd yes
use_rmsd_reference_mol yes
rmsd_reference_filename ../001_structure/bi2852.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 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 ../003_gridbox/grid
gist_score_secondary no
multigrid_score_secondary no
dock3.5_score_secondary no
continuous_score_secondary no
footprint_similarity_score_secondary no
pharmacophore_score_secondary no
hbond_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 ${dock6_path}/parameters/vdw_AMBER_parm99.defn
flex_defn_file ${dock6_path}/parameters/flex.defn
flex_drive_file ${dock6_path}/parameters/flex_drive.tbl
ligand_outfile_prefix 6GJ8_dimer.ligand.min
write_orientations no
num_scored_conformers 1
rank_ligands no
Run energy minimization:
>> ${dock6_path}/bin/dock6 -i min.in -o min.out
To view the RMSD of the minimized ligand conformation, head the generated MOL2 file:
>> head 6GJ8.ligand.min_scored.mol2
The minimized ligand conformation can be visualized in Chimera by opening the output MOL2 file with ViewDock.
[ image ]
Rigid docking
To test orienting without considering ligand flexibility, we perform rigid docking. We use the crystallographic ligand conformation for this test. This is primarily used for testing purposes only, not for screening.
Create the input file for rigid docking called rigid.in:
conformer_search_type rigid
use_internal_energy yes
internal_energy_rep_exp 12
internal_energy_cutoff 100
ligand_atom_file 6GJ8_dimer.ligand.min_scored.mol2
limit_max_ligands no
skip_molecule no
read_mol_solvation no
calculate_rmsd yes
use_rmsd_reference_mol yes
rmsd_reference_filename 6GJ8_dimer.ligand.min_scored.mol2
use_database_filter no
orient_ligand yes
automated_matching yes
receptor_site_file ../002_surface_spheres/selected_spheres.sph
max_orientations 2000
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 ../003_gridbox/grid
gist_score_secondary no
multigrid_score_secondary no
dock3.5_score_secondary no
continuous_score_secondary no
footprint_similarity_score_secondary no
pharmacophore_score_secondary no
hbond_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 ${dock6_path}/parameters/vdw_AMBER_parm99.defn
flex_defn_file ${dock6_path}/parameters/flex.defn
flex_drive_file ${dock6_path}/parameters/flex_drive.tbl
ligand_outfile_prefix rigid.out
write_orientations no
num_scored_conformers 1
rank_ligands no
Run rigid docking:
>> ${dock6_path}/bin/dock6 -i rigid.in -o rigid.out
To view the RMSD of the docked ligand conformation, head the generated MOL2 file:
>> head rigid.out_scored.mol2
Fixed anchor docking
To test the ligand flexibility algorithm without considering orienting, we perform fixed anchor docking. We use the crystallographic ligand placement for this test. This is primarily used for testing purposes only, not for screening, as a ligand in a database will be centered elsewhere, not in the site.
Create the input file for fixed anchor docking called fixed.in:
conformer_search_type flex
write_fragment_libraries no
user_specified_anchor no
limit_max_anchors no
min_anchor_size 5
pruning_use_clustering yes
pruning_max_orients 10000
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 ../001_structure/bi2852.mol2
limit_max_ligands no
skip_molecule no
read_mol_solvation no
calculate_rmsd yes
use_rmsd_reference_mol yes
rmsd_reference_filename ../001_structure/bi2852.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 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 ../003_gridbox/grid
gist_score_secondary no
multigrid_score_secondary no
dock3.5_score_secondary no
continuous_score_secondary no
footprint_similarity_score_secondary no
pharmacophore_score_secondary no
hbond_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
simplex_trans_step 1
simplex_rot_step 0.1
simplex_tors_step 10
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 ${dock6_path}/parameters/vdw_AMBER_parm99.defn
flex_defn_file ${dock6_path}/parameters/flex.defn
flex_drive_file ${dock6_path}/parameters/flex_drive.tbl
ligand_outfile_prefix 6GJ8_fad
write_orientations no
num_scored_conformers 100
write_conformations no
cluster_conformations yes
cluster_rmsd_threshold 2.0
rank_ligands no
Run fixed anchor docking:
>> ${dock6_path}/bin/dock6 -i fixed.in -o fixed.out
To view the RMSD of the top scoring docked ligand conformation, head the generated MOL2 file:
>> head 6GJ8_fad_scored.mol2
Flexible docking
This test combines orienting and ligand flexibility, and is the most applicable to production runs of dock.
Create the input file for flexible docking called flex.in:
conformer_search_type flex
write_fragment_libraries no
user_specified_anchor no
limit_max_anchors no
min_anchor_size 5
pruning_use_clustering yes
pruning_max_orients 5000
pruning_clustering_cutoff 2500
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 6GJ8_dimer.ligand.min_scored.mol2
limit_max_ligands no
skip_molecule no
read_mol_solvation no
calculate_rmsd yes
use_rmsd_reference_mol yes
rmsd_reference_filename 6GJ8_dimer.ligand.min_scored.mol2
use_database_filter no
orient_ligand yes
automated_matching yes
receptor_site_file ../002_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 ../003_gridbox/grid
gist_score_secondary no
multigrid_score_secondary no
dock3.5_score_secondary no
continuous_score_secondary no
footprint_similarity_score_secondary no
pharmacophore_score_secondary no
hbond_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 ${dock6_path}/parameters/vdw_AMBER_parm99.defn
flex_defn_file ${dock6_path}/parameters/flex.defn
flex_drive_file ${dock6_path}/parameters/flex_drive.tbl
ligand_outfile_prefix flex.out
write_orientations no
num_scored_conformers 25
write_conformations no
cluster_conformations yes
cluster_rmsd_threshold 2
rank_ligands no
Run flexible docking:
>> ${dock6_path}/bin/dock6 -i flex.in -o flex.out
To view the RMSD of the top scoring docked ligand conformation, head the generated MOL2 file:
>> head flex.out_scored.mol2
The flex.out_scored.mol2 MOL2 file can be visualized using Chimera's ViewDock.
As seen in the image above, flexible docking to the dimer produces a pose that is almost identical to the crystallographic ligand pose. The introduction of the second protein during docking allows DOCK to better position the methylimidazole group, resulting in a better pose and lower RMSD when docking to the dimer as compared to docking to the monomer.