Generating diversity and securing completeness in algorithmic retrosynthesis

Abstract

Chemical synthesis planning has considerably benefited from advances in the field of machine learning. Neural networks can reliably and accurately predict reactions leading to a given, possibly complex, molecule. In this work we focus on algorithms for assembling such predictions to a full synthesis plan that, starting from simple building blocks, produces a given target molecule, a procedure known as retrosynthesis. Objective functions for this task are hard to define and context-specific. In order to generate a diverse set of synthesis plans for chemists to select from, we capture the concept of diversity in a novel chemical diversity score (CDS). Our experiments show that our algorithm outperforms the algorithm predominantly employed in this domain, Monte-Carlo Tree Search, with respect to diversity in terms of our score as well as time efficiency. SCIENTIFIC CONTRIBUTION: We adapt Depth-First Proof-Number Search (DFPN) (Please refer to https://github.com/Bayer-Group/bayer-retrosynthesis-search for the accompanying source code.) and its variants, which have been applied to retrosynthesis before, to produce a set of solutions, with an explicit focus on diversity. We also make progress on understanding DFPN in terms of completeness, i.e., the ability to find a solution whenever there exists one. DFPN is known to be incomplete, for which we provide a much cleaner example, but we also show that it is complete when reinforced with a threshold-controlling routine from the literature.

Keywords: Chemical diversity score; Computer-Assisted Synthesis Planning (CASP); DFPN; Retrosynthesis.

Fig. 1

The graph G that shows that DFPN is not complete.The molecule nodes and the reaction nodes are depicted as squares and circles, respectively

Fig. 2

Comparison between DFPN* and MCTS: a Chemical diversity score (CDS) for DFPN* and MCTS routes for different search times. DFPN* routes show significantly higher diversity score throughout all data points (Asymptotic Wilcoxon-Mann-Whitney Test, significance level $\alpha=0.01$, colored boxes correspond to the IQR between the first and third quantile, whiskers go up-to/down-to 1.5 times of the upper/lower bound of the IQR). b For smaller search times (60 s/120 s) the mean number of used reaction over all routes and solved molecules by our DFPN* implementation is 3.5/3.9 vs. 2.6/3.5 for the MCTS implementation. For longer search times this trend reverses and the highest median number of used reactions seen for the DFPN* implementation is 5.1 vs 7.3 for the MCTS implementation (colored boxes correspond to the IQR between the first and third quantile, whiskers go up-to/down-to 1.5 times of the upper/lower bound of the IQR)

Fig. 3

Mean multiplicative viability score for DFPN* and MCTS routes for different search times. The viability scores for each reaction of a route are multiplied with each other to give the estimated route viability. The respective mean score from each target molecule is given here. Higher values can be associated with a higher likelihood that the particular route will work in the lab. (Colored boxes correspond to the IQR between the first and third quantile, whiskers go up-to/down-to 1.5 times of the upper/lower bound of the IQR)

Algorithm 1

search (DFPN with TCA)

Fig. 4

The graph G from Simpler counter example without TCA section in vertical orientation. Recall that the nodes in $V_0$ and $V_1$ are depicted as squares and circles, respectively, and that $W_0=\{v_6\},W_1=\mathrm{\varnothing }$

Fig. 5

The explored part of G right after the expansion of $v_0$

Fig. 6

The explored part of G at the second visit to $v_0$

Fig. 7

The explored part of G at the third visit to $v_0$

Fig. 8

The explored part of G at the fourth visit to $v_0$

Fig. 9

Fraction of the cumulative number of samples relative to the full data set (red)

Fig. 10

Top k accuracy and applicability for models trained on transformation classes reweighted according to their size, resulting in the same importance of all classes, and without such reweighting. Reactions found in in-house lab journals are stronger weighted in both cases and for the Applicability metric reactions were checked with our fast filter model

Fig. 11

Percentage of solved molecules for the DFPN* and MCTS algorithms respectively. Shorter search times favor the DFPN* algorithm

Fig. 12

Number of unique molecules per Route. The differences of the medians for DFPN* and MCTS are statistically significant for all computing times (Asymptotic Wilcoxon-Mann–Whitney Test, significance level $\alpha =0.01$). We believe the much broader IQRs in the data generated by the DFPN* are a result of the higher diversity in said data. We assume a superlinear relation between CDS and the number of unique molecules, yielding wider IQRs than those found in the CDS analysis. Due to its mechanism of finding multiple routes via penalizing nodes on the search graph, the DFPN* algorithm is forced to explore a wider variety of molecules. (Colored boxes correspond to the IQR between the first and third quantile, whiskers go up-to/down-to 1.5 times of the upper/lower bound of the IQR, dots represent outliers.)

Fig. 13

Comparison between DFPN* and MCTS with respect to the minimum and maximum route viability: (left) Minimum multiplicative viability score for DFPN* and MCTS routes for different search times. The viability scores for each reaction of a route are multiplied with each other to give the estimated route viability. The respective minimum score from each target molecule (i.e. the assumably worst route found) is distributed around lower values for the DFPN* than for the MCTS implementation for search times 60 and 120 s. For longer search times both algorithms converge to the approximately same level. (Colored boxes correspond to the IQR between the first and third quantile, whiskers go up-to/down-to 1.5 times of the upper/lower bound of the IQR). (right) Maximum multiplicative viability score for DFPN* and MCTS routes for different search times. The viability scores for each reaction of a route are multiplied with each other to give the estimated route viability. The respective maximum score from each target molecule (i.e. the assumably best route found) is distributed around lower values for the DFPN* than for the MCTS implementation for search times 60, 120, and 300 s. At 600 s the distribution is very similar, while for 900 and 1200~s the trend reverses. (Colored boxes correspond to the IQR between the first and third quantile, whiskers go up-to/down-to 1.5 times of the upper/lower bound of the IQR)