In the initial stage of drug discovery, hundreds of thousands of compounds are screened in high-throughput assays to find a small number of hits. To save time and resources, we develop a computational workflow that ranks drug candidates by their binding affinity to the target protein. First, all compounds in a library are docked to the protein. After this filtering step, the drug candidates with the poses scoring best are grouped into pairs, and their relative free energy of binding to the protein is calculated using the accelerated weight histogram method (AWH) for alchemical transformations. Thanks to random moves between λ-windows, a single AWH simulation samples the complete alchemical path connecting two ligands. The probability to visit a λ-window depends on its estimated free energy such that the sampling of the alchemical transformation and the free energy estimate improve each other in a positive feedback loop. To fully integrate docking and free energy calculations, we implemented routines that generate hybrid topologies for drug candidates without human intervention, i.e., they assign force field parameters (AMBER and CHARMM) to the compound and determine the smallest set of atoms to be alchemically transformed for a given compound pair. Being able to convert docking poses into input files for the molecular dynamics simulations in an automated fashion, we extensively test AWH's robustness to (sub-optimal) initial configurations on a benchmark with 13 proteins and approximately 500 ligand pairs for which experimental relative binding free energies are available. In return, free energies obtained from MD for protein-ligand complexes in the PDBbind data set are used as training data for a machine-learning model improving the scoring of docking poses such that both docking and free energy calculations benefit from integration in the workflow.
Read full abstract