SynFlowNet: Design of Diverse and Novel Molecules with Synthesis Constraints

## Abstract

Abstract Generative models see increasing use in computer-aided drug design. However,
while performing well at capturing distributions of molecular motifs, they often produce synthetically inaccessible molecules.
To address this, we introduce SynFlowNet, a GFlowNet model whose action space uses chemical reactions and buyable reactants to sequentially build new molecules. By incorporating forward synthesis as an explicit constraint of the generative mechanism, we aim at bridging the gap between in silico molecular generation and real world synthesis capabilities. We evaluate our approach using synthetic accessibility scores and an independent retrosynthesis tool to assess the synthesizability of our
compounds, and motivate the choice of GFlowNets through considerable improvement in sample diversity compared to baselines.
Additionally, we identify challenges with reaction encodings that can complicate traversal of the MDP in the backward direction.
To address this, we introduce various strategies for learning the GFlowNet backward policy and thus demonstrate how additional constraints can be integrated into the GFlowNet MDP framework. This approach enables our model to successfully identify synthesis pathways for previously unseen molecules.

## 1 Introduction

Designing molecules with targeted biochemical properties is a critical challenge in drug discovery, where computational models could play a significant role to increase efficiency and effectiveness.
Recently, generative models have lead to a renaissance in de novo molecular design (Stanley & Segler, ; Du et al., ).
However, most current de novo design models do not explicitly account for synthetic accessibility. Many approaches operate on SMILES strings (Segler et al., ; Gupta et al., ) or assemble molecules by composing atoms or molecular fragments into a graph (Lewis & Dean, ; Jensen, ; Jin et al., ), leaving no guarantee that the sampled molecules can be synthesized (Stanley & Segler, ; Gao & Coley, ).

To address this limitation, synthetic complexity scores (Ertl & Schuffenhauer, ; Coley et al., ; Liu et al., ) have been proposed as a method to assess molecules and complement generative models with knowledge on synthetic accessibility, however these heuristics and learned metrics are often oversimplified. Alternatively, computer-aided synthesis planning can can be performed as a subsequent step to molecule generation (Segler et al., ; Schwaller et al., ). While synthesis planning can provide principled routes, it can take seconds to minutes to propose pathways, which hampers its integration as a reward function in the generation process (Liu et al., ). Finally, synthetically accessible chemical spaces can be assembled from combinations of reactions involving readily available reactants (Klarich et al., ), but virtual screening on these vast datasets is already impractical today (screening only a fraction of such libraries amounts to thousands of years of compute on a single CPU; Sadybekov et al., ) and these spaces are still growing exponentially.

Given the challenges of synthesizability in atom- and fragment-based generative models, we propose formulating molecule generation in an action space of chemical reactions which naturally enhances synthesizability.
However, simply employing a synthesizable space is not sufficient. To be applicable to drug discovery programs, such a model needs to maximize desired properties and retrieve molecules from distinct modes of the available synthetic space. Generative Flow Networks (GFlowNets; Bengio et al., ) have emerged as a framework capable of generating diverse samples and requiring fewer evaluations of the reward function compared to alternatives like Markov Chain Monte Carlo or Proximal Policy Optimization (Schulman et al., ).

In this work, we introduce SynFlowNet, a GFlowNet specifically trained to generate molecules from documented chemical reactions and purchasable starting materials, thereby constraining exploration to a synthetically accessible chemical space and sampling not only target compounds but also the synthetic routes leading to them. Our main contributions are:

Figure: Figure 1: SynFlowNet allows for synthesis-aware molecule generation. (A) The state space is induced by combining purchasable building blocks (BBs) and chemical reactions. Every final molecule (rectangle box) is associated with a reward. Training trajectories are constructed by sampling the model forward. Our policy $P_{F}(a|s)$ is parameterized as a graph transformer which at each timestep processes the current molecular state $s_{t}$ and outputs a shared embedding which is passed to separate MLP heads to predict the action logits for different action types (5 forward and 3 backward action types). An action $a_{t}$ is then sampled from this hierarchical distribution to transition to the next state via a reaction. Newly added BBs are highlighted in yellow. (B) To allow handling large sets of BBs (up to 200k), we represent them using Morgan fingerprints and compute the probability of sampling a particular BB from the normalised dot product between this representation and the embedding vector of the current state. (C) Finally, when traversing the MDP backwards, to reduce the probability of exiting the MDP defined by our set of reactions and BBs, we train the backward policy to avoid paths that do not terminate in the initial state $s_{0}$.
Refer to caption: x1.png

- •
We train a GFlowNet using an action space defined by documented chemical reactions and purchasable starting materials to generate synthesizable molecules.
- •
We show the advantage of using a reaction-based environment over a fragment-based one in terms of synthesisability, and the benefits of GFlowNets over Reinforcement Learning (RL) in terms of sample diversity on several targets. Additionally, when comparing to other generative methods with synthetically-accessible outputs, we show that SynFlowNet generates more novel candidates w.r.t known drug-like molecules.
- •
We evaluate different action representation alternatives to allow efficient scaling of the action space to up to 200k chemical building blocks.
- •
We identify an inherent problem of employing a reaction-based Markov Decision Process (MDP) with GFlowNets, stemming from the lack of guarantee that backward-constructed trajectories can return to $s_{0}$. To resolve this issue, we propose one of the first attempts at training the backward policy in a GFlowNet with a separate objective from for the forward policy, which we show to correct backward-generated trajectories and improve sample quality and diversity.
- •
We show that our proposed framework can be integrated with target-specific experimental data to inform selection of building blocks, which improves the efficiency of our model.

## 2 Background and Related Work

### 2.1 GFlowNets

GFlowNets (Bengio et al., ) are a class of probabilistic models that learn a stochastic policy to generate objects $x$ through a sequence of actions, with probability proportional to a reward $R(x)$. The sequential construction of objects $x$ can be described as a trajectory $\tau\in\mathcal{T}$ in a directed acyclic graph (DAG) $G=(S,\mathcal{E})$, starting from an initial state $s_{0}$ and using actions $a$ to transition from a state to the next: $s\rightarrow s^{\prime}$. A GFlowNet uses a forward policy $P_{F}(-|s)$, which is a distribution over the children of state $s$, to sample a sequence of actions based on the current states. Similarly, a backward policy $P_{B}(-|s)$ is the distribution over the parents of state $s$, and can be used to calculate probabilities of backward actions, leading from terminal to initial states. The training objective which we adopt in this paper is trajectory balance (Malkin et al., ):

$$ $\mathcal{L}_{\text{TB}}(\tau)=\left(\log\frac{Z_{\theta}\prod_{t=1}^{T}P_{F}(s _{t+1}|s_{t};\theta)}{R(x)\prod_{t=1}^{T}P_{B}(s_{t-1}|s_{t};\theta)}\right)^{2}$ (1) $$

used to learn the forward and backward policies $P_{F}(-|s;\theta)$ and $P_{B}(-|s;\theta)$ parameterized by $\theta$ and to estimate the partition function $Z_{\theta}\approx F(s_{0})=\sum_{\tau\in\mathcal{T}}F(\tau)$.

Bengio et al. () have used GFlowNets to generate molecules with high binding affinity to a protein target by linking fragments to form a junction tree (Jin et al., ). For generative chemistry, the framework was extended to multi-objective optimisation (Jain et al., ; Roy et al., ), where the model was trained to simultaneously optimise for binding affinity to the protein target, Synthetic Accessibility (SA), drug likeness (QED) and molecular weight. Recently, Shen et al. () extended the framework to pocket-conditioned molecular generation, employing the same action space to generate molecules conditioned on different protein targets.

### 2.2 GFlowNets with a trained backward policy

GFlowNets train a forward policy $P_{F}$ to match the backward policy $P_{B}$ according to Eq. 1. The choice of $P_{B}$ therefore impacts the training of GFlowNets and sample quality. Despite this relationship, the choice of backward policy in GFlowNets has attracted limited attention but for a few works discussed here. Malkin et al. () proposed parameterizing the backward policy and training $P_{B}$ and $P_{F}$ simultaneously using the trajectory balance objective (Eq. 1). They also proposed fixing the backward policy to a uniform distribution when modeling the distribution over parents states proves difficult. Mohammadpour et al. () compare a maximum-entropy GFlowNet to GFlowNets with a uniform backward policy. Closer to our work, Jang et al. () propose Pessimistic GFlowNets, which use maximum likelihood over observed trajectories to train $P_{B}$. This ensures that the flow induced by the backward policy is concentrated around observed (training) states, making the model pessimistic about unobserved intermediate states having flow.
In this work, we address the idea of ensuring that backward-generated trajectories belong to the GFlowNet MDP. In doing so, we propose a data-driven solution to a MDP design challenge.

### 2.3 Generative models for molecule design

