eTOX ALLIES
Open Babel 2.3.9
AmberTools 15
Acpype Rev 7828
ParaDockS 1.0.1
GROMACS 4.5.5
scikit-learn 0.17
see below
Within eTOX ALLIES that was developed within our group, we use iterative Linear Interaction Energy (LIE) theory [Stjernschantz & Oostenbrink, 2010 (doi: 10.1016/j.bpj.2010.02.034); Capoferri et al, 2015 (doi: 10.1371/journal.pone.0142232)] to calculate free energies of binding.
Computational details: First, Open Babel was used to perform ligand preparation.
AmberTools was then employed to create ligand topologies and optimize ligand structures, and the topology obtained was converted to GROMACS format using Acpype.
Docking was carried out by ParaDockS (with docking radius of 10 A).
The representative (TIP3P solvated) ligand-protein complex conformations used as input structures for molecular dynamics simulations were obtained after geometric clustering of the obtained docking poses.
All MD simulations including 1 ns production were performed using GROMACS.
LIE parameters (i.e., alpha, beta, and gamma) from training (see below) were used to determine the free energies of binding of the D3R challenge compounds by using iterative LIE theory.
For every compound, the binding free energy was calculated using average ligand-environment interaction energies as obtained from the multiple MD simulations that started from the different poses, using Boltzmann-like weighting as described by Stjernschantz & Oostenbrink, 2010 (doi: 10.1016/j.bpj.2010.02.034) and replicated simulations per binding pose as described by Peric-Hassler et al., 2013 (doi: 10.3390/ijms141224514).
Average interaction energies for unbound ligands were obtained from a duplicated 1 ns production simulation of the ligand solvated in TIP3P water.
Model training: Different models were trained for benzimidazoles (alpha = 0.333, beta = 0.121, gamma = -12.995), isoxazoles (alpha = 0.101, beta = 0.100, gamma = -31.767), and sulphonamides (alpha = 0.141, beta = 0.082, gamma = -33.394).
The training sets of compounds were obtained from literature through ChEMBL [with IC50 data as obtained from doi: 10.1016/j.bmcl.2010.11.039 and doi: 10.1016/j.bmcl.2010.12.123 (benzimidazoles), doi: 10.1016/j.bmcl.2009.03.008 (isoxazoles), and doi: 10.1016/j.bmc.2014.04.014 (sulphonamides)].
Least-square fitting for model training was performed by scikit-learn packages of Python.
We used 3OMK from PDB as protein template structure for training the benzimidazole model, 3FXV for training the isoxazole model and 3BEJ for training the sulphonamide model.
For the D3R challenge compounds, we used the benzimidazole model to predict binding free energies for benzimidazole compounds and the isoxazole model to predict binding free energies for the isoxazole compounds. For the miscellaneous compounds, the model used depends on the protein conformation obtained for the miscellaneous ligand in the crystal structure provided by Roche.
For both the spiros and sulphonamide compounds, we used our sulphonamide model to predict binding free energies.
Note however that we found for most of the sulphonamide and spiro D3R compounds that their prediction fell out of the applicability domain (AD) of both our sulphonamide model and our other LIE models [see Capoferri et al, 2015 (doi: 10.1371/journal.pone.0142232) for a description of the AD analysis].
Therefore we have only reported binding free energies for 5 of the sulphonamide and spiros compounds from the D3R compound set.
For all of the 60 free energy and ranking predictions reported in the csv file in the main submission folder, we have reported our confidence in the predictions (high, medium, or low) in a separate csv file in the Supplementary Information folder.
This confidence was determined based on the AD analysis for every prediction, as described in our Capoferri paper of 2015.
Poses submitted: For every of the 60 ligands for which we provided a free energy prediction, we submitted the docked starting poses for the MD simulations with the highest Boltzmann weight in the main folder, and the poses with lower weights in the Supplementary Information folder.
eTOX ALLIES
ParaDockS 1.0.1
GROMACS 4.5.5
Binding radius = 10 A
Max number of poses to generate = 50
Mutual RMSD = 2 A
The optimized ligand structure is initially rotated by +/- 90 degrees in the x, y, or z direction.
The ligand is subsequently docked into the protein binding site, and maximally 50 poses with mutual RMSD of 2 A are retained for each of the six rotated configurations.
A principal component analysis of the docked poses (represented as heavy-atom coordinates) is performed to reduce the number of variables.
The components explaining more than 5% of the initial variance are retained, and the corresponding scores are used in subsequent k-mean clustering.
An increasing number of clusters is considered in case it would explain at least 5% more of the variance in the score space.
Fitness_function=pScore
Optimizer=PSO
Iterations=100
Docking simulations were conducted using above parameters.
The medoids of the clusters obtained were considered as representative binding poses, used as inputs for MD simulations, and included within this submission.
Note that for every of the ligands for which we provided a free energy prediction, we submitted the docked starting pose for the MD simulations with the highest Boltzmann-like weight (as obtained from MD and LIE) in the main folder, and the poses with lower Boltzmann-like weights in the Supplementary Information folder (cf. LigandScoringProtocol.txt).