An open-source molecular builder and free energy preparation workflow

Abstract

Automated free energy calculations for the prediction of binding free energies of congeneric series of ligands to a protein target are growing in popularity, but building reliable initial binding poses for the ligands is challenging. Here, we introduce the open-source FEgrow workflow for building user-defined congeneric series of ligands in protein binding pockets for input to free energy calculations. For a given ligand core and receptor structure, FEgrow enumerates and optimises the bioactive conformations of the grown functional group(s), making use of hybrid machine learning/molecular mechanics potential energy functions where possible. Low energy structures are optionally scored using the gnina convolutional neural network scoring function, and output for more rigorous protein-ligand binding free energy predictions. We illustrate use of the workflow by building and scoring binding poses for ten congeneric series of ligands bound to targets from a standard, high quality dataset of protein-ligand complexes. Furthermore, we build a set of 13 inhibitors of the SARS-CoV-2 main protease from the literature, and use free energy calculations to retrospectively compute their relative binding free energies. FEgrow is freely available at https://github.com/cole-group/FEgrow, along with a tutorial.

Keywords: Computational chemistry; Structure-based drug design.

Fig. 1. Overview of the FEgrow workflow.

(left) The user specifies the receptor, ligand core, and a list of functional groups, along with their attachment points. (centre) RDKit is used to attach the selected R-group(s) and enumerate the available conformers with a rigid core. (right) Possible bioactive conformers undergo structural optimisation using a hybrid ML/MM potential energy function. The binding affinity is predicted using a convolutional neural network scoring function and molecular properties are assessed. Final structures are output for further free energy based binding affinity assessment.

Fig. 2. Overlay of experimental and predicted protein-ligand benchmark dataset structures.

Crystal structures are shown in yellow and grown compounds in grey. a TYK2 (PDB: 4GIH), b Thrombin (PDB: 2ZFF), c P38 (PDB: 3FLY), d PTP1B with force field optimisation (PDB: 2QBS), e PTP1B using ML/MM optimisation, and f BACE(Hunt) (PDB: 4JPC). Root-mean-square distances (RMSD) between predicted and experimental coordinates of atoms in the built R-groups were calculated using RDKit.

Fig. 3. Structures of the series of cyanophenyl- and uracil-based compounds SARS-CoV-2 main protease (M^pro) inhibitors investigated here.

a Cyanophenyl-based M^pro inhibitors. b X-ray crystal structure of 4 in complex with the protease, with discussed binding pockets labelled. c, d Uracil-based M^pro inhibitors.

Fig. 4. Comparison between experimental and predicted structures of SARS-CoV-2 main protease (M^pro) inhibitors.

Overlay of a 5 and PDBID: 7L11, b 26 and 7L14, c 14 and 7L12, d 21 and 7L13. Crystal structures are coloured in yellow, and modelled binding poses in grey. Root-mean-square distances (RMSD) between predicted and experimental coordinates of atoms in the built R-groups were calculated using RDKit.

Fig. 5. Comparison between free energy calculations and experiment. Binding free energies of 13 analogues of the uracil-based M^pro inhibitors, relative to compound 10.

The error bars indicate one standard error based on least square fitting.