A large number of works have been proposed for generative molecular design (Du et al., ). Early methods employed techniques such as variational autoencoders (VAEs)
(Gómez-Bombarelli et al., ),
deep reinforcement learning (Segler et al., ; Olivecrona et al., ), and generative adversarial networks (GANs) (Guimaraes et al., ; Cao & Kipf, ).
More recently, methods like normalizing flows
(Zang & Wang, )
and diffusion models (Hoogeboom et al., ; Schneuing et al., ; Harris et al., ; Torge et al., ; Igashov et al., ) have gained significant attention due to their ability to model complex distributions more effectively.
Despite these advances, challenges remain, particularly in ensuring that generated molecules adhere to physical and chemical constraints (Stanley & Segler, ; Harris et al., ).

### 2.4 Synthesis-aware molecule generation

The idea of tackling molecule generation and synthesis simultaneously has been investigated early by Vinkers et al. (), who introduced SYNOPSIS, which generates molecules from a starting dataset of available compounds, relying on applying chemical modifications to functional groups and assessing the value of the product with a fitness function. The works of Bradshaw et al. () and Korovina et al. () provided early neural models for one-step synthetic pathways. Bradshaw et al. () generalized this idea as a generative model for synthesis DAGs, which can be optimized in a VAE or RL setup.
Gottipati et al. () used reinforcement learning to generate compounds from reactions and commercially available reactants. Gao et al. () formulate an MDP to model the generation of synthesis trees, which can be optimized with respect to the desired properties of a product molecule. Luo et al. () proposed a model that can project unsynthesizable molecules from existing generative models to synthesizable chemical space by utilizing postfix notations to represent synthesis pathways. Guo & Schwaller () showed that retrosynthesis models can be treated as an oracle in goal-directed molecule generation.
Concurrently to our work, Koziarski et al. () propose a similar GFlowNet-based framework for synthesizable molecular generation but using a different set of reaction templates and making use of a different backward policy.

While our work resembles an RL setting for synthesis-aware molecular generation (Gottipati et al., ; Horwood & Noutahi, ), the key difference lies in the sampling distribution of the learned model. Contrary to RL, the GFlowNet objective is not to generate the single highest-return sequence of actions, but rather to maximise both performance and diversity by sampling terminal states proportionally to their reward. This is especially useful in the context of molecule generation, where we want to explore different modes of the distribution of interest.

## 3 Methods

In this work, we present a framework to train GFlowNets on a Markov Decision Process (MDP) made of molecules obtained from sequences of chemical reactions. Below, we describe how this compositional space of synthesizable molecules is assembled (3.1), we present a method making use of backward policies to palliate imperfect information contained in reaction templates (3.2) and present the model used to navigate that state space and learn the target distribution (3.3).

### 3.1 MDP of chemically accessible space

#### Problem definition

We model synthetic pathways as trajectories in a GFlowNet, starting from purchasable compounds and ending with molecules that are optimized for some desired properties, via a set of permissible reaction templates. At each timestep $t$, the state $s_{t}$ represents the current molecule and stepping forward in the environment consists in building up the molecule by applying new pairs of reactions and reactants until either a termination action is chosen or the path reaches a maximum length. We encode reaction templates using RDKit reaction SMARTS (see Figure A.1).

#### Forward Action Space

We define five types of forward actions: Stop, AddFirstReactant, ReactUni, ReactBi, and AddReactant. ReactUni and ReactBi represent uni-molecular and bi-molecular reactions. The AddReactant action, which is available only after a ReactBi action, represents the choice of reactant for the bi-molecular reaction.
In more detail, each trajectory starts from an empty molecular graph which is followed by a building block sampled from AddFirstReactant. We then continue based on the sampled action type as follows: (a) if the action type is Stop, we reach a terminal state and end the trajectory; (b) if a ReactUni action is sampled, we apply the uni-molecular reaction template to the molecule in state $s$ to yield the product in state $s^{\prime}$; (c) if the action type is ReactBi, the sampled reaction is used as input to an additional MLP, together with the state embedding, to sample a subsequent action of type AddReactant.

#### Backward Action Space

The GFlowNet framework requires us to model travelling backward along trajectories in the state-space. To unfold a reverse trajectory we proceed similarly: (a) if the action-type is a BckReactUni, the action yields the reactant molecule directly; (b) if the action type is BckReactBi, we obtain two reactants, and the molecule that is not a building block becomes the next state (or previous state in the DAG). If the two resulting reactants are both building blocks (which happens at the beginning of the forward trajectory), the molecule that populates the next state is picked with $p=\nicefrac{{1}}{{2}}$ from the two building blocks.
The last action in a backward trajectory is BckRemoveFirstReactant, leading to the empty molecular graph $s_{0}$ (initial state).

#### Masking

Prior to sampling actions in both forward and backward direction, we ensure that the reactions and building blocks to be sampled are compatible with the current state using masks, obtained by checking for substructure match between the reaction template and the reactant (for forward actions) or product (for backward actions) molecules. We also enforce through masking that at least one of the resulting reactants when running a backward reaction is a building block.

#### Data

We use commercially available building blocks (BBs) from , which are small fragments of molecules prepared in bulk to be readily synthesised into candidate molecules.
Reaction templates are obtained from two publicly available template libraries (Button et al., ; Hartenfeller et al., ). After preprocessing, we obtain a total of 105 templates: 13 uni- and 92 bi-molecular reactions respectively (see App. A.2). We also use 12 Enamine REAL reactions in a small number of experiments (see Section 4.6 and App. A.2).

### 3.2 Challenges with Generating Backward Paths

Given the synthesis pointed DAG $\mathcal{G=(S,A)}$, one needs to define a forward probability function $P_{F}$ and a backward probability function $P_{B}$ both consistent with $\mathcal{G}$ (see Eq. 1). Contrary to previous fragment-based or atom-based molecule-generation GFlowNet environments, where any backward action can lead to $s_{0}$ (removing nodes and edges sequentially will lead to an empty graph; Bengio et al., ), defining $P_{B}$ in a reaction-based environment is non-trivial. This is because not every parent state (obtained by applying a reaction template backwards) will ensure that there exists a sequence of actions that leads all the way back to a building block, and therefore $s_{0}$. Note that the masking described in Section 3.1 is insufficient to account for this, as it does not ensure that the state obtained is *further* decomposable into building blocks. Consequently, to maintain a pointed DAG, no flow should be assigned to such a transition. A uniform backward policy, which is a standard choice in GFlowNet literature (Malkin et al., ), will fail at achieving this as it will assign positive flow to every backward action, including those leading to states that are not attainable from forward trajectories initialised in $s_{0}$ (see Figure 1-C). To address this issue, we explore a few training schemes for a parameterized $P_{B}$ that force backward-constructed trajectories to end in $s_{0}$.

#### Training the backward policy

We first explore a training scheme for $P_{B}$ which makes use of forward-generated trajectories. Similarly to Jang et al. (), we train $P_{B}$ using the maximum likelihood objective over trajectories generated from $P_{F}$:

$$ $\mathcal{L}_{B}(\theta_{B})=\mathbb{E}_{\tau\sim P_{F}}[-\log P_{B}(\tau; \theta_{B})].$ (2) $$

In that setting, we (1) generate trajectories using $P_{F}$, (2) update $P_{F}$ according to the trajectory balance objective in Equation 1 and (3) update $P_{B}$ using these same trajectories according to Equation 2 (see see Algorithm 1 in Appendix A.4).

While the maximum-likelihood approach presented above is sufficient to limit $P_{B}$ in allocating flow to paths that do not connect back to $s_{0}$, it restricts the model’s exploration by encouraging $P_{F}$ to collapse on a single path for each terminal molecule. To allow $P_{B}$ to ban erroneous paths while retaining a higher entropy, we also explore the use of policy gradient methods. Specifically, we explore maximising an expected *backwards* reward via REINFORCE (Williams, ), which is suitable for short trajectory environments like our reaction-based MDP:

$$ $J_{B}(\theta_{B})=\mathbb{}{E}_{\tau\sim P_{F},P_{B}}[R_{B}(\tau)]-\alpha H(P_ {B}).$ (3) $$

Here, $H(P_{B})=-\mathbb{E}_{\tau\sim P_{B}}[\log(P_{B}(\tau))]$ is the entropy term and the reward $R_{B}$ is set to 1 for a trajectory that ends in $s_{0}$ and -1 otherwise. In this setting we train the backward policy not only on the trajectories generated by the forward policy, but also on newly generated backward trajectories sampled directly from $P_{B}$ (see Algorithm 2, Appendix A.4).

Interestingly, training the backward policy to navigate back to $s_{0}$ is analogous to the retrosynthesis problem (Corey & Cheng, ; Segler et al., ). We employ this approach of training $P_{B}$ against a different objective than $P_{F}$ to palliate to a design challenge of the MDP, but a similar strategy could also be employed to fold additional preferences over different synthesis routes leading to the same terminal state into the system, for example to take into account the synthesis costs of a particular path. This possibility is left to future work.

### 3.3 Scaling of the reactant space to large number of building blocks

#### Estimating the size of the state space

Figure: Figure 2: Estimated size of state spaces. The full building block set contains 221,181 molecules. $L$ is the maximum trajectory length in the GFlowNet.
Refer to caption: extracted/5932713/figures/state_space_size.png

