Retrosynthesis Prediction via
Search in (Hyper) Graph

Abstract

Predicting reactants from a specified core product stands as a fundamental challenge within organic synthesis, termed retrosynthesis prediction. Recently, semi-template-based methods and graph-edits-based methods have achieved good performance in terms of both interpretability and accuracy. However, due to their mechanisms these methods cannot predict complex reactions, e.g., reactions with multiple reaction center or attaching the same leaving group to more than one atom. In this study we propose a semi-template-based method, the Retrosynthesis via Search in (Hyper) Graph (RetroSiG) framework to alleviate these limitations. In the proposed method, we turn the reaction center identification and the leaving group completion tasks as tasks of searching in the product molecular graph and leaving group hypergraph respectively. As a semi-template-based method RetroSiG has several advantages. First, RetroSiG is able to handle the complex reactions mentioned above by its novel search mechanism. Second, RetroSiG naturally exploits the hypergraph to model the implicit dependencies between leaving groups. Third, RetroSiG makes full use of the prior, i.e., one-hop constraint. It reduces the search space and enhances overall performance. Comprehensive experiments demonstrated that RetroSiG achieved competitive results. Furthermore, we conducted experiments to show the capability of RetroSiG in predicting complex reactions. Ablation experiments verified the efficacy of specific elements, such as the one-hop constraint and the leaving group hypergraph.

I Introduction

The retrosynthesis prediction aims to discern the appropriate reactants necessary for synthesizing a designated product molecule. It was initially formalized by [1] and has since evolved into a foundational issue within the realm of organic chemistry. The challenge lies in the vastness of the search space, encompassing a multitude of potential transformations, as noted by [2]. Thus, there has been a quest for a range of computer algorithms aimed at supporting seasoned chemists in their endeavors. Among them, machine learning based methods have played a vital role and achieved significant progress in this area recently [3].

The template-based and template-free methods have their own set of advantages and disadvantages. The template-based methodologies [3, 4, 5, 6, 7] tackle retrosynthesis prediction as a template retrieval problem. Following template retrieval, these approaches utilize cheminformatics tools like RDKit [8] to construct full reactions. Despite their high interpretability and assurance of molecule validity, they cannot predict reactions outside the template library. In contrast, the template-free approaches leverage deep generative models to produce reactants based on a given product directly. Since molecules can be represented using both graphs and SMILES sequences, existing techniques reframe retrosynthesis as either a sequence-to-sequence problem [9, 10, 11, 12, 13] or a graph-to-sequence problem [14]. Although these generative methods facilitate chemical reasoning within a more expansive reaction space, they lack interpretability.

The semi-template approaches effectively leverage the advantages inherent in generative models along with established chemical knowledge to enhance their methodology. Typical frameworks [15, 16, 17, 18, 19] within this category adhere to a common strategy: they first pinpoint the reaction center and convert the product into synthons. Then, a separate model completes the transformation of these synthons into reactants. These methods offer competitive accuracy and maintain interpretability due to their stage-wise nature. Nonetheless, they continue to face challenges in predicting reactions with multiple reaction center and fail in the scenario where the same leaving group is attached to more than one atom in a molecular graph. This paper refers to the reaction center comprised of multiple bonds or atoms as multiple raction center. Fig. 1 shows a reaction with multiple reaction center and attaching the same leaving group to more than one atom.

Recently, some methodologies [20, 21] have framed retrosynthesis prediction as the task of editing the product molecular graph to the reactant molecular graph. In this paper, we refer to this type of method as the graph edit sequence conditional generative model. The edit operations encompass actions such as Delete-bond, Change-bond, Change-atom, Attach-leaving-group, etc. These techniques extract edit sequences from ground-truth reaction data and then model the conditional probability distribution over these sequences given one product molecular graph. Despite its mitigation of some template-based method limitations, the graph edit sequence conditional generative model remains challenged. It cannot handle attaching the same leaving group to more than one atom in a molecular graph and also underutilizes the prior of chemical rules.

Figure 1: Architecture overview. RretroSiG first identifies reaction center via (a) search in the product molecular graph and then completes leaving groups through (b) search in the leaving group hypergraph. Finally, RetroSiG converts predicted subgraph into reaction center SMARTS and leaving groups SMARTS, subsequently merging them to derive one retrosynthesis template. RetroSiG obtains the reactants by applying the merged template to the given product. At state $t$, the red highlighted part represents the explored nodes, and the green nodes denote the action space after applying the one-hop constraint.

