Digitizing chemical discovery with a Bayesian explorer for interpreting reactivity data

Current attempts to automate discovery take an ad hoc approach to the problem of finding novelty (10, 37). The search in chemical space is often guided by predicted reactivity, but the desired outcome is often making a discovery, i.e., reactivity that appears unlikely according to previous knowledge. We sought to create a framework wherein expert chemists can describe their theories—including any bias from their training and experience—quantitatively as a probabilistic model. This quantitative description makes it possible to define discoveries formally in terms of the state of beliefs before and after a set of experiments.

As a very simple example of a probabilistic theory of chemistry, we postulated each compound to be capable of having one or more abstract properties to varying degrees, indicated by a number ranging continuously between 0 and 1. The assignment of these abstract properties could be equally interpreted as partitioning the compounds in chemical space into overlapping fuzzy sets (38). A prior distribution is also selected for the reactivity between each set of compounds (Fig. 2A). Combining these two distributions gives the joint probability distribution for compounds α and β to belong to mutually reactive sets A and B and react as a result (Fig. 2B). This formulation can be extended to 3 and 4 component reactions (detailed mathematical formulation in SI Appendix).

Before using our system in an experimental setup, we simulated its behavior when given artificial reactivity data. This validation step served to ensure that the output reflected the basic intuition underlying our model. We first looked at the Diels–Alder reaction, notable for motivating the formulation of pericyclic mechanisms. As the Diels–Alder reaction could not be explained using the existing properties known at the time, e.g., acid and base, the labels diene and dienophile had to be devised in order to describe the structural features of its participants. We looked at reactions within a small chemical space of molecules known at the time to test whether the model would be able to make the same deduction (Fig. 3A). A typical sample from the model after revealing the reactivity data in Fig. 3A is shown in Fig. 3B. Cyclopentadiene is seen to possess two mutually reactive properties distinct from those of all other compounds. We were encouraged by this finding since the two properties in question have direct analogues in organic chemistry, i.e., diene and dienophile. In line with the use of a conservative prior distribution for properties, only the minimum (namely four) needed to explain the reactivity observations were used, in accordance with Occam’s razor.

The assignment of molecules to sets does not preclude the use of molecular structures. It is possible to create an alternative probabilistic model that reasons about the chemical structure of molecules by representing each molecule as a vector indicating the presence or absence of certain structural features. These bit string representations, known as molecular fingerprints, are a well-studied subject within the field of cheminformatics (39), and there are several widely used algorithms available (40, 41), notably extended connectivity fingerprints based on the Morgan algorithm (42) and substructure query sets such as MACCS keys (43). Given a suitable bit vector representation, simply adapting the model to assign memberships to each fingerprint bit instead of each molecule enables reasoning about reactivity in terms of structural motifs (see SI Appendix for implementation details).

To validate the use of probabilistic reasoning about structural features using structural fingerprints, we constructed an artificial chemical space consisting of 36 two- and 84 three-component combinations between the compounds in Fig. 3C (see SI Appendix for observations in this dataset). The model was seeded with the outcome of the 36 binary reactions as initial knowledge and tasked with exploring the chemical space by randomly picking one of the remaining experiments at each step. Following this exploration phase, the model was able to infer a three-component reactivity mode, the Morgan fingerprint bits for which correspond to the isocyanide, carboxylic acid, and carbonyl motifs, i.e., the Passerini reaction. Likewise, a binary reactivity mode was inferred related to the carboxylic acid and amine groups (Fig. 3C).

By tracking the likelihood of observations—that is, how probable or unsurprising each observation is, whether reactive or nonreactive—as the model explores this chemical space, it is possible to pinpoint when anomalous observations are made and at what point these observations are interpreted as a discovery rather than anomaly (Fig. 3D). Once the empirical outcome of an experiment is known, it is possible to compare the likelihood before and after, that is, a priori versus a posteriori. The top graph in Fig. 3D shows the model’s degree of surprise to the outcomes before (indicated by the logarithm of observation likelihood) any observations have been revealed to it. As the exploration progresses (middle plot), the model revises its beliefs, accepting the Passerini reaction as predictably reactive rather than anomalous. At the end of the exploration, no observation is indicated as anomalous, meaning the model’s final interpretation is consistent with all outcomes.

Probabilistic interpretation is most versatile in combination with a robotic chemistry platform, such as the Chemputer-based setup (Fig. 4A). This platform is composed of a set of syringe pumps and valves for liquid handling, dispensing chemicals, moving the reaction mixtures, and cleaning the robot. It is capable of mixing reagents from a pool of up to 24 stock solutions to prepare reactions in 20 concurrent experiments kept under inert atmosphere (Fig. 4D).

To control the platform, we also expanded χDL⁵ (a programming language for digital execution of chemistry) to permit a nonlinear sequence of operations, collection of online analytical data such as ¹H NMR, and real-time decision-making based on the outcome of previous reactions. The current system implements real-time flow benchtop NMR, mass spectrometer, and HPLC. All instruments are remotely controlled using Python libraries that we developed to interface with manufacturer APIs using χDL steps (e.g., “Acquire HPLC”, “Shim NMR”) (Fig. 4B). Reaction temperature can also be adjusted using remote-controlled hotplates, allowing each reaction to proceed at constant temperature. χDL can perform the sequence of operations required for reaction setup and analysis in parallel, so the platform is capable of analyzing the contents of one reactor every hour (Fig. 4C).

