Molecular docking with VEGA ZZ and Fred


1 Introduction
2 What's you need
3 FTase download
4 Protein preparation
5 Creation of the input files for NAMD
6 Run the NAMD minimization
7 Protein refinement
8 Ligand build and conformational search
9 Definition of the docking region
10 Molecular docking
11 Analysis of the docking results
12 Tips & Ticks


1 Introduction

VEGA ZZ, Fred and NAMD are very efficient tools to perform a complete molecular docking analysis. This tutorial explains in details how was performed the docking between the farnesyltransferase enzyme (FTase) and a set of inhibitors as published in the paper: Bolchi C., Pallavicini M., Rusconi C., Diomede L., Ferri N., Corsini A., Fumagalli L., Pedretti A., Vistoli G., Valoti E., "Peptidomimetic inhibitors of farnesyltransferase with high in vitro activity and significant cellular potency", Bioorg. Med. Chem. Lett., 17(22), 6192-6 (2007).


2 What's you need


3. FTase download

You can download the FTase complex (1JCQ) structure using the PDB Web interface or the tool integrated in VEGA ZZ:


4 Protein preparation

The zinc ion must be unbonded to all atoms in order to assign  the correct atom types. It results bonded to the thiolate group of the L739,750, to the ASP 297 B carboxylate, to the CYS 299 B thiolate group, to the ring nitrogen NE2 of the HIS 362 B. To remove the wrong bonds, follow these steps:

Adding the hydrogens, the thiolate groups of the L739,750 and the CYS 299 B will be protonated, but they must be without hydrogens:

The carbonyl of both C-terminal residues (HIS 367 A and GLU 424 B) are in aldehydic form and must be edited to carboxylate:

The sucrose and the acetate are far to the binding site and so can be removed:

Select all atoms (View Select All) and remove the labels (View Label atom Off).
In order to optimize the crystal structure of the complex, it will performed an energy minimization, using the CHARMM22_LIG force field keeping fixed the protein backbone. The CHARMM22_LIG contains more parameters than the standard CHARMM22_PROT for proteins.

Before the energy minimization, the atom constraints must be applied. We want to fix the protein backbone and the zinc ion, in order to avoid too many modifications in the experimental structure.


5 Creation of the input files for NAMD

NAMD requires some files: the PDB and the PSF file of the molecule, one or more parameter files and a command file defining the condition of the calculation.

numsteps                50000
minimization            on
dielectric              1.0
coordinates             1JCQ.pdb
outputname              1JCQ
outputEnergies          5000
binaryoutput            no
DCDFreq                 5000
restartFreq             5000
structure               1JCQ.psf
paraTypeCharmm          on
parameters              par_all22_na.inp
parameters              par_all22_prot.inp
parameters              par_all22_vega.inp
exclude                 scaled1-4
1-4scaling              1.0
switching               on
switchdist              8.0
cutoff                  12.0
pairlistdist            13.5
margin                  0.0
stepspercycle           20
fixedAtoms		on
fixedAtomsCol		B

Save the file with the 1JCQ_min.namd name. This file allows to perform a 50000 steps conjugate gradients minimization, saving the output (coordinates and restart files) every 5000 iterations. For more information about the parameters, please consult the NAMD User Guide.


6 Run the NAMD minimization

namd2 1JCQ_min.namd > 1JCQ_min.out

If you have more than one CPU installed, you can speed-up  the calculation specifying the total number of CPUs:

namd2 +p2 1JCQ_min.namd > 1JCQ_min.out

In this case the PC has two CPUs and both are used for the calculation (+p2 option). The dual core (e.g. Athlon X2, Pentium D, Core 2 Duo, etc) and the Pentium 4 with hyperthreading should use the +p2 option.

At the end of the calculation two files are generated: the 1JCQ.dcd (trajectory file) and the 1JCQ.out. The first one is a binary file that can't be opened with a text editor. It contains the atom coordinates of each saved frame (10 frames, because one frame every 5000 was saved). The second one is a text files containing the output messages generated by NAMD and the energy information. We need to save the lowest energy structure:


7 Protein refinement