We leverage the properties of GFlowNets (Bengio et al., ) to estimate the size of the state space induced by our action space. Specifically, noting that GFlowNets learn $\log Z=\log\sum_{x}R(x)$, we train a model with $R=1$ for all terminal states to estimate their total count. We do so for different numbers of building blocks and different maximum trajectory lengths ($L$), and find that SynFlowNet using $L=3$ and 10k building blocks matches the size of the Enamine REAL space (), and that the size of the space quickly increases with the number of building blocks. We use our full set of 105 reactions. Note how reaction-constrained models considerably limit the exploration of the chemical space, with a fragments GFlowNet exploring a space $\sim$10 orders of magnitude larger.

#### Scaling method

To ensure synthetic accessibility of all samples, our model is inherently constrained by the initial set of available building blocks (BBs). To cover a large chemical space, it is thus crucial to use an extensive BB collection. The model must demonstrate scalability to accommodate larger sets of BBs, both in terms of training efficiency and overall performance. To do so, we change the representation of the BBs and their selection mechanism, as shown in Figure 1B.
Instead of the weight matrix of the mapping from hidden units to logits associated with BBs to be randomly initialized, it is *fixed* to be the matrix of binary Morgan fingerprints (Rogers & Hahn, ).

## 4 Results

Our results support a number of claims: (i) a reaction-based MDP greatly improves the synthesisability of generated molecules (Sec. 4.1), (ii) employing GFlowNets enables much more diverse molecule sampling over RL (Sec. 4.5), (iii) learning a backward policy results in higher-reward molecules and enables finding retrosynthetic pathways for molecules belonging to our state-space (Sec. 4.4), and (iv) our method of sampling molecules based on chemical fingerprints allows for efficient scaling to large chemical spaces (Sec. 4.5). Finally, we show that this approach also allows to improve results for particular programs by curating the set of building blocks based on target-specific experimental data (Sec 4.6). Unless otherwise specified, all experiments below use a backward policy trained with maximum likelihood, and Morgan fingerprint embeddings for our action space, with a library of 10,000 Enamine building blocks and 105 reaction templates.

### 4.1 Reaction-based MDP design

Figure: Figure 3: Example of synthesis pathway generated by SynFlowNet. The reward here is the sEH binding affinity proxy. Right shows the Vina docked pose and molecule metrics.
Refer to caption: extracted/5932713/figures/synflownet_synthesis_pymol.png

Figure: Figure 4: Comparison across MDPs. We evaluate four GFlowNet models: SynFlowNet is trained with an action space of chemical reactions and maximum trajectory lengths (L) of 3 and 4, FragGFN and FragGFN SA are trained with an action space of fragments, however the latter also optimises for synthetic accessibility (SA), besides the sEH binding proxy reward. SynFlowNet molecules are achieving higher binding scores and better synthesisability.
Refer to caption: extracted/5932713/figures/syn_vs_frag_vs_frag_sa_with_len3.png

We compare GFlowNets using two action spaces – fragments vs. reactions – and report reward, diversity and synthesizability of the generated samples.
We retrain the GFlowNet model proposed by Bengio et al. () with fragments derived from our building blocks set. This ensures that we use similar chemical spaces and that we get as close as possible to a fair comparison between the two MDPs. In the fragment space, it is common practice to optimise for synthetic accessibility scores to improve synthesizability from the model. We therefore train two versions of the FragGFN model: one using the sEH binding proxy as reward function, and another that optimizes for both sEH binding and synthetic accessibility (SA) score. SynFlowNet was only trained with the sEH proxy as reward.

The results in Fig. 4 show that our MDP design achieves great improvement in terms of synthesisability, while preserving high rewards and diversity.
We see that within the reaction environment, MDP design choices can further influence sample quality: a maximum trajectory length of 3 achieves better SA scores and AiZynthFinder (Genheden et al., ) retrosynthesis success (62%) compared to a maximum trajectory length of 4 (40%)(^1^11The reactions available to AiZynthFinder are the US patent office (USPTO) set (), which are different to the reactants used by SynFlowNet, so one should not expect a 100% success score.).
Note that synthesisability metrics are correlated with molecule size (Skoraczyński et al., ), and that molecules assembled from longer synthesis routes are naturally larger (see Table A.2). Both FragGFN models score 0% with AiZynthFinder. Diversity, measured as average Tanimoto distances between molecular fingerprints, is preserved from a fragments to a reaction environment, albeit less so for shorter trajectories. Figure 3 is an example of output from SynFlowNet, returning also the synthesis pathway extracted from the MDP.

### 4.2 GFlowNets as samplers of chemical space

Figure: Figure 5: Comparison between GFlowNet and RL. SynFlowNet (GFlowNet with synthesis actions) discovers more modes compared to entropy-regularised RL trained on the same state & action space.
Refer to caption: extracted/5932713/figures/comparison_to_baselines.png

While the works of Gottipati et al. () and Horwood & Noutahi () are closest to our method, where RL is paired with a synthesis action space, their code is not publicly available. We reproduce their setting by pairing an RL algorithm with our MDP. We train our model with a soft Q-learning (Haarnoja et al., ) objective and compare to SynFlowNet on three different reward functions: sEH, GSK3$\beta$ and DRD2 (see Appendix A.1).

In Fig. 5, we show the reward distribution of generated samples after equal numbers of training steps. We notice that soft Q-learning collapses to sampling a narrower distribution of high-reward molecules for all targets. This is also shown in the number of Bemis-Murcko scaffolds counted for molecules with rewards above a certain threshold, where soft Q-learning quickly collapses (note the logarithmic scale of the y-axis).

### 4.3 Further comparison to baselines

We continue by comparing SynFlowNet to strong baselines from the literature. First, REINVENT uses a policy-gradient method to tune an RNN pre-trained to generate SMILES strings (Blaschke et al., ; Loeffler et al., ). REINVENT has been shown to outperform many models in terms of sample efficiency (Gao et al., ) and proposing realistic 3D molecules upon docking (Ciepliński et al., ). Second, SyntheMol shares a similar MDP design with SynFlowNet, but is a search-based method using Monte Carlo Tree Search (Swanson et al., ). Contrary to our model, synthesis pathways in SyntheMol contain one reaction, and its estimated state space has 30 billion molecules, reached from the Enamine building block set and 13 reactions (Swanson et al., ).

The comparison is summarized in Figure 6. SynFlowNet achieves comparable rewards to REINVENT, and better sample diversity.
Although REINVENT does not use explicit knowledge of synthesizability in its generation process, it implicitly optimizes for it by being a likelihood model overfit on a curated ChEMBL dataset (Zdrazil et al., ). This is reflected in our last sub-plot, where we look at the maximum similarity to ChEMBL molecules, as a proxy for novelty. This proxy is strong enough, as ChEMBL contains one of the largest collection of diverse drug-like molecules to date (Zdrazil et al., ).
While REINVENT seems to be competitive in terms of synthesizability and reward, it generates molecules with poor novelty, staying close to its pretraining distribution. Notably, SynFlowNet achieves high novelty for high-scoring molecules. Overall, SyntheMol scores comparably to SynFlowNet, however SynFlowNet achieves a better reward/diversity trade-off and SynFlowNet’s capabilities are extended to explore larger state spaces (see Fig. 2). We repeat the comparison on other targets, and again find that SynFlowNet outperforms SyntheMol (see App. B).

Figure: Figure 6: SynFlowNet is competitive against other popular models from the literature. SynFlowNet achieves a good balance between reward optimisation, diversity, synthesizability and novelty (assessed by maximum similarity to ChEMBL molecules). REINVENT stays close to its pretraining distribution, harming the novelty of the proposed molecules. SynFlowNet is closest to ideal.
Refer to caption: extracted/5932713/figures/synflownet_vs_reinvent_vs_synthemol.png

### 4.4 Improved MDP consistency through trained backward policy

We evaluate the effectiveness of our proposed training schemes for the backward policy $P_{B}$ by measuring (1) whether $P_{B}$ can ensure that backward-constructed trajectories reliably find a trajectory back to the initial state $s_{0}$ and (2) whether it brings any benefit to the forward policy. In Table 1 we compare a fixed, uniform backward policy to three versions of a parameterized backward policy: a free policy, updated w.r.t. the trajectory balance loss, (see Eq. 1), a policy trained with maximum likelihood on the forward-generated trajectories and a policy that is allowed to explore the backward action space and trained with REINFORCE to find paths leading back to $s_{0}$. We observe that both the maximum likelihood and REINFORCE policies succeed in ensuring that backward flow is not lost outside of the MDP, as they manage to construct trajectories that start from terminal states sampled from $P_{F}$ all the way to $s_{0}$. We refer to such routes as solved routes (train) in Table 1, as the terminal states have been visited by the GFlowNet during training. We also test the ability of the trained policies to retrieve synthesis routes for molecules which have not been visited during training ($test$). These were obtained from a random sampler using our set of reactions and building blocks. While we see small differences in the number of high-reward modes discovered for different backward policies, we see that a free policy consistently fails to be competitive with the rest. Overall, the maximum likelihood and REINFORCE policies prove effective in guiding $P_{F}$ to high reward modes, while ensuring MDP consistency.

