- Launch MD simulation following the recipes in Subheadings
3.2and 3.3. The procedures for the complex are identical as for
the receptor alone. In the two initial energy minimization steps
of the procedure, the ligand should be considered as part of the
protein, i.e., it should be restrained. Given that ligands are very
small compared with the receptor, there will be no noticeable
differences in MD execution times. The above steps should be
repeated for all ligands obtained from the docking procedure. - The final validation of structures comes from the analysis of
MD trajectories of the studied complexes. A stable complex
should be characterized by a low root mean square deviation
(RMSD) value for the ligand throughout the entire trajectory.
However, the RMSD criterion is not sufficient to conclude
about the stability of the complex. In spite of thermal fluctua-
tions, the ligand should remain in the same position within the
binding site, showing persistent intermolecular contacts. The
information about the amplitude of atomic thermal motions
can be obtained from the root mean square fluctuations
(RMSF). The intermolecular contacts can be determined
from the analysis of atomic proximities, which also permit to
observe the formation of intermolecular hydrogen bonds. The
analysis makes use of the ptraj (or cpptraj) utility, which can be
run with the command “ptraj complex.prmtop ptraj.in”, where
the contents of the input script ptraj.in depend on the profile of
the analysis. To obtain the values of RMSD and RMSF, the
following example can be used:
trajin complex5.mdcrd 1 999999
strip :292-999999
center :1-290 mass origin
image origin center familiar
rms first mass out RMSD-rec.txt time 10 :1-290@C,N,CA
rms first mass out RMSD-lig.txt time 10 :291 nofit
atomicfluct out RMSF-lig.txt :291In this example we read all frames from the input trajectory,
then discard all components except the complex (in this example
the receptor has 290 residues and the ligand has only one residue,
#291). We superpose each frame on the first one using the mass-
weighted receptor’s backbone and calculate RMSD for the recep-
tor. Then, using the previous superposition, we compute the
ligand’s RMSD from the dispersion of all of its atoms. The time
step (10 ps) is meant to set the units on theX-axis of the generated
plot. This graph will tell us if the gravity center of the ligand
remains in the same position, but not whether the pose is stable.
If ligand moves within the binding site, its atoms will show high
values of RMSF. The last line of the above script produces a file
Molecular Dynamics in Virtual Screening 167