As a starting point for validating the probabilistic approach to discovery in an experimental setting, we studied the system’s behavior while operating in a chemical space containing a set of historically significant chemical reactions (44–51). The goal was to verify whether the Oracle is able to identify these discoveries and quantify their significance without referring to prior information from the literature. To this end, nine landmark name reactions (listed in Table 1) from a range of periods in the history of synthetic chemistry were considered, and a corresponding set of 11 compounds, which participate in these reactions were selected to use as reagents in our robotic platform (Table 1). We deliberately selected a chemical space with known discoveries so that the oracle’s interpretation could be directly compared and validated against current understanding of these reactions.

A simplification used so far and commonly encountered in systems interfacing robotic platforms with machine-learning algorithms is the use of binary reactivity observations (52, 53)—that is, the outcome of each experiment is simply described as reactive or nonreactive (0 or 1). This restriction prevents reasoning in situations where the precise mode of reactivity is of interest. For instance, with reaction outcomes stored independently as binary labels, it is not possible to know whether A + B + C: 1 signifies a three-component given A + B: 1 also. The Bayesian Oracle is not inherently limited to binary observations, so in the next stage, we devised a multibit vector that allows the description of different types of reactivity. Specifically, we binned the HPLC-DAD retention times of the reactants (prior to the reaction) and final reaction mixture into a set of regions and compared. The presence of new peaks in each region is recorded as a reactivity vector used as the outcome (Table 1 and SI Appendix).

We used our system to perform and interpret the outcome of 550 reactions between two, three, or four reactants. At the start of each iteration, the probabilistic model was conditioned on the outcome of all previous reactions, and the most “disruptive” combination—that is, the combination whose outcome would most radically change expected reactivity for the remaining unexplored portion of the chemical space—was selected to perform next. To understand how these historical reactions were discovered and interpreted by the model, we asked whether they were perceived as unexpected when discovered—that is, at the point during exploration when first observed—and whether they seemed justified in retrospect once all reactions had been performed.

Examining the relative likelihood of each discovery (Fig. 5) reveals that the model can make inferences about related reactions. Based on the prior distributions used in this case, all entries are initially seen as highly unlikely. The model quickly learns about the Heck, Buchwald–Hartwig, and Sonogashira coupling reactions following the discovery of the Suzuki reaction at step 24. The Wittig reaction shows a similar dependence on the Wittig–Horner reaction (discovered at step 12); only the first to be discovered is found surprising.

The theories of chemistry demonstrated so far have been simplified to illustrate the types of insight enabled by automated Bayesian interpretation, but we envision that the greatest utility will be possible within a workflow where theories of much greater detail can be defined and evaluated by domain experts (54). To facilitate and formalize this workflow, we have created Delphi, a platform for hosting and interrogating Bayesian theories, Fig. 6A. An arbitrary number of hypotheses can be deposited in the Delphi, which assigns them unique identifiers so they can be iteratively refined and derivatized. The same set of results can be interpreted under multiple theories, providing a quantitative and objective means of evaluating the relative merit of rival theories (55–57). As a demonstration, we reinterpreted the robot’s experimental findings under an alternative theory that links structural motifs in participating reagent molecules (represented by their MACCS keys) to the number of unique product HPLC peaks in each reaction (rather than their location). The resulting interpretations are compared side by side in Fig. 6 B and C.

Founded upon recent progress in laboratory automation and online analytics, the long march toward digitizing chemical discovery has passed a number of landmark developments involving varying degrees of reliance upon and interpretability by human chemists. With this work, we demonstrate that eliminating expert input is not a necessary condition for removing hidden bias or using modern hardware to reason about reactivity in large chemical spaces. Present probabilistic methods still present computational challenges when exploring hypothesis spaces parameterized by a large number of latent dimensions or discrete parameters. Defining bespoke probabilistic models to represent domain knowledge also requires familiarity with probability theory and the inference method used. We expect the steady improvement of both inference algorithms and probabilistic programming languages to lower the computational as well as cognitive barriers to wider adoption of the probabilistic paradigm. Meanwhile, progress toward a shared standard for describing the experimental space (inputs) and predictions (outputs) will allow systems like Delphi to act as repositories for reusable expert chemical knowledge, facilitating reproducible collaboration on discovery campaigns carried out around the globe.

We thank Matthew Craven for help with implementing closed-loop experiments in χDL and Dario Cambié for help with improving the reliability of our closed-loop nuclear magnetic resonance (NMR) setup.We gratefully acknowledge financial support from the EPSRC (Grant Nos. EP/L023652/1, EP/R020914/1, EP/S030603/1, EP/R01308X/1, EP/S017046/1, and EP/S019472/1), the ERC (Project No. 670467 SMART-POM), the EC (Project No. 766975 MADONNA), and DARPA (Project Nos. W911NF-18- 2-0036, W911NF-17-1-0316, and HR001119S0003).