MD refinement of docking poses


Gromacs 4.6.2, Amber 14, Acpype

System Preparation Parameters

Assumed pH 7.4
GAFF force field for ligands
FF14SB force field for proteins
TIP3P water model
Octahedral solvation box extending 8 Angstrom outside complex
Particle mesh Ewald with a cutoff of 10 Angstrom

System Preparation Method

The antechamber module of Amber 14 was used to parametrise the ligands with GAFF force field. The tLeap module integrated with Amber 14 was used to generate complex topologies where the protein was treated with ff14sb
force field. The complexes were immersed into a truncated orthorhombic box with TIP3P water molecules such that no protein atom was within 8 Angstrom of the box edge. Acpype was used to convert Amber generated topology and coordinates to Gromacs compatible files.

Pose Prediction Parameters

Standard MD simulation
Length 50 ns
Clustering using 2 Angstrom cutoff

Pose Prediction Method

Starting poses were generated by Autodock Vina (for details see "Multi-target Docking using Autodock Vina" submission), with docking attempted for all FXR structures in the PDB. For each ligand, the complex predicted to have the lowest binding free energy was selected.

Energy minimization was performed with steepest-descent for 200 steps. NPT equlibration (1 ns) was performed by using a position restraint on the protein C-alpha atoms and ligand heavy atoms. The temperature was kept at 310 K with a velocity rescaling thermostat and the pressure was kept at 1 bar with a berendsen barostat.

The restraints were then removed and the simulation was continued as long as the ligand stayed close to its starting conformation, to verify the stability of the pose. More specifically, the simulation was terminated when the ligand heavy-atom RMSD exceeded 0.25 nm from a reference configuration taken as the average over the first nanosecond of unrestrained simulation (to allow for relaxation with the new force field compared to docking), or when the simulation length reached 51 ns. All simulations were performed using Gromacs 4.6.2.

Each simulation (sampled every 20 ps) was clustered using the g_cluster utility (part of Gromacs 4.6.2), using the "Gromos" clustering method and a cutoff of 2 Angstrom. Clusters with less than 10 structures were discarded. Only for 5 ligands (1, 5, 12, 15, and 16), two clusters were obtained, whereas for the rest of the ligands, a single cluster was obtained. In the cases with two clusters, the largest cluster was ranked as the best pose and the other cluster was ranked as the second best pose (the scores were arbitrarily set to 1 and 2, assuming that a higher score denotes a better binding). The central structure of each cluster was submitted.