Motivated by the aforementioned constraints, we introduce the RetroSiG (Retrosynthesis via Search in (Hyper) Graph) framework (Fig. 1)—a semi-template-based approach—aimed at facilitating single-step retrosynthesis prediction. There are two phases in our RetroSiG, i.e., (a) search in the product molecular graph and (b) search in the leaving group hypergraph. (a) (Fig. 1a): We regard the reaction center identification as the search in the product molecular graph. A crucial observation dictates that the single or multiple reaction center should constitute a node-induced subgraph within the molecular product graph [22]. Therefore, RetroSiG utilizes a reinforcement learning agent with the graph encoder to explore the appropriate node-induced subgraph in the specific product molecular graph, identifying it as the predicted reaction center. At each stage, it contemplates the selection of one node within the molecular product graph as an action, incorporating it into the explored node-induced subgraph. As the reaction center graph maintains connectivity, we implement a one-hop constraint to streamline the search space and optimize the overall performance. (b) (Fig. 1b): We regard the leaving group completion as the search in the leaving group hypergraph. Since the unique leaving group patterns derived from the training dataset can cover 99.7% of the test set [17], it is feasible to predict leaving group patterns from the training leaving group corpus. Implicit dependencies exist among various leaving groups. Naturally, we can model these dependencies by creating a leaving group hypergraph. Each node denotes one unique leaving group pattern, while a hyperedge connects the corresponding leaving groups that co-occur within a template (reaction). Therefore, RetroSiG also utilizes a reinforcement learning agent with the hypergraph encoder to explore the appropriate node-induced subgraph in the leaving group hypergraph, selecting it as the resulting leaving groups. At each step, it contemplates the selection of one node within the leaving group hypergraph as an action, incorporating it into the explored node-induced subgraph. In the evaluation set, we observe that 98.9% of the leaving groups that co-occur within a template are contained within the same hyperedges. Thus, we make full use of the prior and impose a one-hop constraint to reduce the search space. Finally, we convert two predicted node-induced subgraphs into reaction center SMARTS and leaving groups SMARTS, subsequently merging them to derive one retrosynthesis template. We obtain the reactants by applying the retrosynthesis template to the given product [22]. Comprehensive experiments demonstrate that RetroSiG achieves a competitive result. Furthermore, we conduct experiments to show the capability of RetroSiG in predicting complex reactions. Ablation experiments confirm the efficacy of specific elements, such as the one-hop constraint and the leaving group hypergraph. In summary, we outline our primary contributions as follows:

• We propose, a semi-template-based method, the RetroSiG framework to conduct single-step retrosynthesis prediction. It has the capability to predicting reactions with multiple reaction center and attaching the same leaving group to more than one atom.

• The key novelty of RetroSiG is viewing the reaction center identification and leaving group completion as the search in the product molecular graph and the leaving group hypergraph respectively. The hypergraph effectively represents the implicit dependencies among leaving groups, and the search in the graph makes full use of the prior, i.e., one-hop constraint.

• Comprehensive experiments demonstrate that RetroSiG achieves a competitive result. Furthermore, we conduct experiments to show the capability of RetroSiG in predicting complex reactions. Ablation experiments confirm the efficacy of specific elements.

II Related Work

We summarize the comparison of existing single-step retrosynthesis in Table I. Compared with TB, RetroSiG does not need to perform Subgraph Matching (NP-Hard) since we can determine the exact location of the reaction center by the explored node-induced subgraph. Also, the specific position leads to the unique solution after applying the template, and thus RetroSiG does not need the additional model to ranking Reactants. RetroSiG can obtain Explicit Template and have the capability of predicting complex reactions with multiple reaction center (RC) or attaching the same leaving group (LG) to more than one atom.

TABLE I: Comparison of different baselines in seven dimensions.

	Interpretability	Subgraph Matching	Ranking Reactants	Explicit Template	Multi-RC	Attaching the Same LG to More Than One Atom	Extrapolation Ability
TB	✔	✔	✔	✔	✔	✔	✗
TF	✗	✗	✗	✗	✔	✔	✔
STB	✔	✗	✗	✔	✗	✗	✔
Graph2Edits	✔	✗	✗	✔	✔	✗	✔
RetroSiG	✔	✗	✗	✔	✔	✔	✔

II-A Graph Edit Sequence Conditional Generative Methods

MEGAN [20] and Graph2Edits [21] express a chemical reaction as a sequence of graph edits $\boldsymbol{E}=\left(e_{0},\cdots,e_{t},\cdots,e_{T}\right)$ that transform the product molecular graph into the reactant molecular graph. Here, the edits $e$ encompass various actions such as Delete-bond, Change-bond, Change-atom, Attach-leaving-group, and so forth. The probability of a sequence of graph edits $\boldsymbol{E}=\left(e_{0},\cdots,e_{t},\cdots,e_{T}\right)$ can be factorized to the probabilities conditional on the product molecule $\boldsymbol{p}$ and $\boldsymbol{E}_{<t}=\left(e_{0},\cdots,e_{t-1}\right)$, i.e.,

$$p(\boldsymbol{E}\mid\boldsymbol{p})=\prod_{t=1}^{T}p\left(e_{t}\mid\boldsymbol{p},\boldsymbol{E}_{<t}\right),T\geq 1$$ (1)

When $T=0$, $\boldsymbol{E}=\left(e_{0}\right)$ and $p(e_{0}\mid\boldsymbol{p})=1$ because $e_{0}$ is a fixed start token. These techniques extract edit sequences from ground-truth reaction data and then model the conditional probability distribution over these edit sequences given one product molecular graph. MEGAN [20] was the initial attempt to represent reactions as edits sequences for retrosynthesis prediction. However, this approach faced difficulties in reactant generation due to the atomic-level addition operations. It particularly struggled with reactions involving the attachment of large leaving groups. Graph2Edits [21] presents retrosynthesis as a process of reasoning involving reactions between products, intermediates, and reactants, achieved through a sequence of graph modifications. Importantly, in contrast to the MEGAN model, Graph2Edits substitutes the ”add-atom” actions with the attachment of substructures, reducing the number of generation steps and thereby enhancing reactant generation efficiency. These methods acquire knowledge of reaction transformation rules to a certain degree, thus improving their interpretability and generalization capabilities in complex reactions. However, It cannot predict reactions with attaching the same leaving group to more than one atom in a molecular graph and also underutilizes the prior of chemical rules.