**Table 1: Effect of different training paradigms for $\boldsymbol{P_{B}}$. Training the GFN ${P_{B}}$ ensures that backward-constructed trajectories belong to our MDP and can marginally improve ${P_{F}}$. % of solved routes refers to the ability of $P_{B}$ to retrieve synthesis routes for on-policy (train) and off-policy (test) molecules reachable through our state space. Test molecules have not been seen during training. # of high reward ($R>0.9$) modes from $P_{F}$ is reported out of $\sim 500,000$ samples seen during training.**
| $P_{B}$ policy | % of solved routes (train) | % of solved routes (test) | # of high reward modes from $P_{F}$ |
| --- | --- | --- | --- |
| Uniform | $11.0\pm 3.7\%$ | $11.0\pm 4.1\%$ | $47,515.0\pm 11,264.7$ |
| Free | $67.3\pm 3.7\%$ | $1.0\pm 0.8\%$ | $6,754.3\pm 4,980.1$ |
| MaxLikelihood | $99.3\pm 0.5\%$ | $32.3\pm 7.3\%$ | $37,708.2\pm 13,992.6$ |
| REINFORCE | $100.0\pm 0.0\%$ | $44.3\pm 2.6\%$ | $55,387.6\pm 28,886.3$ |

### 4.5 Scaling to the entire Enamine set

Figure: Figure 7: Ablation over action-selection mechanism. Average reward of samples (A), utilization of building blocks (B), and training times (C) for both models on different subsets of building blocks. Embeddings enable more efficient learning (D), as they help the agent navigate the set of available building blocks and tend to exploit more relevant BB clusters (F) rather than explore new clusters (E). For each experiment, three models trained with different random seeds were used.
Refer to caption: x2.png

In this section, we study the ability of SynFlowNet to scale to large sets of building blocks using the method proposed in Section 3.3 as opposed to using a learnable embedding for each building block.
First, we compare the performance of our models using the softmax action-selection mechanism and the fingerprints-based softmax on the full Enamine set of $221,181$ building blocks and randomly sampled subsets containing $75\%$, $50\%$, $25\%$, $10\%$, and $3\%$ of the full set. As shown in Figure 7A, increasing the number of building blocks in general negatively affects the quality of the sampled molecules in terms of reward. However, the model with fingerprints consistently outperforms the former and makes the reward degradation effect almost negligible. Next, we study whether models are able to use the blocks they are given. As shown in Figure 7B, the model with fingerprints consistently uses more unique building blocks than the baseline model. Notably, the baseline model is only able to use a small subset of building blocks to produce high-quality samples ($R>0.8$). In contrast, fingerprints enable the use of a 5-fold larger set of blocks for the same quality level, thus significantly surpassing the baseline samples in terms of diversity. Finally, we measure the average time of the forward pass and average time of forward pass with auxiliary RDKit operations like template matching. Figure 7C shows that increasing the number of building blocks slows down training, and using fingerprints introduces additional overhead.

As shown in Figure 7D, fingerprints enable more efficient learning on a large building block set. Indeed, instead of memorising all available BBs, the model learns to navigate in the space of embeddings and to maximise dot product with the relevant building blocks. To illustrate this, we clustered all available building blocks based on Tanimoto similarity. Next, we studied how both models explore new building blocks within and across these clusters. As shown in Figure 7E-F, the former model tends to explore more clusters and less unique building blocks within clusters. The model with fingerprints, on the contrary, is able to focus on clusters that maximise the reward.

### 4.6 Guiding SynFlowNet with experimental fragment screens

Figure: Figure 8: Effect of curating the building block set. (A) We adapt SynFlowNet’s building block library based on experimentally validated fragments for a given target. (B) The curated set improves Mpro rewards over random building blocks.
Refer to caption: extracted/5932713/figures/mpro_main_text.png

One advantage of our framework is that building block sets can be used to specialize the model for a particular target.
Fragment-based Drug Discovery (FBDD) (Thomas et al., ) is a major strategy to increase the efficiency of drug discovery campaigns. In FBDD, instead of screening large chemical spaces, a low number of small compounds (called fragments) are screened and experimentally validated to bind, and are then linked in silico to make full molecules.

We hypothesise that experimental data from x-ray fragment screens (Murray & Rees, ) can guide and enhance SynFlowNet’s capabilities. As a proof-of-concept, we focus on the SARS-CoV-2 main protease (Mpro), leveraging strucural data from the COVID Moonshot project (Boby et al., ). We compile a building block set from Enamine with high similarity to fragments confirmed by x-ray crystallography (see Fig. 8A). We find that biasing the building block library towards fragments with known protein-ligand complementarity increases the reward over randomly selected molecules (see Fig. 8B). Further details on methods and results are given in App. Section B.3.

## 5 Conclusion

In this work, we discuss the application of GFlowNets to de novo molecular design paired with forward synthesis. We demonstrate that an action space of chemical reactions is an effective way of enforcing synthesisability, and that pairing it with GFlowNets excels in terms of diversity. Comparisons to state-of-the-art baselines emphasised that SynFlowNet explores novel regions of the chemical space. We also proposed a novel paradigm for training the backward policy in the GFlowNet and in doing so we improved and validated the correctness of our MDP design. Furthermore, we studied the building block exploration and exploitation mechanisms of SynFlowNet, showing efficient scaling to using hundreds of thousands of building blocks. Finally, as a proof of concept for the adaptability of SynFlowNet to real drug discovery programs, we showed that the framework can be specialised to target-specific molecule generation by making use of experimental fragment screens.

## Acknowledgements

We thank Enamine for providing the building blocks and for useful conversations. We thank Daniel Cutting, Austin Tripp and Jeff Guo for discussions regarding baselines. Miruna Cretu is supported by the SynTech CDT at the University of Cambridge and part of the work was completed during an internship at Valence Labs.

## Appendix A Extended methods

### A.1 Reward functions

We train SynFlowNet for a number of reward functions/targets. In some cases (e.g. when benchmarking against the Fragment-based GFlowNet), we also perform multiple objective optimisation with SA score.

#### sEH

Our main reward function is defined as the normalized negative binding energy as predicted by a pretrained proxy model, available from Bengio et al. () and trained on molecules docked with AutoDockVina (Trott & Olson, ) for the sEH (soluble epoxide hydrolase) protein target, a well studied protein which plays part in respiratory and heart disease (Imig & Hammock, ). The proxy model, which utilizes the weights from Bengio et al. (), was trained using a message-passing neural network (MPNN) (Gilmer et al., ) that processes atom graphs as input. Details of the model architecture are available in Bengio et al. (). It was trained on a dataset of 300,000 randomly generated molecules, achieving a test mean squared error (MSE) of 0.6. Note that the reward scale in our results differs from the original GFlowNet publication, with rewards adjusted by a factor of 1/8 in our analysis.

#### gsk3 β 𝛽 \beta italic_β and DRD2

We also employ two oracle functions from the PMO Gao et al. () benchmark, which provide machine learning proxies trained fit to experimental data to predict the bioactivities against their corresponding disease targets. The two targets we use here are gsk3$\beta$ Li et al. () and dopamine receptor D2 (DRD2) (Olivecrona et al., ).

#### Easy adoption to other targets using GPU-accelerated Vina docking

Finally, we wish for users to rapidly be able to adapt SynFlowNet to learn binding for their target of interest without having to retrain a new proxy or relying on slow docking simulations. We accomplish this using the new GPU-accelerated Vina-GPU 2.1 docking algorithm (Tang et al., ; Alhossary et al., ). For our experiment, we use the PDB:6W63 (Mesecar et al., ) structure for the SARS-CoV-2 main protease (Mpro) target and PDB:2XJX (Murray et al., ) for Heat Shock Protein 90 (HSP90) target. Receptors are prepared for docking with prepare_receptor4.py and the center of the docking is defined as the center of mass for the ligand with a size of 25
$\mathrm{\SIUnitSymbolAngstrom}$ in accordance with previous work (Buttenschoen et al., ).

### A.2 Chemical reactions

Figure: Figure A.1: Example of reaction templates (Hartenfeller et al., ). The templates act as rules, which match any molecule that left side (before the arrow) as a subgraph. The matching part is then transformed into the right hand side of the rule. “R” represents any group, “A” represents an aromatic atom. Note the implicit reagent used in reaction #2.
Refer to caption: extracted/5932713/figures/templates.png

#### Reaction pre-processing and action masking

Reaction templates pre-processing was necessary to ensure that the templates could be run backwards.
All templates containing wildcards (*) in the reactant SMARTS were duplicated with replacements for all atoms that the wildcard substituted for.

When running experiments with a REINFORCE backward policy (which requires to sample trajectories backward on-policy), we enforced additional masking that forward reactions are sampled if and only if they are reversible (i.e. when applying the template backward on the product, the same reactants are obtained as the ones utilised by the template in the forward direction).

