A Unified View of Deep Learning for Reaction and Retrosynthesis Prediction: Current Status and Future Challenges

1 Introduction

Drug discovery is crucial to human healthcare, but the process is notoriously labor-intensive and costly. As Eroom's law suggests [Scannell et al., 2012], the exploration of new drugs becomes increasingly slower and more expensive over time. Therefore, it is natural and significant to leverage machine learning techniques to accelerate the drug discovery process. In recent years, the use of deep learning approaches to enhance different stages of drug discovery has become prevalent due to the rise of deep learning. Among these stages, reaction prediction and retrosynthesis prediction are two fundamental steps that can benefit from deep learning tools.

In actual production environments, chemists often aim to design synthesis routes that can lead to target molecules through a series of chemical reactions. One common strategy is to decompose target molecules into simpler precursor structures that can be synthesized more easily, a process known as retrosynthetic analysis. Automating this planning process using deep learning is crucial for discovering and optimizing synthesis routes. Retrosynthetic planning involves two subtasks: multi-step retrosynthetic planning and singlestep retrosynthesis prediction. In this survey, we focus on the latter, as multi-step planning is usually approached as a searching problem, which is fundamentally different from the conditional structured prediction problem posed by reaction prediction and single-step retrosynthesis prediction. For the remainder of this survey, we will refer to single-step retrosynthesis prediction as simply retrosynthesis prediction. Reaction prediction is another key task in organic synthesis analysis. A robust reaction prediction model can provide insight into the underlying mechanisms of biochemical reactions, as well as generate virtual reactions to expand the database for retrosynthetic planning. In summary, reaction and retrosynthesis prediction are interrelated and can enhance one another.

∗This work is done when Ziqiao Meng worked as an intern in Tencent AI Lab.

†Corresponding authors: Peilin Zhao, and Irwin King.

Several surveys exist on reaction and retrosynthesis prediction. Engkvist [Engkvist et al., 2018] provides an overview of various computational approaches for reaction prediction, ranging from quantum computation to deep learning tools. However, this survey lacks detail for each approach and provides limited coverage of deep learning-based methods. On the other hand, [Dong et al., 2021] reviews some deep learning-based retrosynthetic planning methods and datasets, but it falls behind the state-of-the-art research on retrosynthesis prediction. It is worth noting that DualTF [Sun et al., 2021] provides a unified framework for retrosynthesis prediction from the perspective of energy-based models, but it does not offer a unified understanding of reaction and retrosynthesis, nor does it provide specific discussions on the limitations and major challenges of each approach. Overall, the existing literature lacks a comprehensive and unified understanding of advanced reaction and retrosynthesis prediction models.

Compared to previously mentioned surveys and works, our survey unifies the formulation of reaction and retrosynthesis prediction for the first time. We systematically discuss the strengths and weaknesses of each approach from different perspectives for both problems. Additionally, our survey presents novel challenges and limitations that were not explicitly stated in previous works. Finally, based on the current status, we list several future directions for further enhancement, along with detailed analysis.

2 Preliminary and ProblemFormulation

Reaction prediction and retrosynthesis prediction are dual tasks of each other. They are also known as forward reaction prediction and backward reaction prediction, respectively. To apply deep learning, both tasks are formulated as conditional generation tasks. In this section, we first introduce the problem formulation for both tasks. Then, we provide an introduction to the basic background knowledge required for understanding different methods.

Molecular Formulation. A chemical molecule M can be represented in two major data formats: the SMILES string and the molecular graph. (1) For the SMILES format, a molecular structure M is described as a sequence of characters mi such that M := m1m2...mL, where L denotes the total length of the string. The sequence represents a spanning tree of the 2D molecular structure, and each character mi denotes a structural element such as an atom element, chemical bond, branching notation, and so on. (2) A molecule can also be abstracted as an undirected graph G = {V, E}, where V = {v1, .., vn} denotes the set of n atoms and E = {e1, .., em} denotes the set of m edges. Each node is associated with a feature vector hi ∈ Rd containing atomic information such as aromaticity and electric charge. Then, we have a feature matrix H ∈ Rn×d containing allatom information. An adjacency matrix A ∈ Rn×n×c describes the topological structures of M, where Aijk indicating the presence or absence of a chemical bond of type k between atom i and atom j. Multiple molecules can be easily represented by the above two formats. For SMILES format, multiple SMILES strings can be concatenated by a full stop "." into one single SMILES sequence. For molecular graph, a set of molecules is regarded as a single disconnected graph with each molecule as an independent connected component.