II-B Template-based Methods (TB)

The template serves as a reaction representation that’s both more efficient and interpretable [22]. After establishing a library of reaction templates, the methods involve comparing a target molecule against these templates. Upon finding a match, the corresponding template is utilized to transform product molecules into reactant molecules. Numerous research efforts have introduced various strategies for prioritizing templates. In prior studies, [3] utilized molecular fingerprint similarity between the intended product and compounds in the dataset to prioritize potential templates. Meanwhile, [5] employed Neuralsym, a hybrid neural-symbolic model, to glean insights for a multi-class classification task focused on template selection. Additionally, [4] interpreted chemical insights within reaction templates as logical rules, utilizing graph embeddings to derive the conditional joint probability of these rules and reactants. Furthermore, [6] evaluated pertinent local templates within forecasted reaction centers of a target molecule. It also accounted for the non-local impacts of chemical reactions by integrating global reactivity attention. [7, 23] composes templates by selecting and annotating molecule subgraphs from training templates. Despite their interpretability, template-based methods have restricted applicability due to their incapability to forecast reactions beyond the template library. Moreover, their suitability diminishes when dealing with extensive template sets due to the challenges inherent in subgraph matching.

II-C Template-free Methods (TF)

Template-free Methods utilize deep generative models to generate reactant molecules directly. Prior research [9, 10, 11, 12, 24, 13] has employed the Transformer architecture [25], reconceptualizing the issue as a sequence-to-sequence translation task, transforming from product molecules to corresponding reactants. [14] replaces the initial sequence encoder with a graph encoder to guarantee the permutation invariance of SMILES representations. Although these generative methods facilitate chemical reasoning within a more expansive reaction space, they lack interpretability.

II-D Semi-template-based Methods (STB)

Semi-template-based approaches amalgamate the advantages of generative models with additional chemical knowledge. Considering that chemical reactions generally entail altering a limited segment of the molecular structure, most existing research approaches have divided retrosynthesis into a two-stage procedure: first, identifying the reaction center utilizing a graph neural network to generate synthons via molecular editing [26], and subsequently converting these synthons into reactants employing methods like a graph generative model [16], a Transformer [16, 18], or a subgraph selection model [17]. These methods offer competitive accuracy and maintain interpretability due to their stage-wise nature. Nonetheless, they continue to face challenges in predicting reactions with multiple reaction center. They cannot effectively address scenarios where the same leaving group is attached to more than one atom in a molecular graph [27, 28, 29, 30].

III Methodology

III-A Preliminaries

III-A1 Edge Graph Attention Network Layer

Given $\mathbf{H}^{(t)}=[\boldsymbol{h}_{1}^{(t)};\boldsymbol{h}_{2}^{(t)};\cdots;\boldsymbol{h}_{n}^{(t)}]\in R^{n\times d}$ and $\mathbf{f}_{ij}^{(t)}\in R^{1\times d}$ are the node embedding matrix and the embedding of edge $(i,j)$ at the $t$-th layer repectively, where $\boldsymbol{h}_{i}^{(t)}\in R^{1\times d}$ is the node-level embedding for node $i$ of the graph and is also the $i$-th row of $\mathbf{H}^{(t)}$, $d$ is the dimension of node-level embedding and $n$ is the number of nodes. $\boldsymbol{h}_{i}^{(0)}$ and $\boldsymbol{f}_{ij}^{(0)}$ are the initial features of node and edge in molecular product graph $\mathcal{G}_{p}$. EGAT injects the graph structure into the attention mechanism by performing masked attention, namely it only computes $\alpha_{ij}$ for nodes $j\in\mathcal{N}_{i}$, where $\mathcal{N}_{i}$ is the first-order neighbors of node $i$ in the graph:

$$\boldsymbol{f}_{ij}^{(t+1)} =\text{ LeakyReLU }\left(\left[\boldsymbol{h}_{i}^{(t)}\mathbf{W}\left\|\boldsymbol{f}_{ij}^{(t)}\right\|\boldsymbol{h}_{j}^{(t)}\mathbf{W}\right]A\right),$$ (2)

$$e_{ij} =\mathbf{a}\cdot{\boldsymbol{f}_{ij}^{(t+1)}}^{T},\alpha_{ij}=\frac{\exp\left(e_{ij}\right)}{\sum_{k\in\mathcal{N}_{i}}\exp\left(e_{ik}\right)},$$

where $e_{ij}\in R$ and $\alpha_{ij}\in R$ are a non-normalized attention coefficient and a normalized attention coefficient representing the weight of message aggregated from node $j$ to node $i$ respectively in the $t$-th layer of EGAT, and $\|$ is the concatenation operation. Besides,$\mathbf{W}\in R^{d\times d}$, $A\in R^{3d\times d}$ and $\mathbf{a}\in R^{1\times d}$ are learnable parameters in the $t$-th layer.