Note that reaction templates can be uni-molecular (employing one reactant, Fig. A.1 #2) or bi-molecular (employing two reactants, Fig. A.1 #1).

#### Enamine REAL Space Reactions

For the case study in Section 4.6, we used reactions from to stay close to an experimental setup where the molecules generated would be readily purchasable from Enamine.
Enamine assembled a molecular space called the REAL space, which is a vast catalogue of 48B purchasable compounds(^2^22[https://enamine.net](https://enamine.net)). We use a subset of Enamine reactions available from Swanson et al. () which produces 93.9% of the REAL space(^3^33Note: this dataset contains a single tri-molecular reaction which has been removed here for simplicity).

### A.3 Model and training

A graph neural network based on a graph transformer architecture (Yun et al., ) is used to parameterize the forward and backward policies. The model’s action space is defined using separate MLPs for each action type (see Figure 1). Our model is trained in an online fashion, meaning that it learns exclusively from trajectories sampled from the GFlowNet policy, without relying on an external dataset of trajectories or a set of target molecules. Note however that this framework is compatible with offline training, which makes use of such datasets as starting point for exploring the molecular space.

#### SynFlowNet Training

We adapt the framework from Bengio et al. () to train a GFlowNet sampler over a space of synthesisable molecules, which are assembled from an action space of chemical reactions and reactants. A graph neural network with a graph transformer architecture (Yun et al., ) is used to produce a state-conditional distribution over the actions. A state is represented as a molecular graph in which nodes contain atom features. Edge attributes are bond type and the indices of the atoms which are its attachment points. This representation is augmented with a fully-connected virtual node, which is an embedding of the conditional encoding of the desired sampling temperature, obtained using an MLP. The sampling temperature is controlled by a temperature parameter $\beta$, which also plays a role in reward modulation, allowing for exponential scaling of the rewards (by making rewards received during training equal to $R^{\beta}$). We experimented with sampling $\beta$ from multiple distributions, and use a constant distribution in the reported results in this paper. We used a thermometer encoding of the temperature (Buckman et al., ).

The model is trained using the trajectory balance objective (Malkin et al., ) and thus is parameterized by forward and backward action distributions $P_{F}$ and $P_{B}$ and an estimation of the partition function $Z=\sum_{\tau\in\mathcal{T}}F(\tau)$.

For the state space estimation experiment, we used $R=1$, $10\,000$ training steps and varying trajectory lengths.
For the rest of the SynFlowNet experiments with SeH proxy as reward, we used $\beta=32$, 5000 training steps and varying trajectory lengths (see main text).
For training with the DRD2 and gsk3$\beta$ targets as reward, we use a maximum trajectory length of 4 and a reward exponent $\beta=16$.
For the backward policy training, see App. A.4. The rest of the hyperparameters are fixed and are presented in Table A.1.

#### FragGFN training

We obtained fragments and their attachment points from Enamine building blocks by following the protocol provided by Jin et al. (Jin et al., ). The model trained with the sEH reward was trained with the default implementation in Yoshua Bengio (), and an optimised reward exponent $\beta=64$. The model optimising for synthetic accessibility, as well as sEH binding, was trained with a reward obtained by multiplying the two scores.

#### Soft Q-learning training

We implemented a version of Soft Q-learning (SQL) (Haarnoja et al., ), an energy-based policy learning method, that operates on SynFlowNet’s MDP. We optimised the training procedure by performing both manual and grid searches across several values of entropy regularisation parameter $\alpha$ and reward scaling parameter $\beta$. For $\beta$, contrarily to GFlowNets, since the method does not only learn a policy but tries to estimate the Q-values of each actions directly, we found that using large values of reward scaling such as $\beta=64$, which are common for GFlowNets, would destabilise the algorithm and had to be lowered to $\beta=4$ or $\beta=2$. For $\alpha$, we tried several values that would strike the best tradeoff between allowing the model to find high-performing molecules while maximising diversity and avoiding to collapse the agent’s distribution on only a few modes, with our best model using $\alpha=0.01$.

#### REINVENT training

We benchmarked our approach against REINVENT4 (Loeffler et al., )(^4^44We used the code available at [https://github.com/MolecularAI/REINVENT4](https://github.com/MolecularAI/REINVENT4)). Following the setup described in their methodology, we fine-tuned the REINVENT prior model using reinforcement learning. We change little else other than the reward function used to train the model. The training was conducted with the default batch size of 100 over 3,000 training steps, resulting in a total of 300,000 oracle calls, which is consistent with the number of oracle calls used in SynFlowNet experiments.

#### SyntheMol generation

We run SyntheMol’s (Swanson et al., ) Monte Carlo Tree Search (MCTS) algorithm using their standard sets of $139\,493$ building blocks and 13 chemical reactions from the Enamine REAL space, but replace the provided bioactivity prediction models with the reward functions used in this study (Sec. A.1).
We use the publicly available code repository(^5^55[https://github.com/swansonk14/SyntheMol](https://github.com/swansonk14/SyntheMol)) and perform $50\,000$ rollouts.
Other hyperparameters were kept at default settings, including a maximum of 1 reaction and exploration parameter $c=10$.
Before the search, building block scores were pre-computed using the target reward functions.
For the final selection, we sample 1000 molecules randomly from all returned molecules.

**Table A.1: Hyperparameters used in our SynFlowNet training pipelines.**
| Hyperparameters | Values |
| --- | --- |
| Batch size | 64 |
| Number of GNN layers | 4 |
| GNN node embedding size | 128 |
| Graph transformer heads | 2 |
| Learning rate ($P_{F}$) | $10^{-4}$ |
| Learning rate ($P_{B}$) | $10^{-4}$ |
| Learning rate ($Z$) | $10^{-3}$ |
| max_nodes | 9 |

### A.4 Backward Policy Training Algorithms

For the maximum likelihood backward policy, at each training iteration, $P_{B}$ is updated to minimize the maximum likelihood loss over trajectories generated in that batch from $P_{F}$. Similarly for a REINFORCE $P_{B}$, we train $P_{F}$ as above, but maintain a replay buffer. To train $P_{B}$, we sample terminal states from the buffer and sample trajectories backwards from $P_{B}$, which are used in a REINFORCE update. To improve training stability, we also used trajectories generated from $P_{F}$ (in a 1:1 ratio) to update $P_{B}$.

Contrary to the rest of the models (see App. A.3), the backward policy models were trained for $8000$ steps, with reward exponent $\beta=64$ and $\text{max_len}=5$. The results in Section 4.4 are reported for a REINFORCE loss with an entropy multiplier term $\alpha=1.0$.

Figure: Algorithm 1 Training of Maximum Likelihood Backward Policy for GFlowNets

Figure: Algorithm 2 Training of REINFORCE Backward Policy for GFlowNets

### A.5 SynFlowNet scaling

For Section 4.5, to compute clusters of building blocks, we used BitBIRCH (Jung et al., ), a recent adaptation of Balanced Iterative Reducing and Clustering using Hierarchies (BIRCH) algorithm (Zhang et al., ), recently proposed for efficient clustering of large molecular libraries.

## Appendix B Extended Results

### B.1 Comparison to baselines

Table A.2 contains metrics for all the models and baselines we have run. Figure A.2 contains results for molecules generated with SyntheMol for the DRD2 target, compared to SynFlowNet. SynFlowNet is superior in finding high-reward molecules. We tried various hyperparameters for SyntheMol to optimise for the gsk3$\beta$ target and were unsuccessful – no nonzero rewards were found.

**Table A.2: Comparison to baselines. Results obtained by averaging over 1000 random molecules sampled from the trained models. Standard errors obtained from training using 3 seeds. Due to high computational cost, AiZynthFinder scores are computed over 100 random samples.**
| Method | sEH proxy ($\uparrow$) | Diversity ($\uparrow$) | SA ($\downarrow$) | AiZynth ($\uparrow$) | QED ($\uparrow$) | Mol. weight ($\downarrow$) | ChEMBL Similarity ($\downarrow$) |
| --- | --- | --- | --- | --- | --- | --- | --- |
| REINVENT | $0.91\pm 0.01$ | $0.68\pm 0.02$ | $2.19\pm 0.02$ | $0.95$ | $0.57\pm 0.04$ | $429.78\pm 23.02$ | $0.64\pm 0.02$ |
| FragGFN | $0.77\pm 0.01$ | $0.83\pm 0.01$ | $6.28\pm 0.02$ | $0.00$ | $0.30\pm 0.01$ | $724.62\pm 32.39$ | $0.24\pm 0.09$ |
| FragGFN SA | $0.70\pm 0.01$ | $0.83\pm 0.01$ | $5.45\pm 0.05$ | $0.00$ | $0.29\pm 0.01$ | $683.31\pm 59.92$ | $0.21\pm 0.01$ |
| SyntheMol | $0.64\pm 0.01$ | $0.86\pm 0.01$ | $3.08\pm 0.01$ | $0.82$ | $0.63\pm 0.01$ | $412.24\pm 0.98$ | $0.49\pm 0.01$ |
| Soft Q-Learning | $0.80\pm 0.07$ | $0.42\pm 0.04$ | $2.63\pm 0.39$ | $0.96$ | $0.39\pm 0.02$ | $408.02\pm 12.76$ | $0.52\pm 0.02$ |
| SynFlowNet ($L=3$) | $0.92\pm 0.01$ | $0.79\pm 0.01$ | $2.92\pm 0.10$ | $0.65$ | $0.59\pm 0.02$ | $365.23\pm 2.42$ | $0.43\pm 0.01$ |
| SynFlowNet SA ($L=3$) | $0.94\pm 0.01$ | $0.75\pm 0.02$ | $2.67\pm 0.03$ | $0.93$ | $0.68\pm 0.01$ | $358.27\pm 3.52$ | $0.48\pm 0.01$ |
| SynFlowNet ($L=4$) | $0.88\pm 0.01$ | $0.82\pm 0.01$ | $3.54\pm 0.03$ | $0.40$ | $0.27\pm 0.01$ | $557.49\pm 8.60$ | $0.38\pm 0.01$ |
| SynFlowNet QED | $0.86\pm 0.03$ | $0.81\pm 0.03$ | $4.02\pm 0.26$ | 0.55 | $0.74\pm 0.04$ | $398.50\pm 8.84$ | $0.38\pm 0.01$ |

Figure: Figure A.2: Comparison between SynFlowNet and SyntheMol with DRD2 as reward. SyntheMol struggles to find high-reward molecules. We also run SyntheMol with GSK3$\beta$ and were not able to optimise it. Results for SynFlowNet with both DRD2 and GSK3$\beta$ are reported in Figure 5.
Refer to caption: extracted/5932713/figures/drd2_synflownet_vs_synthemol.png

Figure: Figure A.3: (LHS) Effects of different parameterizations of the backward policy. For REINFORCE, ($\alpha$) represents the multiplier to the entropy term. (RHS) REINFORCE is trained on-policy and the mean rewards of the sampled backward trajectories are plotted. Note that the backward reward is 1 if the backward-constructed trajectory ends in $s_{0}$ and -1 otherwise. A baseline for Uniform policy is shown, obtained by sampling 100 random trajectories backward from terminal objects in the DAG. Results averaged over 8 seeds.
Refer to caption: extracted/5932713/figures/bck_policies_modes.png

### B.2 Example Trajectories

Figure A.4 shows examples of molecules and synthesis pathways generated from SynFlowNet.

Figure: Figure A.4: SynFlowNet-generated molecules and synthesis pathways for SeH.
Refer to caption: extracted/5932713/figures/synthesis_pathways_app.png

### B.3 Curated building blocks sets given target data

Here we provide further details on the methods and rationale for the fragment-based building block curation strategy outlined in 4.6. To enable high molecule purchaseability in a real world drug discovery campaign, we use the Enamine reactions in all these experiments. GPU-accelerated Vina docking was used for reward computation (see App. A.1).

#### Experimental fragment extration pipeline

Fragment-based drug discovery (FBDD) simplifies the search space by focusing on smaller molecular fragments rather than designing or screening full molecules, enabling a more efficient exploration of potential binding interactions (Murray & Rees, ; Thomas et al., ). To identify molecules bound within the same pocket across different protein structures, we performed a sequence-based search using MMSeqs2 (Steinegger & Söding, ) across the Protein Data Bank (PDB). We retained only hits with a sequence identity of 90% or greater and an alignment overlap exceeding 80% to ensure high similarity to the reference target. The resulting PDB structures were structurally aligned(^6^66We use the superimpose_homologs function from biotite: [www.biotite-python.org](https://www.biotite-python.org/latest/index.html)) according to the ligand-bound chain with the reference structure. Any ligand with at least one atom within 2 Å of the reference ligand in the binding site was selected for further analysis.

Figure: Figure A.5: Curation of HSP90 fragments and building blocks. Red line shows threshold at which we denote a molecule a ‘fragment’ and was choosen to max the upper limit of the typical BB size.
Refer to caption: x3.png

#### Selection of small molecule fragments

To isolate small molecular fragments suitable for FBDD, we filtered out ligands containing more than 25 atoms, as this atom count aligns with typical building block molecule sizes. This is done for 2 reasons; (i) we can pretend that we are in a drug discovery campaign where a fragment screen has just been conducted and there is no leaked information from ‘full’ molecules and (ii) this threshold aligns well with the typical Enamine building block in terms of size (Figure A.5).

#### Curation of Enamine Building Blocks

For every experimentally validated fragment, we selected the top 100 closest building block molecules from the Enamine library based on molecular similarity. After removing duplicate SMILES entries, this process resulted in thousands of curated building blocks for each target.

#### Example fragment screens

We study two targets. The first is Heat Shock Protein 90 (HSP90), for which there are a large number of ligand bound structures, mostly thanks to FBDD campaigns conducted in industry (Murray et al., ; Woodhead et al., ). The second is the SARS-CoV-2 Main Protease (Mpro), for which there is a large amount of structural data and in particular fragments from a fragment screen performed by the COVID Moonshot open-science initiative (Boby et al., ).
In the case of HSP90, using PDB entry 2XJX (Woodhead et al., ) as the reference structure, our ligand extraction pipeline identified 92 molecules. Filtering for small molecular fragments of fewer than 25 atoms further refined this set for analysis. The methodology was repeated for Mpro using PDB entry 7GAW, where we extracted 138 molecules.

Examples of molecules designed for Mpro using our curated building block set are shown in Figure A.6.

#### Alignment of SynFlowNet designs with experiments

Murray et al. () performed a fragment screen of HSP90 and identified two distinct lead classes of binding modes: (i) an aminopyrimidine class that formed hydrogen bonds with an asparagine residue and several conserved water molecules within the pocket and (ii) a phenolic class that primarily binds through water-mediated hydrogen bonds networks. These complicated water-dependent interaction networks are extremely challenging to model computationally and is a clear limitation of Vina docking (which does not take into account waters). Figure A.7 shows that by guiding SynFlowNet in this way, we can design molecules and poses that fit with these experimentally validated interactions.

Figure: Figure A.6: Example molecules generated by SynFlowNet for the Mpro case study. Green molecules are those generated by SynFlowNet. Magenta shows the molecule from the reference structure. ‘Sim’ is the Tanimoto similarity between a designed molecule and the reference structure.
Refer to caption: extracted/5932713/figures/mpro_pymol.png

Figure: Figure A.7: SynFlowNet designs molecules that align with real world fragment experiments. Murray et al. () performed a fragment screen for the HSP90 target and identified two classes of binding modes for the target: those built around (A) phenolic compounds (e.g. PDB:5M4H) and (B) aminopyrimidine compounds (e.g. PDB:2YEE). SynFlowNet trained on a target-curated building block set was able to consistently design compounds based on these experimentally validated binding modes, meaning we can be more certain that they are likely to bind well in reality. Green shows SynFlowNet-designed molecule and magenta are fragments from x-ray crystallography experiments.
Refer to caption: extracted/5932713/figures/hsp_pymol.png

## Peer Reviews

confidence: 4
contribution: 3
initial_rating: 8
presentation: 4
questions: * I could not find how $s_0$ is initialized, nor how important this is (i.e., how different initializations here impact the performance of SynFlowNet). Can the authors clarify? Does it correspond to a molecular fingerprint/embedding of the target molecule, is it random noise, or something else?
* One of the contributions of this paper is stated as "we propose one of the first attempts at training the backward policy in a GFlowNet with a separate objective from for the forward policy, which we show to correct backward-generated trajectories and improve sample quality and diversity" but I am not sure I got that. Which experiment was it exactly that demonstrated the new backward policy led to improved sample quality/diversity? Are these results in the paper? If so this needs to be clarified.
* How sample efficient is SynFlowNet? I would have liked to see an analysis of this, which is really important for molecular generation tasks as we are almost always operating in the low data regime and data/sample efficiency is key.
* How does SynFlowNet perform when optimizing for targets that are highly symmetric, e.g., molecules with a Cn rotation axis, where n$\ge$2? Has this been evaluated? This is something that is a known weakness of other synthesizability-constrained molecular generation frameworks (e.g., SynNet). 
* What are the failure modes of SynFlowNet? This should be more thoughtfully analyzed, as all models have failure modes and it is important to understand them to help future research build on the results of SynFlowNet.
* How challenging would it be to expand SynFlowNet beyond bi-molecular reactions? I am not suggesting this should have been done here, but asking more from curiousity and wondering if this is a future direction the authors are aiming to explore in future work.

Minor:
* Are the points in Figure 6 (last panel) shuffled for plotting? It was not clear to me if there are red points under the blue ones or vice-versa. This should be clarified as the results from this plot are informing the conclusions.
* The references are poorly formatted, this should be fixed.
* gsk3$\beta$ should be GSK3$\beta$ (weird capitalization in the Appendix).
soundness: 4
strengths: * The background and related work section was strong, and covered the previous work in this domain well.
* The paper, including the appendix, was well-structured and easy to follow.
* There are more strengths, but I mention them below amongst some of the specific questions I had, to put the questions in context.
summary: It was a pleasure to review this paper on synthesizability-constrained molecular generation using GFlowNets. Overall, I felt the paper was strong and focuses on an important application area of generative models, it was well-written, and generally builds well on previous work in the field and introduces some novelty. In my review below, I focus on aspects and unresolved questions that, if addressed, could enhance the paper. I summarize my main concerns as: (1) in its current state, the work is not reproducible, and (2) the evaluation metrics and comparisons to prior work used in evaluating SynFlowNet are not always the most relevant for molecular design, and do not tell us much about how this method would actually perform in a realistic setting. There is room for improvement here. I think this paper has the potential to be a really strong paper if it addresses the concerns below.
weaknesses: * I searched throughout the paper and in the submission, and could not find any code to accompany the paper. Given that the authors lamented the lack of publicly-available code for two previous works in their paper, I would have expected a greater emphasis on publishing the code here, as it would also help with the reproducibility of this work. To me this is unacceptable for a computer science conference submission in 2024, so please let me know if I simply missed it, and I can raise my score (I ranked it below the acceptance threshold because of this, and would update my score if the code was made available as it was overall a good paper).
* Throughout the paper, the authors repeat that one of the key advantages of SynFlowNet is its ability to generate more diverse candidates than other methods (namely, REINVENT and SyntheMol). However, I argue that diversity is only a meaningful metric in molecular generative models if the types of molecules being generated are from the accessible molecular space; if the generated molecules are diverse but not accessible, then the model is not very practically useful. That being said, it was good the the authors evaluated synthesizability using AiZynthFinder, but an "easy" experiment that the authors could have also done that would add a lot of value is to consider molecular rediscovery benchmarks. That is, can the model also find known, previously-synthesized binders; celecoxib rediscovery is a widely-used standard benchmark for this, but there are also more challenging targets. Such rediscovery metrics are a nice alternative to experimental validation which is not necessarily possible on the review timeline. To summarize, I would like to know that SynFlowNet can generate diverse and *known* synthesizable molecules, not just diverse *predicted* synthesizable molecules, because at the end of the day CASP tools are not a replacement for experiment.
* I found the lack of comparison to other methods for synthesizability-constrained molecular generation a weakness. The authors mentioned that the models of Gottipati and Horwood & Noutahi were not publicly available, and thus focused on comparison to SyntheMol and REINVENT. While these comparisons were interesting and certainly well-motivated, I would have also liked to know how SynFlowNet compares to other synthesizability-constrained molecular generation models which *do* have publicly available codebases, like MoleculeChef or SynNet. As it stands, we do not know which model is best, because SyntheMol has not done any comparisons to those models either.
* In the Appendix, the authors state that 300K oracle calls were used both to train REINVENT and SynFlowNet, however, this is an extreme number of oracle calls; can the authors justify this? REINVENT demonstrates good performance after only 10K oracle calls for most objectives, so are 300K oracle calls simply what is needed for good performance in SynFlowNet? I would like to see the analogous experiments and comparisons/results with only 10K oracle calls, which is the standard number of oracle calls used for comparing molecular generative models.

---

confidence: 4
contribution: 4
initial_rating: 6
presentation: 3
questions: 1. In Table 1, while the two backward policies, the Maxlikelihood and REINFORCE, show relatively high ratios of the solved routes in both the training and test sets, their number of high reward modes is similar to that of the Uniform Policy. This result is rather disappointing, considering the objective of the generative model is to sample more high-reward modes. In this context, the authors need to discuss the effect of the proposed backward policies on the performance of the forward policy in more detail.
2. In Figure 4, lower rewards achieved by FragGFN and FragGFN SA seem questionable because FragGFN covers a much larger chemical space, as seen from Figure 2, and so is more likely to find molecules with higher rewards. Has the poor performance of FragGFN originated from insufficient training? Or is it because of the use of SeH proxy as a reward, which may have a severe generalization problem for unseen molecules? This should be justified in the revision. 
3. In Section 4.1, what is the rewarding strategy used for multi-objective optimization (FragGFN SA, SynFlowNet SA, SynFlowNet QED)? The possible reward structures could be Multi-objective GFlowNet (MOGFN) or multiplication-based rewards.
4. Regarding the metric Reward/Novelty trade-off, it is possible that improving the novelty for SynFlowNet may be simply due to the use of the Enamine building blocks, while REINVENT has been trained with ChEMBL. Thus, the high similarity of the molecules generated by REINVENT to ChEMBL is not unexpected and would be reduced if training it with the Enamine REAL. The authors need to address this issue.
soundness: 3
strengths: 1. The use of the GFlowNet improves both diversity and high reward modes compared to other optimization methods such as reinforcement learning.
2. The proposed backward policy ensured that backward-constructed trajectories belong to the GFlowNet MDP. 
3. Various benchmark and ablation studies were performed to ensure the strategies of the proposed framework.
summary: The authors proposed a molecular generative model that explicitly considers the synthesizability of generated molecules. Specifically, they formulated a molecular design problem into a Markov decision process (MDP) that successively adds a purchasable building block (BB) with a selected reaction to the previous state. A GFlowNet has been adopted to tackle the MDP. Namely, SynFlowNet shows better performance in improving the diversity and high reward modes simultaneously, whereas conventional optimization approaches like reinforcement learning only focus on the highest reward modes. To handle more BBs, they are represented with binary Morgan fingerprints, and the probability of sampling a particular BB is computed as a function of its fingerprint and the embedding vector of the current state. In addition,  new backward policies have been tested to ensure that the backward trajectories can return to existing building blocks. They also performed a case study that leverages experimental data for more efficient drug design as a practical application. Overall, the manuscript is well-written, but several issues need to be addressed.
weaknesses: 1. The conclusions of some experiments are arguable (see the questions below), which need to be justfied probably with additional experiments for clarity.  
2. Increasing the number of BBs slows down training linearly that may cause a computational burden when dealing with a large amount of BBs.

---

confidence: 4
contribution: 2
initial_rating: 6
presentation: 2
questions: NA
soundness: 2
strengths: - Focus on Synthetic Accessibility: SynFlowNet adopts a framework that effectively addresses a key limitation in generative models for drug design by ensuring that generated molecules are synthetically feasible, which is crucial for practical drug development applications.
- Diverse and Optimized Molecule Sampling: Leveraging GFlowNets, SynFlowNet provides a balanced sampling process that maximizes both reward and structural diversity, effectively mitigating the typical mode collapse encountered in reinforcement learning-based models.
- Validation through Synthetic Accessibility Metrics: In silico experiments demonstrate SynFlowNet’s advantage over fragment-based approaches, particularly in terms of synthesizability, underscoring the practical benefits of synthesis-based molecule generation for drug discovery.
summary: The paper presents SynFlowNet, a model designed to improve the generation of novel drug-like molecules while ensuring their synthetic accessibility. SynFlowNet constructs molecules sequentially using commercially available reactants and documented chemical reactions, effectively constraining the generation process within a synthesizable chemical space. The model integrates a forward synthesis constraint, evaluates synthetic accessibility through external retrosynthesis tools, and optimizes diversity and reward (e.g., binding affinity) without relying on reinforcement learning, which tends to sample fewer, high-reward molecules. SynFlowNet introduces a backward policy for navigating reaction trajectories, improving synthetic feasibility and sampling accuracy. By utilizing molecular fingerprints and scalable building blocks, SynFlowNet can explore vast chemical spaces efficiently, making it adaptable for specific drug targets when provided with experimental fragment screening data, such as in SARS-CoV-2 protease inhibitors.
weaknesses: Limited Testing on Optimization Tasks: While SynFlowNet showcases improvements in SeH proxy, its optimization abilities are not fully tested against established molecular optimization benchmarks. Incorporating PMO benchmark tasks and directly comparing results with other generative baselines, such as those in Gao et al. (2022) [1], could provide a more comprehensive evaluation of SynFlowNet’s optimization performance.

### Reference
[1] Gao, Wenhao, et al. "Sample efficiency matters: a benchmark for practical molecular optimization." Advances in neural information processing systems 35 (2022): 21342-21357.

---

confidence: 4
contribution: 3
initial_rating: 5
presentation: 3
questions: 1. **Page 2, Figure 1B:**
    - As I understand it, the building block selection probability depends on not only the state embedding but also the selected reaction. (Page 4 Line 213). However, the Figure 1B indicates that the selection probability is calculated only from the state information $s_k$. I suggest revising the figure and including the series of equations for the neural network to clarify this process.
    - Also, I was wondering if dot product operation between Morgan fingerprint (binary vector) and state embedding can satisfy the GFlowNet Objective. When the 4-bit fingerprints for building blocks $b_1$, $b_2$, and $b_3$ are $[1,0,0,0]$, $[0,1,0,0]$, and $[1,1,0,0]$, respectively, the probability of $b_3$ can be derived from those of $b_1$ and $b_2$:
      $$\log P(b_3|s) = \frac{\log P(b_1|s) + \log P(b_2|s)}{\sqrt{2}}$$
      I think that this relationship can introduce bias in the forward transition probability.
3. **Page 7, Footnote:**
Is it typo?: '… which are different to the **reactants** used by SynFlowNet …'. As I know, AiZynthFinder allows user to use different reactants. I think 'reactants' should be replaced to 'reactions'.
4. **Page 9, Figure 7C:**
I suggest changing the title "Performance" to "Runtime" or a similar term.
5. **Page 19, Table A.1 (Hyperparameters)**:
Could you clarify the meaning of *max_nodes*? In the standard GFlowNet implementation, it refers to the maximum number of heavy atoms in molecules, but the generated molecules illustrated in the paper have more than 9 heavy atoms.
soundness: 2
strengths: 1. Novel Backward Policy of GFlowNet training: SynFlowNet proposes a novel learning strategy for backward policy to ensure the directed acyclic graph and facilitate synthesizable molecules.
2. Connection to Wet-Lab Experiments: The case study in Section 4.6 provides a clear demonstration of how to connect wet-lab experiments and *in silico* drug design, highlighting the potential real-world applicability of the method.
3. Comprehensive Ablation Studies: The paper provides ablation studies for many design choices.
summary: This paper presents SynFlowNet, a practical approach for drug design that incorporates reaction-based Markov Decision Processes (MDPs) into a Generative Flow Network (GFlowNet) to generate diverse molecules with guaranteed synthetic pathways.
Notably, authors propose a novel backward policy training strategy to ensure the presence of backward trajectories.
Furthermore, the authors investigate the potential of SynFlowNet to be integrated with wet experiments such as fragment-based crystallography.
weaknesses: **Overall Comment**

I acknowledge the novelty of this work in successfully adapting the GFlowNet framework for synthesizable molecular design, but there are several critical aspects of the current manuscript that need to be addressed.
I think this manuscript may be suitable for acceptance after revision, and I’m positive about changing the scores during the rebuttal period.

*\* The issues are sorted in the order they appeared in the manuscript.*

**Issue 1: Lack of novelty / appropriate citations for the major contribution.** Page 2, Line 103, Contribution 3

One major concern is that the manuscript does not contain appropriate citations for this major contribution.
To use 200k building blocks, the authors used action embedding [1], which has already been employed in existing synthesis-oriented molecular generative models [2-4].
In particular, Koziarski et al. [4] introduced action embedding into GFlowNet and investigated its effectiveness in large action spaces.
While the authors use a different fingerprint representation (Morgan fingerprints instead of MACCS fingerprints), Gottipati et al. [2] have already evaluated various representations, including both fingerprints.
Authors should highlight the novel differences from existing approaches or provide appropriate citations.

**Issue 2: Lack of experimental details.** Section 4.1 and Section 4.5
- In Section 4.1, what is the rewarding strategy used for multi-objective optimization (FragGFN SA, SynFlowNet SA, SynFlowNet QED)? The possible reward structures could be Multi-objective GFlowNet (MOGFN) [5] or multiplication-based rewards [6].
- In Section 4.5, I couldn't find the experimental details. The related section is on Page 19, Section A.5, but it does not provide crucial information about this experiment. I suggest to include experimental details and mention it like Section 4.6.
  - Reward Function: What reward function was used in the ablation studies? Is it sEH?
  - Average Times: Clarify what *average times* refers to–e.g., runtime per batch or runtime per molecule.
  - For Figure 7 D-F, entire building block set is used?

**Issue 3: Fairness of evaluation metric.** Page 8, Figure 6

The improvement in the reward/diversity trade-off is impressive, but I have questions about the metric, **Reward/Novelty trade-off**.
Regarding the Reward/Novelty trade-off, it is possible that the improvement in novelty is due to the use of the Enamine building blocks rather than the effectiveness of well-structured generative model architecture.
In this case, the REINVENT's high similarity to ChEMBL would be greatly reduced by simply retraining it on Enamine REAL or ZINC.

For a fair comparison in terms of novelty, the comparison should be performed on the same exploration space:
- Train REINVENT on Enamine REAL.
- Use ChEMBL-like building blocks, such as blocks obtained from the retrosynthetic decomposition [3] of ChEMBL molecules.

**Issue 4: Limited comparative study.** Sections 4.2 and 4.3

The comparison between the proposed method and baselines was performed only on proxy models.
I recognize that proxy models such as sEH [7] are widely used to evaluate the effectiveness of optimization strategies in the optimization community [8], so I think this is a minor issue.
However, I have concerns about the appropriateness of evaluating and comparing the performance of discovering novel molecules based on proxy models trained on a restricted and known chemical space.

**Issue 5: Restricted scalability to a large building block set.** section 4.5

I definitely agree with the argument for employing an extensive building block (BB) set to explore the vast chemical space.
I believe that it can mitigate the reward-diversity trade-off by expanding the sample space.
However, I concern that the current experimental results do not seem to provide sufficient evidence about this argument and the scalability of SynFlowNet.

The authors demonstrated that a larger BB set yielded a higher selection of unique BBs (Figure 7B), but this has a negative impact on the average reward (Figure 7A) and speed (Figure 7C).
Consequently, the current findings seem to emphasize the negative effects of the usage of extensive BB sets on the quality and throughput for generation.
To justify the significance of utilization of unique BBs, I suggest including structural diversity metrics for generated samples, such as the number of unique Bemis-Murcko scaffolds or the average pairwise Tanimoto distance.

---
**Reference**
1. Dulac-Arnold, Gabriel, et al. "Deep reinforcement learning in large discrete action spaces." arXiv preprint arXiv:1512.07679 (2015).
2. Gottipati, Sai Krishna, et al. "Learning to navigate the synthetically accessible chemical space using reinforcement learning." International conference on machine learning. PMLR, 2020.
3. Seo, Seonghwan, Jaechang Lim, and Woo Youn Kim. "Molecular generative model via retrosynthetically prepared chemical building block assembly." Advanced Science 10.8 (2023): 2206674.
4. Koziarski, Michał, et al. "RGFN: Synthesizable Molecular Generation Using GFlowNets." ICML'24 Workshop ML for Life and Material Science: From Theory to Industry Applications (2024).
5. Jain, Moksh, et al. "Multi-objective gflownets." International conference on machine learning. PMLR, 2023.
6. Lee, Seul, Jaehyeong Jo, and Sung Ju Hwang. "Exploring chemical space with score-based out-of-distribution generation." International Conference on Machine Learning. PMLR, 2023.
7. Bengio, Emmanuel, et al. "Flow network based generative models for non-iterative diverse candidate generation." Advances in Neural Information Processing Systems 34 (2021): 27381-27394.
8. Gao, Wenhao, et al. "Sample efficiency matters: a benchmark for practical molecular optimization." Advances in neural information processing systems 35 (2022): 21342-21357.

Meta-review: In this work, authors introduce SynFlowNet, a generative model that combines GFlowNets with reaction-based Markov Decision Processes to generate synthetically accessible drug-like molecules. The key technical innovation is incorporating forward synthesis constraints directly into the generative process by using commercially available reactants and documented chemical reactions. The model introduces novel backward policy training strategies to ensure generated molecules have valid synthetic routes.

The paper demonstrates several important scientific findings. First, SynFlowNet achieves better diversity in generated molecules compared to baseline methods like REINVENT while maintaining high rewards, as shown through extensive experiments with both proxy models and docking scores (Reviewer f18i). Second, the model scales effectively to large building block sets of over 200k compounds, though with expected computational overhead (Reviewer 6JHQ). Third, the backward policy training innovations help ensure generated molecules have valid synthetic routes back to building blocks, addressing a key challenge in synthesis-constrained generation (Reviewers AADH and 6JHQ).
The paper's key strengths include: 1) A well-engineered solution to the critical problem of ensuring synthetic accessibility in molecular generation, validated through multiple metrics and external retrosynthesis tools, 2) Comprehensive empirical evaluation including ablation studies, comparison to strong baselines, and demonstration on practical drug discovery tasks, 3) Clear practical relevance demonstrated through integration with experimental fragment screening data.

The initial submission had some limitations around reproducibility and evaluation completeness, but these were thoroughly addressed during review. The authors added code, conducted additional experiments on molecular rediscovery tasks and sample efficiency (addressing Reviewer f18i's concerns), and validated the novelty claims by testing with ChEMBL-derived building blocks (addressing Reviewers 6JHQ and AADH's concerns). The revised evaluation provides strong evidence for the method's capabilities.

The authors were highly responsive during review, strengthening the evaluation to conclusively demonstrate the benefits of their approach. While there remain some scaling challenges with large building block sets, the current capabilities are sufficient for many practical applications. The combination of technical innovation, thorough evaluation, and clear practical relevance makes this a strong contribution to the field. I recommend acceptance because the paper makes a significant technical contribution to synthesizable molecular generation while maintaining high standards of empirical validation, an important problem that has received lesser attention in the molecule generation community.

Decision: Accept (Spotlight)