Definition 1 (Reaction Prediction) Given a set of N reactant molecules {MRi }Ni=1 , the target is to predict the set of M possible product molecules {MPi }Mi=1 .

Definition 2 (Retrosynthesis Prediction) Given a set of M product molecules {MPi }Mi=1 , the target is to predict a set of N reactant molecules {MRi }Ni=1 that can lead to {MPi }Mi=1 .

Note that in real implementations, M = 1 because only the main product is recorded in public benchmark datasets. Unfortunately, this issue makes retrosynthesis prediction more difficult than reaction prediction because M < N requires retrosynthesis to attach newly appeared atoms, leading to a much larger combinatorial search space. In general, reaction prediction aims to model the conditional probability distribution P(MP |{MRi }Ni=1 ), while retrosynthesis prediction aims to model the distribution P({MRi }Ni=1 |MP ).

Reaction Center & Reaction Template. In reaction prediction, a reaction center C is a subset of atom pairs C = {(vi, vj)} ⊆ V × V that change bond types when a chemical reaction occurs. In retrosynthesis prediction, a reaction center C is defined as a subset of existing bonds C = {ei} ⊆ E that can be modified to obtain simpler structures. A reaction template pool T is a set of reaction subgraph rules derived from a large chemical reaction database. A reaction template T ∈ T is an extracted subgraph pattern from the corresponding reaction center C. Specifically, T := tR1 + tR2 + ...tRN → tP , where tRi denotes the subgraph pattern inside the ith reactant, and tP denotes the subgraph pattern inside the product.

Atom-mapping. Both reaction and retrosynthesis prediction follow the atom-mapping principle. This principle states that each atom in the reactants/products has exactly one corresponding atom in the products/reactants. This fundamental one-to-one mapping relation physically constrains reaction space and determines that chemical reactions are mainly about bond breaking and bond formation.

Synthon and Leaving Group. A target molecule can be decomposed into a set of synthons S = {G S i }Ni=1 , which are simpler precursor substructures that can constitute the target molecule with few additional bond connections. Note that the atom set V covered by S is exactly the same as the target molecule. Leaving group L = {G L i }Ni=1 is a group of atoms or substructures in original reactants that do not appear in the target molecule after a reaction occurs. In short, synthon S and leaving group L can form reactants {GR i }Ni=1 .

Evaluation Metrics. Both reaction and retrosynthesis prediction adopt top-k accuracies to evaluate the model performance. Top-k accuracy is the percentage of reactions that have the ground-truth product in the top-k predicted sets of molecules. As long as the top-k predicted products include the ground-truth main product, it would be counted as a correct prediction. Usually, the k is in range {1, 2, 3, 5, 10}.

3 Deep Learning on Reaction and Retrosynthesis Prediction

In this section, we discuss different approaches by classifying them into four categories which are template-based methods, sequence-based and graph-based autoregressive models, graph-based two-stage models, and graph-based nonautoregressive models. Specifically, we introduce their learning mechanism design decisions, weaknesses, and strengths.

3.1 Template-basedMethods

Template-based (TB) methods are mainly leveraging a pool of reaction templates to deduce possible reaction centers. TB methods are mimicking how human experts conduct chemical reasoning. Assume we have a pool of reaction templates T =

Figure 2: This figure illustrates template-based models, sequence-based autoregressive models, and graph-based autoregressive models.

{T1, T2, ...}, then the inference process of TB methods aims at picking the best template as follows:

where M P is a pool of candidate product molecules. Eq. (1) describes the steps of choosing the best template T ∗ and deducing the possible product ˆMP according to the given T ∗ . The critical part is matching molecules with templates based on some similarity measure. Specifically, P(Ti|{MRi }Ni=1 ) and P(MPi |T ∗ ) are evaluated in the following way:

where sim(·, ·) denotes a similarity measure function. It can be a simple inner product between molecular embedding and template embedding obtained from deep neural networks, or a more complicated similarity function evaluated by some subgraph matching algorithms. TB retrosynthesis has the same learning mechanism but in a reverse direction:

where M R denotes the pool of candidate reactant molecules, T ∗ := tR1 + tR2 + ...tRN → tP . The second and third equation in Eq. (3) indicates that each predicted reactant molecule ˆMRi matches with the subgraph template tRi with highest score.

Advantages. (1) TBmethods are reliable since they are using extracted human knowledge, which can always provide good interpretations for predictions. (2) The training and inference process of TB methods are relatively simple, which is easy for domain experts to manipulate.

Disadvantages. (1) The performance of TBmethods highly relies on the scale of template database. Therefore, the template database must be updated frequently, which is apparently very expensive. (2) TB methods have poor generalization to out-of-domain unseen reactions. (3) Templates are extracted local subgraph rules while ignoring global information over the reaction. Consequently, TB methods fail to capture global information interactions and easily make false predictions based on local rules. The illustration of TB methods is shown in Figure 2.

3.2 Sequence-based and Graph-based Autoregressive Models

Sequence-based autoregressive (SAR) models are widely adopted in both forward and backward prediction. It regards both problems as the neural machine translation problem. For reaction prediction, the input is SMILES strings of reactants MR := mR1 mR2 ...mRL1 with length L1 and output is SMILES strings of products MP := mP1 mP2 ...mPL2 with length L2. The input source and output target are flipped for retrosynthesis prediction. Specifically, SAR models are estimating the following conditional probability distribution:

where P(mPi |mP<i ,MR) and P(mRi |mR<i ,MP ) are approximated by the transformer [Vaswani et al., 2017] model. Each generation step does a greedy search over token space and select the optimal token. To generate top-k candidates, we only need to conduct beam search over greedy results. Graphbased autoregressive (GAR) models have similar learning mechanisms but with different generation sequence definitions. GAR models first define an action space π including several edit actions like atom addition/deletion, bond addition/deletion, and termination. A sequence of these actions would transform reactants/products into their corresponding products/reactants. Thus, they are estimating the following probability distribution:

where GR,<i and G P,<i denote the states of edited reactants and products in previous i - 1 steps respectively. The sampling process is exactly the same as that of SAR models. Each step i selects the optimal action π∗ i and applies the optimal action to transform reactants GR,i-1 and products G P,i-1 to next state GR,i and G P,i respectively. The optimal action would be "terminate" when all generation steps are predicted to be over. Note that both Eq. (4) and Eq. (5) can be simplified by Markov assumptions.

Advantages. (1) Autoregressive models do not need atommapping information. (2) Autoregressive modeling has very natural sampling process, such as beam search. (3) Sequencebased modeling can directly utilize some mature techniques from natural language processing.

Disadvantages. (1) Autoregressive modeling can only generate predictions step-by-step, which is very inefficient. (2) Autoregressive modeling requires pre-defined generation orders. However, the molecular generation order is ambiguous. (3) Sequence-based modeling requires data augmentation techniques to improve performance. The illustration of both SAR models and GAR models are shown in Figure 2.

3.3 Graph-based Two-Stage Models

Reaction. For graph-based two-stage models in reaction prediction, they split reaction prediction into two stages which are reaction center identification stage and candidate ranking stage. For reaction center identification, it aims at selecting atom pairs assigned with high reactivity scores:

where s(·) denotes any function that outputs a score in range (0, 1). Based on identified reaction centers ˆ C, a set of candidate products ˆGP = { ˆ G P 1 , ˆ G P 2 , ...} will be enumerated in a combinatorial way through hand-written rules or reaction templates. Then the second stage is learning to rank candidate products in ˆGP . To rank generated products, each pair ({G R i }Ni=1 , ˆ G P i ) should have a score evaluated as follows:

where f(·, ·) can be complex neural networks and σ(·) denotes the sigmoid function.

Retrosynthesis. For graph-based two-stage models in retrosynthesis prediction, they split retrosynthesis into two stages which are edit prediction stage and synthon completion stage. For edit prediction, they select some existing edges with high predicted reactivity scores to decompose:

After obtaining predicted edit centers ˆ C, we decompose the target molecule by breaking the predicted edit centers, which will results in a set of synthons ˆ S = { ˆ G S i }Ni=1 , (N = k + 1). A basic assumption for this approach is that the number of synthons is the same as the number of reactants. For synthon completion stage, it models the following distribution:

The above Eq. (9) indicates that the synthon completion is equivalent to attaching the leaving group G L i to the corresponding synthon G S i . Leaving group attachment can be either formulated as an autoregressive conditional generation problem or a classification problem with a predefined dictionary. With predicted leaving group ˆ G L i , each predicted reactant ˆ GR i can be recovered by attaching ˆ G L i to the corresponding synthon ˆ G S i .

Advantages. (1) The inference process of graph-based twostage models is very similar to the deduction process of chemists. (2) Graph-based two-stage models split a hard task into two simpler tasks that are easier to tackle.

Disadvantages. (1) Some two-stage models like WLDN [Jin et al., 2017] require costly hand-crafted combinatorial enumerations. (2) The overall performance of two-stage models is restricted by the bottleneck of each stage. The illustration of the above methods are shown in Figure 3.

3.4 Graph-based Non-autoregressive Models

Only reaction prediction has a graph-based nonautoregressive model achieved by NERF [Bi et al., 2021] and currently no non-autoregressive models for retrosynthesis. NERF reformulates the problem into electron redistribution modeling. Previous approaches adopt three dimensional adjacency matrix A ∈ Rn×n×c while NERF adopts twodimensional adjacency matrix by converting one-hot type encoding to scalar values such that A ∈ Rn×n . Formally, Aij is a scalar value in range [0, 3], representing the number of shared electrons (number of bonds) between atom i and atom j. Note that aromatic bond is represented as Aij = 1 and atoms i, j will be marked as aromatic atoms. The key idea of NERF is predicting ∆ˆA by combining self-attention mappings. NERF adopts the conditional variational autoencoder [Sohn et al., 2015] (CVAE) architecture to approximate P(G P |G R) by introducing a latent variable z. Instead of directly maximizing the log-likelihood logP(G P |G R), CVAE maximize its evidence-lower bound (ELBO):

where q(z|G P , GR) is the reaction encoder with reaction (G R, G P ) as input and low-dimensional representation hz as output, P(G P |G R, z) is product decoder with reactants GR and latent embedding hz as input, P(z|G R) denotes the prior distribution of latent variable z. KL term is minimizing the gap between q(z|G P , GR) and P(z|G R). The backbone network architectures of q(z|G P , GR) are combinations of GNNs and transformers. With this architecture, G R and G P will be projected to reactants embedding hR and products embedding hP respectively. A cross attention layer is used for mapping hR to latent hz with hP as teacher forcing during training. The conditional latent embedding is derived such that ˆhz = hR + hz . Then apply self-attention mechanism on ˆhz to derive two electron redistribution matrix W + and W - for bond increase and bond decrease respectively. Then ∆ˆA=W + -W - and ˆAP =AR +∆ˆA.

Advantages. (1) Non-autoregressive models enable parallel sampling, which result in much faster sampling speed compared to autoregressive ones. (2) Non-autoregressive models have already achieved the best class of top-1 accuracy, which demonstrates that the non-autoregressive decoder is very powerful for reaction modeling. (3) Non-autoregressive models do not need pre-defined generation order.

Disadvantages. (1) The uncertainty modeling of nonautoregressive models is very tricky. The top-k sampling process is not as natural as beam search in autoregressive modeling. (2) Non-autoregressive models rely on the atommapping information while this also requires additional algorithms for alignment. The illustration of NERF learning mechanism is shown in Figure 3.

4 Limitations and Challenges

In this section, we discuss some important limitations and challenges existing in current solutions.

4.1 Side Products

Reaction. Side products are missing in public USPTO benchmark datasets, which results in incomplete supervision signals. Particularly for non-autoregressive modeling, missing side products makes electron redistribution matrix violate the conservation rule and thereby the reaction space is not fully constrained. How to complete and infer this missing information is a fundamental challenge for reaction prediction.