EGAT [31] employs multi-head attention to stabilize the learning process of self-attention, similar to Transformer [25]. If there are $K$ heads, $K$ independent attention mechanisms execute the Eq. 2, and then their features are concatenated:

$$\boldsymbol{h}_{i}^{(t+1)} =\operatorname{MLP}\left(\|_{k=1}^{K}\sigma\left(\sum_{j\in\mathcal{N}_{i}}\alpha_{ij}^{k}\boldsymbol{h}_{j}^{(t)}\mathbf{W}^{k}\right)\right)$$ (3)

where $\|$ represents concatenation, $\alpha_{ij}^{k}$ are normalized attention coefficients computed by the $k$-th learnable $\mathbf{W}^{k}\in R^{d\times d}$, $A^{k}\in R^{3d\times d}$ and $\mathbf{a}^{k}\in R^{1\times d}$ following Eq. 2. Besides, $\operatorname{MLP}$ denotes multi-perceptron.

III-A2 Hyper Graph Neural Network Layer

Specifically, for the $l$-th layer in HGNN [32], it takes hypergraph $\mathcal{G}_{hg}$’s incidence matrix $\mathrm{H}\in R^{n\times e}$ and hidden representation matrix $\mathrm{X}^{l}\in R^{n\times d^{l}}$ as input, where $n$ and $e$ are the number of nodes and hyper-edges respectively. Then, the node representations in next layer will be computed as follows:

$$\mathbf{X}^{l+1}=\sigma\left(\mathbf{D}_{v}^{-\frac{1}{2}}\mathbf{HWD}_{e}^{-1}\mathbf{H}^{\top}\mathbf{D}_{v}^{-\frac{1}{2}}\mathbf{X}^{l}\boldsymbol{\Theta}^{l}\right)$$ (4)

where $\sigma(\cdot)$ is the nonlinear activation function. $\mathbf{D}_{v}\in R^{n\times n},\mathbf{D}_{e}\in R^{e\times e},\mathbf{W}\in R^{e\times e}$ are the diagonal node degree, edge degree and edge weight matrices, respectively. $\mathbf{\Theta}^{l}\in R^{d^{l}\times d^{l+1}}$ is a trainable parameter matrix.

III-B Search in the Product Molecular Graph

In this study, RDChiral [22] is employed to identify the super general reaction center. Consequently, it ensures that whether it’s a single or multiple reaction center, it must be a node-induced subgraph within the molecular product graph. A product molecular graph $\mathcal{G}_{p}=(\mathcal{V}_{p},\mathcal{E}_{p})$ is represented as a set of $|\mathcal{V}_{p}|$ nodes (atoms) and a set of $|\mathcal{E}_{p}|$ edges (bonds). The reaction center graph $\mathcal{G}_{rc}=(\mathcal{V}_{rc},\mathcal{E}_{rc})$ is a node-induced subgraph of the product molecular graph $\mathcal{G}_{p}$ such that $\mathcal{V}_{rc}\subseteq\mathcal{V}_{p}$ and $\mathcal{E}_{rc}=\left\{(u,v)\mid u,v\in\mathcal{V}_{rc},(u,v)\in\mathcal{E}_{p}\right\}$. Given the product molecular graph $\mathcal{G}_{p}$, reaction center identification aims to detect the corresponding reaction center graph $\mathcal{G}_{rc}$, i.e. the node set $\mathcal{V}_{rc}$ of $\mathcal{G}_{rc}$. Thus, we regard the reaction center identification as the search in the product molecular graph. At each stage, RetroSiG deliberates on selecting a node within the molecular product graph as an action, incorporating it into the explored node-induced subgraph. Ultimately, we transform the identified subgraph into the SMARTS of reaction center.

State Space: RetroSiG has a discrete state space $\mathcal{S}$. Each state $s_{t}\in S$ is represented as $s_{t}=\{\mathcal{G}_{p},\hat{\mathcal{V}}_{rc}^{t}\}$, where $\hat{\mathcal{V}}_{rc}^{t}$ denotes the explored node set at step $t$. Notably, $\hat{\mathcal{V}}_{rc}^{t=0}=\emptyset$.

Action Space: RetroSiG has two types of actions in its action space $\mathcal{A}$: (1) selecting one node $v$ in the product molecular graph $\mathcal{G}_{p}$, and (2) stop action, (i.e., doing nothing), denoted as STOP. Since the reaction center graph $\mathcal{G}_{rc}$ is connected, we impose a one-hop constraint on the first type of action. This constraint serves to reduce the search space and enhance overall performance. In other words, the action $a_{t}$ is to select one node $v$ from the first-order neighbour of the current explored subgraph $\hat{\mathcal{G}}_{rc}^{t}$ or a STOP action implying stop of the exploration process. At step $t$, action space is represented as $\mathcal{A}^{t>0}=\{v|v\in\operatorname{ONE-HOP}(\hat{\mathcal{G}}_{rc}^{t})\}\cup\{\text{STOP}\}$. Notably, $\mathcal{A}^{0}=\{v|v\in\mathcal{V}_{p}\}$.

