Our hit identification workflow combines physics-based cheminformatics methods together with novel machine learning algorithms. We employ a fragment-based virtual screening with significant speed-ups from our novel pharmacophore matching algorithm. Secondly, we enrich the pool of the potential hits with de novo generated drug-like candidates. These candidates are then ranked and refined using sequential binding affinity estimation techniques of increasing accuracy. Initial rounds of selection are done using docking scores. Subsequently, the most promising candidates are scored using MD-based methods, including MM-PBSA and FEP with QM-enhanced force fields. Additionally, we employ a Bayesian optimization process using a probabilistic machine learning model to score compounds. This model is iteratively trained on scores derived from MD and guides the selection for subsequent rounds of MD simulation.

We leverage a fragment combination approach to search a vast combinatorial space of reagent combinations while reducing the need for computationally expensive methods. We construct several pharmacophore models based on the target site, which we use for generating complementary pharmacophore embeddings. Each pharmacophore is divided into non-overlapping sub-components. Fragment-sized compounds are matched against each pharmacophore subcomponent. Hits from complementary pharmacophore sets are fed into a reaction predictor which determines whether the fragments can be combined, based on both chemical and geometrical constraints. If the combined compound maintains compatibility with the full pharmacophore this compound is passed for scoring.

We use Enamine (building blocks) fragment libraries, as these are combined by Enamine to derive the full Enamine REAL database. This way we yield compounds that are available from Enamine. Furthermore, this approach enables us to effectively explore the full 4.5bn compound space while only having to run pharmacophore matching on ~300,000 fragment compounds. This approach explores a combinatorial reaction product space that grows quadratically (in case of two-component reactions) or cubically (in case of three-component reactions) while computational costs are linear with the number of fragments.

Fragment combinations that pass the pharmacophore model are further scored by consensus docking [Autodock Vina, PLANTS, Schrodinger Glide, Molecular Operating Environment]. For the most promising pool of candidates, we compute the absolute free energy of the binding using the molecular mechanics Poisson−Boltzmann surface area (MM-PBSA) [AmberTools and Gromacs software]. Free energy perturbation calculations will be employed to choose between top-ranking ligands with similar structures and to evaluate minor modifications to the final structures.

One of the aforementioned pipelines of molecules into the binding affinity estimation funnel is provided by a variational autoencoder (VAE) trained on SMILES string representations from the Enamine Hit Locator Library. By training on this library of ~460k molecules we are able to generate novel molecules with similar characteristics to those present in Enamine. In addition to the novelty advantage, the VAE also allows us to do molecular sampling in the neural network’s latent space. This can act as a molecular scaffold similarity search and, additionally, allow interpolation between molecular structures. Once we have estimated binding affinities for molecules from any of the pipelines, we will sample the highest scoring molecules and retrain the VAE using the high-scoring subset. This has the effect of tuning the VAE’s generative distribution towards molecular structures that have higher predicted binding affinity, thereby increasing the sampling efficiency.

A key aspect of our strategy is ensuring efficient sampling of the 4.5bn compound space of the Enamine REAL database, while also ensuring that the most computationally intensive methods (MM-PBSA & FEP) are reserved for the most promising candidate molecules. Using cluster compute resources we are able to dock several million molecules: this includes carrying out ensemble docking against multiple conformations of the receptor.

The ~4 million curated set of molecules that we dock and score will be composed of the following:

• Pharmacophore matches generated from combining Enamine building blocks
• Molecules generated by the VAE
• ~3m molecules from the Enamine high-throughput screening library
• ~460k molecules from the Enamine hit locator library
• ~200k molecules from the Enamine building blocks library

Using this combination of building-block combinatorial models, ML generative models, and targeted libraries we have probed most of the relevant molecular space of the 4.5bn compound Enamine REAL database while reducing docking computations by multiple orders of magnitude.

To ensure that we maximise the information gained from the most computationally intensive MD-based scoring methods, we will follow the Bayesian optimization strategy as outlined in [Hernández-Lobato et al. 2017]. We train a Bayesian graph neural network (BGNN) to estimate the MD-derived binding affinity, we then use the predictions of the BGNN to select which molecules we should run MD scoring on next, the network is retrained with the new data and the cycle is repeated.

To ensure a diversity of structures in the final submission, we will cluster compounds based on structure and submit one sample from each cluster. FEP will be used to decide between structurally similar molecules.

• Molecular Operating Environment (MOE)
• Glide (Schrödinger)
• Q-Chem
• Gaussian

• AutoDock Vina
• Protein-Ligand ANT System (PLANTS)
• GROMACS
• Dock 3.7 (Kuntz Group UCSF)