Retrosynthesis. Retrosynthesis directly use reaction data in USPTO and thereby all single-step analysis only contains one single outcome. Missing side products in the outcome side may not affect the general process of retrosynthesis prediction. However, side products still provide important information about leaving groups.

4.2 Limitations in Dataset

Reaction. The USPTO-479K dataset has two major issues. First, reaction types are very imbalanced. Indicated by Bi [Bi et al., 2021], reactions with linear topology dominate reaction types while few reactions with cyclic topology exist in dataset. How to learn transferable knowledge from rare reactions is an important lesson for reaction modeling. Second, in real applications, same set of reactants can result in different products under different physical conditions, which is called multi-modality in reaction prediction. Multi-modality can provide rich information for conditional generative models to generate valid and diverse candidate products. However, most reactions in UPSTO-MIT are one-to-one mapping, which means that same set of reactants can only result in a unique major product.

Retrosynthesis. The USPTO-50K dataset has two major limitations. First, the size of USPTO-50K is not large enough, which only contains 50k backward reactions. Considering many recent approaches only have slight numerical differences in top-k accuracies, the current small-scale dataset is not adequate for testing model capability. Second, the current dataset will bias edit predictions and leaving group selections. Most backward reactions only have one single edit while very few have multiple edits, which results in poor prediction accuracy in multiple-edit cases. Additionally, leaving group distribution is very imbalanced, which makes graphbased models tend to select a few frequently occurring atoms.

4.3 Limitations in Evaluation

Reaction. [Kovacs et al., 2021] points out that the random split of USPTO-MIT is not sufficient enough for evaluations since most of the reaction types in the testing set are already covered in the training set. Therefore, harder dataset split like scaffold splits and time splits should be taken into account for future evaluations of reaction predictions. Scaffold split is testing whether the model can generalize well under out-of-distribution settings, in which the training data distribution is very different from testing data distribution. Time split is splitting reactions in the order of discovery time. This split aims at testing whether the model can truly discover new chemical reactions or not.

Retrosynthesis. Current evaluation metrics are not sufficient to evaluate retrosynthesis predictions since they require exact matching between predicted reactants and ground-truth reactants. This will lead to low diversity of generations. Although some predicted reactants are distinct from given ground-truth reactants, they are acceptable solutions for synthesis choices. Since the single-step retrosynthesis model is finally served for multi-step retrosynthetic planning purposes, low diversity for each step would make the whole planning monotonous and inflexible.

4.4 Challenges in Non-autoregressive Modeling

Reaction. Non-autoregressive models [Bi et al., 2021; Meng et al., 2023] achieve the state-of-the-art top-1 accuracy among template-free approaches, which demonstrates its high potential in reaction modeling. However, the top-k accuracies of these models are not comparable to that of autoregressive ones. These experimental results indicate that NERF has a powerful encoder-decoder architecture while its uncertainty modeling is not successful. They adopt the CVAE framework to inject uncertainty modeling. However, training CVAE will easily encounter the posterior collapse issue. In fact, the KL divergence loss term gets close to zero with few iterations during training. Also, imposing uncertainty on latent can lead to the uncontrollable generation of the following electron redistribution prediction, which will lower the validity of predicted products. To sum up, uncertainty estimation is an important challenge for non-autoregressive models. Non-autoregressive modeling with correct uncertainty estimation will revolutionize reaction prediction research.

Retrosynthesis. Non-autoregressive retrosynthesis model is very tricky since the total number of atoms in reactants must be pre-determined before attaching leaving groups. However, the distribution over the total number of reactant atoms is biased by the benchmark dataset. Therefore, to generate reactants in a non-autoregressive manner, enough blank entries should be reserved for predicted leaving groups. Unfortunately, if there are too many blank entries in the 3D adjacency matrix and atom feature matrix, the task would become a very tough imbalanced classification task since only a few blank entries will be filled up with predicted atoms and bonds.

4.5 Challenges in Generalization