Figure 2: a. Policy network architecture in the search in the product molecular graph. b. Policy network architecture in the search in the leaving group hypergraph.

Reward Function: At intermediate step $t$, reward function $\mathcal{R}$ gives $s_{t}$ a reward 0. At the terminal step $T$, if the predicted $\hat{\mathcal{V}}_{rc}^{T}$ is identical to the ground-truth $\mathcal{V}_{rc}$, $\mathcal{R}$ gives $s_{t}$ a reward 1. Otherwise, it gives $s_{t}$ a reward of 0.

Policy Net (Fig. 2a): We use the EGAT [31] to obtain node-level embedding. Each atom $u$ has an initial feature vector $\mathbf{x}_{u}$. Each bond $(u,v)$ has a feature vector $\mathbf{x}_{uv}$. The specifics regarding the initial features can be located in Table II and Table III. For ease of reference, we symbolize the encoding process as $\mathrm{EGAT}(\cdot)$ and elaborate on the architectural specifics in Section III-A. The EGAT computes atom representations $\left\{\mathbf{c}_{u}^{t}\mid u\in\mathcal{V}_{p}\right\}$ via

In order to incorporate information about the explored node set $\hat{\mathcal{V}}_{rc}^{t}$ at step $t$, we add a learnable embedding $\mathbf{x}_{explored}$ to the initial feature of the explored node. $\mathbb{I}[u\in\hat{\mathcal{V}}_{rc}^{t}]\mapsto\{0,1\}$ be the predicate that indicates whether atom $u$ is the element of $\hat{\mathcal{V}}_{rc}^{t}$. Here, $\mathcal{N}(u)$ refers to the neighboring atoms of atom $u$. The graph representation $\mathbf{c}_{\mathcal{G}}^{t}$ is formed by aggregating atom representations, i.e. $\mathbf{c}_{\mathcal{G}}^{t}=\sum_{u\in\mathcal{V}_{p}}\mathbf{c}_{u}^{t}$. Finally, we calculate the logits $s_{u}^{t}\in\mathbb{R}^{1},s_{s}^{t}\in\mathbb{R}^{1}$ for selecting one atom $u$ and STOP action at each step $t$ through the fully connected layers:

$$s_{u}^{t} =\boldsymbol{W}_{b}\left(\sigma\left(\boldsymbol{W}_{a}\mathrm{c}_{u}^{t}+\boldsymbol{B}_{a}\right)\right)-\mathbb{I}\left[u\notin\mathcal{A}^{t}\right]\infty,$$ (6)

$$s_{s}^{t} =\boldsymbol{W}_{b}\left(\sigma\left(\boldsymbol{W}_{a}\mathrm{c}_{\mathcal{G}}^{t}+\boldsymbol{B}_{a}\right)\right)-\mathbb{I}\left[\mathrm{STOP}\notin\mathcal{A}^{t}\right]\infty,$$

where $\boldsymbol{W}_{a}$ and $\boldsymbol{W}_{b}$ are the weights and $\boldsymbol{B}_{a}$ is the bias. $\mathbb{I}[u\notin\mathcal{A}^{t}]\,(\mathbb{I}[\text{STOP}\notin\mathcal{A}^{t}])\mapsto\{0,1\}$ be the predicate that indicates whether atom $u$ (STOP) is not in the action space $\mathcal{A}^{t}$ at step $t$. We add the $-\infty$ to the logit of the action not appearing in the action space $\mathcal{A}^{t}$ at step $t$. It ensures that the probability of an action not appearing in the action space equals 0 after $\mathrm{SOFTMAX}$ operation.

TABLE II: Atom Features used in EGAT. All features are one-hot encoding.

Feature	Description	Size
Atom type	Type of an atom by atomic number.	100
Total degree	Degree of an atom including Hs.	6
Explicit valence	Explicit valence of an atom.	6
Implicit valence	Explicit valence of an atom.	6
Hybridization	sp, sp2, sp3, sp3d, or sp3d2.	5
# Hs	Number of bonded Hydrogen atom.	5
Formal charge	Integer electronic charge assigned to atom.	5
Aromaticity	Whether atom is part of aromatic system.	1
In ring	Whether an atom is in ring	1

TABLE III: Bond features used in EGAT. All features are one-hot encoding.

Feature	Description	Size
Bond type	Single, double, triple, or aromatic.	4
Conjugation	Whether the bond is conjugated.	1
In ring	Whether the bond is part of a ring.	1
Stereo	None, any, E/Z or cis/trans.	6
Direction	The direction of the bond.	3

III-C Search in the Leaving Group Hypergrah

In this paper, we use RDChiral [22] to extract the super general leaving groups from the training set and obtain unique leaving group patterns. Since the unique leaving group patterns derived from the training dataset can cover 99.7% of the test set [17], it is feasible to predict leaving group patterns from the training leaving group corpus. Implicit dependencies exist among various leaving groups. Naturally, we build a hypergraph to model these dependencies. Each node denotes one unique leaving group pattern, while a hyperedge connects the corresponding leaving groups that co-occur within a template. A leaving group hypergraph $\mathcal{G}_{hg}=(\mathcal{V}_{hg},\mathcal{E}_{hg})$ is represented as a set of $|\mathcal{V}_{hg}|$ nodes (unique leaving group) and a set of $|\mathcal{E}_{hg}|$ hyperedges (co-occurrence). Given the product molecular graph $\mathcal{G}_{p}$, leaving group completion aims to identify the corresponding leaving groups, i.e. the sub-node set $\mathcal{V}_{lg}$ of $\mathcal{V}_{hg}$. Thus, we regard the leaving group completion as the search in the leaving group hypergraph. RetroSiG considers choosing one node in the leaving group hypergraph and adding it to the explored node set as an action at each step. In the end, we transform the identified nodes into the SMARTS of leaving groups.