The protein isn't ready for the docking because the ligands (FPP and L739,750) are already inside and should be removed to create enough space in the active site to dock the new ligands. In order to enlarge a little bit the cavity, the co-crystallization water molecules included in a sphere of 8 radius centred at the zinc ion coordinates, will be removed (1006, 1165, 1252, 1425, 1674, 1676, 1677,1678, 1680, 1684, 1686) because in a biological environment, a ligand could replace a water molecule and in the molecular docking is impossible to simulate it.

The protein is now ready for the docking.


8 Ligand build and conformational search

In this section will be explained how one ligand (compound 4) reported in the paper was built and how was performed the conformational analysis. This step is required by the docking software  to consider the ligand flexibility because both ligand and protein are kept rigid.

Compound 4

The subsequent conformational analysis will be performed by AMMP applying the SP4 force field parameters that are dynamically assigned starting from the bond types (single, partial double, double and triple). For this reason, the bond order must be correctly assigned, but looking the ligand structure you can see that the methionine carboxylate group has the two oxygens bonded by single and double bonds instead of two partial double bonds.

At this step, the ligand structure need an energy minimization in order to refine the 3D structure. The assignment of charges and potentials is not required because it was done automatically during the 3D conversion.

Now we perform the conformational search based on the Boltzman jump algorithm because the systematic search is too expensive in time terms due to the nine flexible torsions of the molecule.

You must define the torsion angles that will be changed during the conformational analysis:

The Boltzmann jump approach generates the conformations randomly without checking if they are similar. For this reason, the cluster analysis must be performed:

Due to the stochastic nature of the Boltzmann jump approach, running two conformational analysis of the same molecule, you could obtain different results (e.g. different number of clusters).


9 Definition of the docking region

As explained above, the docking study will be performed by OpenEye's Fred. This software requires a special file called box defining the FTase space region in which the ligand will be docked. We assume that the ligand, developed as antagonist, could occupy the catalytic pocket in which the FPP and the L739,750 are placed in the crystal structure. For this reason, the box can be defined by FPP, L739,750 atoms and the zinc ions (known to play an important role in the catalysis). To increment the docking possibilities, the box will be enlarged of 8 in all X, Y, Z dimensions (see the next section).


10 Molecular docking

To run the docking, you must have the following files in the same directory:

File name Description
1JCQ_dock.pdb FTase refined structure without the co-crystallized ligands (FPP and L739,750) and the water molecules inside a sphere of 8 radius centred on the zinc ion.
box.mol2 Atoms subset defining the FTase region in which the ligand will be docked.
compound4_clust.mol2 Ligand conformations used to simulate its flexibility.

You should have the OpenEye installation directory in the search path otherwise you must type the full path to run Fred.

fred2 -pro 1JCQ_dock.pdb -box box.mol2 -addbox 8 -dbase compound4_clust.mol2 -oformat mol2
      -conftest none -prefix 1JCQ-compound4 -refine lig_mmff


Option Argument Description
File name of the target structure in which the ligand will be docked.
File name of the box file defining the target region in which the ligand will be docked.
Expand the box size of 8 for each dimension.
Database containing all structures/conformations that will be docked.
Save the docked structures in mol2 multi model files.
Disable the database check to find redundant structures. If you don't specify this switch, one conformation only is docked.
Prefix used naming the output files.
Minimize each complex using the MMFF force field.

For a description of the Fred theory and for a more exhaustive description of the command parameters, consult the software documentation.

Fred generates a *_docked.mol2 and a *_scores.txt file for each scoring function:

File Score
ChemGauss 2
Consensus. This is not a real score function, because it's the summary of all scoring functions. The Consensus of a compound/conformation is calculated as sum of the rank position of each score function in the score list (ChemGauss 2 + ChemScore + Plp + ScreenScore + ShapeGauss). The first molecule in the consensus list is the overall best docked one considering all scoring methods.

The .mol2 files contain the ligand poses with the same order of the corresponding .txt score files starting from the best to the worst complex.


11 Analysis of the docking results

The analysis of the docking results will be performed by VEGA ZZ. As first step, you need to choose the docking pose and to assembly the complex because the mol2 files contains the ligand only without the FTase structure.

To put the ligand inside the FTase structure you must copy & paste it.

The view is a little bit chaotic because all atoms are shown. To simplify it, you could select the residues closer to the ligand:

Another analysis possibility, is the use of the VEGA ZZ Evaluation of interaction tool, that requires the atomic charges and the CVFF force field correctly assigned:


12 Tips & Tricks