There are two common challenges in generalization for reaction and retrosynthesis prediction. The first challenge is outof-distribution prediction. Out-of-distribution mainly evaluates whether the model trained on main reaction types can generalize to rare reaction types. This imbalance issue naturally exists in chemical reactions since we can expect that some reactions easily occur while some reactions rarely occur. Therefore, it is important to design stable reaction and retrosynthesis prediction models under different distribution splits. The second challenge comes from the low-level representation quality. Both problems require low-level molecular and reaction representation learning techniques. Therefore, it is significant to derive powerful molecular and reaction representations with strong generalization to related tasks.

5 Conclusion and Future Directions

From the above discussions, we know that current approaches still have some major shortcomings although they have achieved promising results. Therefore, we list several future directions for further refinement of current solutions.

5.1 3DMolecular Information

A chemical molecule is inherently a point cloud consisting of a set of 3D points with Cartesian coordinates. However, 3D molecular information has not been used in the literature while sequence and graph representation of molecules have been widely explored. 3D position vectors provide important complementary distance information for each pair of atoms. The relative pairwise distance can be very different in 3D Euclidean geometry compared to 2D molecular graphs. For example, atom v1 and atom v2 may be distant from each other in a non-euclidean molecular graph while they might be close to each other in 3D Euclidean space. This is particularly useful for reaction center ranking. Therefore, effectively incorporating 3D molecular information into modeling can facilitate the more accurate reaction and retrosynthesis predictions.

5.2 Diverse Benchmark Dataset and New Evaluation Metrics

For reaction prediction, the size of the current USPTO-479K dataset is large enough. However, its modality and diversity are still not sufficient. New benchmark datasets should include more reaction types and more complex reactions. Also, different dataset splits, such as scaffold split and time split, should be included for cross-validation. For retrosynthesis prediction, the USPTO-50K dataset is at small-scale. A new benchmark dataset should be a large-scale dataset containing at least 100K samples. Furthermore, future benchmark datasets should contain more target molecules with multiple edits. In addition, a new evaluation metric for retrosynthesis is necessary and urgent. FusionRetro [Liu et al., 2023] attempts to evaluate single-step retrosynthesis models in the context of multi-step planning. More diverse evaluation metrics can be designed in the future.

5.3 Non-autoregressive Modeling

Just as mentioned in the previous section, the uncertainty estimation of non-autoregressive models is still inaccurate. Thus, it is crucial to explore more effective uncertainty modeling solutions for non-autoregressive models. We believe it will be a ground-breaking moment if the top-k accuracies of nonautoregressive models can be comparable to that of autoregressive ones. Non-autoregressive retrosynthesis prediction has not even been explored in current works. Due to the small size of leaving group in dataset, the significance of inference speed have been neglected in current research. We expect that the inference efficiency will be emphasized when analyzing complicated synthesis routes.

5.4 Self-supervised Learning

As previously discussed in limitations, reaction and retrosynthesis prediction are suffering from some issues brought by the dataset itself. Leveraging self-supervised learning (SSL) strategies on unlabeled datasets to overcome these limitations is a naturally motivated direction to explore. Although various SSL strategies have been presented for general molecular representation learning, currently seldom specific SSL strategies are designed for reaction and retrosynthesis prediction. PMSR [Jiang et al., 2021] explores SSL strategies for retrosynthesis but it fails to achieve a strong performance. Therefore, exploring more powerful SSL strategies for both problems is a feasible direction.

5.5 Multi-task Learning

Combining both tasks with other related tasks is a promising direction. DualTF [Sun et al., 2021] gives a first attempt by presenting a dual model to show that reaction and retrosynthesis prediction model can reinforce each other. FusionRetro [Liu et al., 2023] and GNN-Retro [Han et al., 2022] leverages reaction contexts and molecular contexts respectively to improve retrosynthetic analysis. Lu and Zhang [Lu and Zhang, 2022] show that multi-task learning can further improve reaction prediction accuracy. Intuitively, reaction classification, reaction yield prediction, and many other tasks have close connections with reaction and retrosynthesis predictions. How to leverage these tasks together is crucial for further improvements on both problems.

Acknowledgments

The work described here was partially supported by grants from the National Key Research and Development Program of China (No. 2018AAA0100204) and from the Research Grants Council of the Hong Kong Special Administrative Region, China (CUHK 14222922, RGC GRF, No. 2151185).