State Space: Similar to search in the product molecular graph, each state $s_{t}\in S$ is represented as $s_{t}=\{\mathcal{G}_{p},\mathcal{G}_{hg},\hat{\mathcal{V}}_{lg}^{t}\}$, where $\hat{\mathcal{V}}_{lg}^{t}$ denotes the explored node set at step $t$. Notably, $\hat{\mathcal{V}}_{lg}^{t=0}=\emptyset$.

Action Space: Here, RetroSiG also has two types of actions in its action space $\mathcal{A}$: (1) selecting one node $v$ in the leaving group hypergraph $\mathcal{G}_{hg}$, and (2) stop action, (i.e., doing nothing), denoted as STOP. In the evaluation set, 98.9% of the leaving groups that co-occur within a template are contained within the $\mathcal{E}_{hg}$. It implies that, during inference, the prediction is likely to be incorrect if the chosen leaving group originates from different hyperedges. Thus, we impose a one-hop constraint on the first type of action. This constraint serves to reduce the search space and enhance overall performance. In other words, the action $a_{t}$ is to select one node $v$ from the first-order neighbour of the current explored node set $\hat{\mathcal{V}}_{lg}^{t}$ or a STOP action. At step $t$, action space is represented as $\mathcal{A}^{t>0}=\{v|v\in\operatorname{ONE-HOP}(\hat{\mathcal{V}}_{lg}^{t}\}\cup\{\text{STOP}\}$. Notably, $\mathcal{A}^{0}=\{v|v\in\mathcal{V}_{hg}\}$.

Reward Function: At intermediate step $t$, reward function $\mathcal{R}$ gives $s_{t}$ a reward 0. At the terminal step $T$, if the predicted $\hat{\mathcal{V}}_{lg}^{T}$ is identical to the ground-truth $\mathcal{V}_{lg}$, $\mathcal{R}$ gives $s_{t}$ a reward 1. Otherwise, it gives $s_{t}$ a reward of 0.

Policy Net (Fig. 2b): We use the HGNN [32] to obtain node-level embedding. We assign each atom $u$ in the hypergraph a randomly initialized learnable embedding $\mathbf{f}_{u}$. To simplify, we represent the encoding process as $\mathrm{HGNN}(\cdot)$ and elaborate on the architectural specifics in Section III-A. The HGNN computes atom representations $\left\{\mathbf{h}_{u}^{t}\mid u\in\mathcal{V}_{hg}\right\}$ via

$$\left\{\mathbf{h}_{u}^{t}\right\}=\operatorname{HGNN}\left(\mathcal{G}_{hg},\left\{\mathbf{f}_{u}+\mathbb{I}[u\in\hat{\mathcal{V}}_{lg}^{t}]\mathbf{f}_{explored}\right\}\right).$$ (7)

In order to incorporate information about the explored node set $\hat{\mathcal{V}}_{lg}^{t}$ at step $t$, we add a learnable embedding $\mathbf{f}_{explored}$ to the initial embedding of the explored node. $\mathbb{I}[u\in\hat{\mathcal{V}}_{lg}^{t}]\mapsto\{0,1\}$ be the predicate that indicates whether atom $u$ is the element of $\hat{\mathcal{V}}_{lg}^{t}$. The graph representation $\mathbf{h}_{\mathcal{G}}^{t}$ is an aggregation of atom representations, i.e. $\mathbf{h}_{\mathcal{G}}^{t}=\sum_{u\in\mathcal{V}_{hg}}\mathbf{h}_{u}^{t}$. We also use the EGAT [31] to obtain graph-level embedding $\mathbf{c}_{\mathcal{G}}$ of the product molecular graph $\mathcal{G}_{p}$. Finally, we calculate the logits $s_{u}^{t}\in\mathbb{R}^{1},s_{s}^{t}\in\mathbb{R}^{1}$ for selecting one atom $u$ and STOP action at each step $t$ through the fully connected layers:

$$s_{u}^{t} =\boldsymbol{W}_{b}\left(\sigma\left(\boldsymbol{W}_{a}\left(\mathbf{h}_{u}^{t}+\mathbf{c}_{\mathcal{G}}\right)+\boldsymbol{B}_{a}\right)\right)-\mathbb{I}\left[u\notin\mathcal{A}^{t}\right]\infty,$$ (8)

$$s_{s}^{t} =\boldsymbol{W}_{b}\left(\sigma\left(\boldsymbol{W}_{a}\left(\mathbf{h}_{\mathcal{G}}^{t}+\mathbf{c}_{\mathcal{G}}\right)+\boldsymbol{B}_{a}\right)\right)-\mathbb{I}\left[\mathrm{STOP}\notin\mathcal{A}^{t}\right]\infty,$$

where $\mathbb{I}[u\notin\mathcal{A}^{t}]\,(\mathbb{I}[\text{STOP}\notin\mathcal{A}^{t}])\mapsto\{0,1\}$ be the predicate that indicates whether atom $u$ (STOP) is not in the action space $\mathcal{A}^{t}$ at step $t$. The $-\infty$ operation ensures that the probability of an action not appearing in the action space equals 0 after $\mathrm{SOFTMAX}$ operation.

III-D Training and Inference

The goal of policy network aims to maximize the expected reward, enhance search quality across episodes and it is updated using Proximal Policy Optimization (PPO) [33]:

$$L^{CLIP}(\theta)=\hat{E}_{t}[\min(r_{t}(\theta)\hat{A}_{t},\operatorname{clip}\left(r_{t}(\theta),1-\epsilon,1+\epsilon\right)\hat{A}_{t})].$$ (9)

Here, $\hat{E}_{t}$ represents the expected value at timestep t, $r_{t}(\theta)$ is the ratio of the new policy and the old policy, and $\hat{A}_{t}$ is the estimated advantage at timestep t.

Inference is executed via beam search [34] employing a log-likelihood scoring function. For a beam size $k$, we can obtain two sets of $k$ candidate solutions from two stages. Finally, we have $k^{2}$ candidate solutions. From the $k^{2}$ possibilities, we select $k$ candidate solutions with the highest cumulative log-likelihood score as the top-k prediction results.

IV Evaluation

IV-A Data

We utilize the well-established retrosynthesis benchmark dataset, USPTO-50K [35], as the basis for evaluating our approach. This dataset comprises 50,016 reactions with accurate atom-mapping, categorized into ten distinct reaction types. Our methodology aligns with the data partition detailed in [3], wherein we allocate 40,000 reactions for training, 5,000 for validation, and 5,000 for testing. To mitigate potential information leakage inherent in the USPTO-50K dataset, as highlighted in previous studies [17, 15], we apply additional steps. This involves canonicalizing the product SMILES and reassigning the atom-mapping to the reactant atoms, following the methodology outlined in [17].

IV-B Evaluation Metric

In line with previous studies, we adopt top-k exact match accuracy as our evaluation metric. In our experiments, we utilized different values of k—1, 3, 5, and 10—for comparative analysis. This metric assesses whether the predicted set of reactants precisely matches the ground truth reactants.

IV-C Baseline

We contrast the predictive outcomes of RetroSiG against various template-based, template-free, and semi-template-based methodologies. Notably, we introduce semi-template-based methods including G2G [16], RetroXpert [15], RetroPrime [18], MEGAN [20], GraphRetro [17], and Graph2Edits [21] as primary baselines due to their exceptional performance. We also considered template-based models like Retrosim [3], Neuralsym [5], GLN [4], LocalRetro [6], and template-free models such as SCROP [10], Augmented Transformer [11], GTA [12], Graph2SMILES [14], Dual-TF [36], Retroformer [13] as robust baseline models for comparison.

IV-D Implementation Details

[37], DeepHypergraph (DHG), and PyTorch [38]. In the EGAT framework, we configure four identical four-head attentive layers with a hidden dimension of 512, while the HGNN structure consists of two layers with a hidden dimension of 512. Embedding sizes across the model are uniformly set at 512. We use $\operatorname{ELU}(x)=\alpha(\exp(x)-1)$ for $x\leq 0$ and $x$ for $x>0$ as our activation function where $\alpha=1$ as our activation function. Our experiments are executed on a machine equipped with an Intel Xeon 4114 CPU and a Nvidia Titan GPU. Training parameters include a learning rate of 0.0001, 10,000,000 iterations, and the utilization of the Adam optimizer [39]. At intervals of 1000 iterations, we store checkpoints to identify the best-performing ones based on evaluation criteria.

IV-E Performance

Main Results: Table IV displays the top-k exact match accuracy outcomes obtained from the USPTO-50k benchmark. When the reaction class is unknown, our method achieves a 54.9% top-1 accuracy which is competitive and only 0.2% lower than the SOTA (Graph2Edits) result of 55.1%. However, for k equals 3 and 5, RetroSiG outperforms the performance of Graph2Edits. RetroSiG also beats LocalRetro by a margin of 1.5% and 0.1% respectively in top-1 and top-3 accuracy. With the reaction class given, our method achieves a 66.5% top-1 accuracy which is only 0.6% lower than the SOTA (Graph2Edits). For k equals 3, 5 and 10, RetroSiG outperforms the performance of Graph2Edits. We attribute it to the solid prior knowledge carried by the one-hop mechanism.

TABLE IV: Top-$k$ accuracy for retrosynthesis prediction on USPTO-50K.

Model	Top-k accuracy (%)
	Reaction class unknown				Reaction class known
	k=1	3	5	10	1	3	5	10
	Template-Based Methods
Retrosim	37.3	54.7	63.3	74.1	52.9	73.8	81.2	88.1
Neuralsym	44.4	65.3	72.4	78.9	55.3	76.0	81.4	85.1
GLN	52.5	69.0	75.6	83.7	64.2	79.1	85.2	90.0
LocalRetro	53.4	77.5	85.9	92.4	63.9	86.8	92.4	96.3
	Template-Free Methods
SCROP	43.7	60.0	65.2	68.7	59.0	74.8	78.1	81.1
Aug.Transformer	53.2	-	80.5	85.2	-	-	-	-
GTA	51.1	67.6	74.8	81.6	-	-	-	-
Graph2SMILES	52.9	66.5	70.0	72.9	-	-	-	-
Dual-TF	53.6	70.7	74.6	77.0	65.7	81.9	84.7	85.9
	Semi-Template-Based Methods
G2G	48.9	67.6	72.5	75.5	61.0	81.3	86.0	88.7
RetroXpert	50.4	61.1	62.3	63.4	62.1	75.8	78.5	80.9
RetroPrime	51.4	70.8	74.0	76.1	64.8	81.6	85.0	86.9
MEGAN	48.1	70.7	78.4	86.1	60.7	82.0	87.5	91.6
GraphRetro	53.7	68.3	72.2	75.5	63.9	81.5	85.2	88.1
Graph2Edits	55.1	77.3	83.4	89.4	67.1	87.5	91.5	93.8
RetroSiG	54.9	77.6	84.1	89.0	66.5	87.9	92.0	94.1

Figure 3: Analysis of complex samples.

Complex Sample Results: We investigate the performance effect of some complex reactions in the USPTO-50k, including reactions with multiple reaction center and attaching the same leaving group to more than one atom. As is shown in Fig. 3a, Most reactions with 1, 2 atoms in the reaction center account for the majority, which have 1121 (22.3%) and 3604 (72.0%) cases respectively. However, the reactions with 3, 4, 5, 6, and $\geq$7 account for a small proportion, which have 184 (3.7%), 29 (0.6%), 24 (0.5%), 36 (0.7%) and 9 (0.2%) cases respectively. Multiple reaction center consists of three or more atoms, and 92.6% (261 / 282) of the samples with multiple reaction center are with attaching the same leaving group to more than one atom. Thus, we report the top-10 accuracy of the number of atoms 3, 4, 5, 6, and $\geq$7 in Fig. 3b to investigate the performance of complex reactions. We observe that our model’s performance remains stable even as the number of atoms increases. It demonstrates that our model can handle these complex reactions.

TABLE V: Ablation Study.

	Hypergraph Structure	One-Hop Constraint(rc)	One-Hop Constraint(lg)	Top-1(%)
1⃝	✔	✔	✔	54.9
2⃝	✔	✔	✗	54.6
3⃝	✔	✗	✔	53.5
4⃝	✔	✗	✗	53.4
5⃝	✗	✗	✗	52.8

Ablation Study: In Fig. V, Hypergraph-Structure ✔denotes that we use HGNN and hypergraph structure to encode each node (leaving group) embedding. One-Hop Constraint ✔denotes that we use a one-hop constraint. 1⃝ outperforms 2⃝ and 3⃝, showing that one-hop constraint is effective. 2⃝ outperforms 3⃝, demonstrating that one-hop constraint in reaction center identification contributes more to performance improvement. 4⃝ outperforms 5⃝, implying that hypergraph effectively models dependencies between leaving groups.

Figure 4: Visualization of Example Predictions.

IV-F Visualization of Example Predictions

In Fig. 4, we illustrate the model predictions alongside the corresponding ground truth for three specific cases. Fig. 4a shows an example with attaching the same leaving group to more than one atom. Moreover, the model identifies both the reaction center and leaving groups correctly. Fig. 4b, the model identifies the reaction center correctly, while the predicted leaving groups are incorrect. We hypothesize this is because carboxyl (COOH) is more common in the training set than acyl chloride (COCl). However, the prediction results of top-1 are also feasible chemically. Fig. 4c, the model incorrectly identifies the reaction center and consequently the leaving group, resulting in top-1 prediction being inconsistent with ground truth. Top-1 prediction is an open-loop operation, which is not necessarily harmful from the perspective of multi-step retrosynthesis. Interestingly, the reaction of top-1 prediction is to attach more than one leaving group to one synthon. Previous semi-template-base methods [17] could not predict such complex samples.

V Discussion and Conclusion

Limitations: One common drawback among all semi-template-based and graph-edit methods is their heavy dependence on atom-mapping details between reactants and products. Incorrect mapping will extract the wrong template. It will provide incorrect reaction center labels and leaving group labels during model training. We look forward to the continuous improvement of the atom-mapping algorithm and the extraction/application template algorithm of Rdchiral [22] in the community, which is very beneficial to the RetroSiG framework.

Future Works: There are three directions for future work. End2End: We need to design an end-to-end architecture to connect the two search processes so that the signal at the final step can update both stages. Decision Transformer : Our current policy net consists of a simple (hyper) graph encoder and MLP; hence, we must design a more powerful policy net. The ability of Transformer [25] has recently been proven in various fields, especially in large language models [40]. Naturally, Decision Transformer [41] is a choice as a policy net. Termination Reward: Upon completion, the environment currently incentivizes the agent solely when its prediction aligns precisely with the ground truth. However, there exist various valid methods to synthesize a product. It is crucial to devise a more reasonable reward mechanism at termination step.