PharmacoMatch: Efficient 3D Pharmacophore Screening via Neural Subgraph Matching
Authors: Daniel Rose
Oliver Wieder
Thomas Seidel
Thierry Langer
Venue: ICLR 2025 Poster (2025)
Keywords: Contrastive Representation Learning
Neural Subgraph Matching
Virtual Screening
Pharmacophore Modeling
Primary area: learning on graphs and other geometries & topologies
TL;DR: PharmacoMatch, a new contrastive learning approach, accelerates pharmacophore screening by encoding query-target relationships in the embedding space.

Abstract: The increasing size of screening libraries poses a significant challenge for the development of virtual screening methods for drug discovery, necessitating a re-evaluation of traditional approaches in the era of big data. Although 3D pharmacophore screening remains a prevalent technique, its application to very large datasets is limited by the computational cost associated with matching query pharmacophores to database molecules. In this study, we introduce PharmacoMatch, a novel contrastive learning approach based on neural subgraph matching. Our method reinterprets pharmacophore screening as an approximate subgraph matching problem and enables efficient querying of conformational databases by encoding query-target relationships in the embedding space. We conduct comprehensive investigations of the learned representations and evaluate PharmacoMatch as pre-screening tool in a zero-shot setting. We demonstrate significantly shorter runtimes and comparable performance metrics to existing solutions, providing a promising speed-up for screening very large datasets.

Published as a conference paper at ICLR 2025
PHARMACOMATCH! EFFICIENT 3D PHARMACOPHORE SCREENING VIA NEURAL SUBGRAPH MATCHING
Daniel Rose!.2, Oliver Wieder Thomas Seidel & Thierry Langer
'Department of Pharmaceutical Sciences, Division of Pharmaceutical Chemistry, Faculty of Life Sciences, University of Vienna, 1090 Vienna, Austria 2 Christian Doppler Laboratory for Molecular Informatics in the Biosciences, Department of Pharmaceutical Sciences, University of Vienna, 1090 Vienna, Austria Vienna Doctoral School of Pharmaceutical, Nutritional and Sport Sciences, University of Vienna, 1090 Vienna, Austria
Email: rdse wiever,ieder, thomas.seidel thierry.langer @ @univie.ac.at
ABSTRACT
The increasing size of screening libraries poses significant challenge for the development of virtual screening methods for drug discovery, necessitating evaluation of traditional approaches in the era of big data. Although 3D phar- macophore screening remains prevalent technique, its application very large datasets is limited by the computational cost associated with matching query phar- macophores to database molecules. In this study, we introduce PharmacoMatch, a novel contrastive learning approach based on neural subgraph matching. Our method reinterprets pharmacophore screening as an approximate subgraph match- ing problem and enables efficient querying of conformational databases by encod- ing query-targe target relationships in the embedding space. We conduct comprehen- sive investigations of the learned representations and evaluate PharmacoMatch as pre-screening tool in zero-shot setting. We demonstrate significantly shorter runtimes and comparable performance metrics to existing solutions, providing a promising speed-up for screening very large datasets.
INTRODUCTION
The immense scale of the chemical space, covering over 1060 small organic molecules (Virshup al., 2013), makes the identification of molecules that bind to protein target particularly challeng- ing. Virtual screening methods are important tools in computer-aided drug discovery, addressing this complexity and supporting medicinal chemists navigate large molecular databases search of potential hit compounds (Sliwoski et al., 2014). Among these methods, established approach 3D pharmacophore screening. Pharmacophores represent an ensemble steric and electronic molecular features that is necessary to (Wermuth et al., 1998). A pharmacophore query can, for instance, be generated from the interac- tion profile of ligand-protein complex and then aligned positionally with the three-dimensional conformations compounds in a database (Wolber & Langer, 2005). These molecules are sub- sequently ranked based on their agreement with the pharmacophore query, with those exhibiting similar pharmacophoric patterns retrieved as potential hit compounds (Wolber et al., 2006).
optimal interaction with specifie biological target
An important development of the last years has been the emergence of make-on-demand libraries, such as Enamine REAL (Shivanyuk et al., 2007). These libraries contain billions of commercially available compounds that can be rapidly synthesized and are continuously expanding, driven by advances in synthetic accessibility (Ilanos et al., 2019). Screening larger libraries increases the chances of identifying hits, but it also brings the trade-off of longer screening times. Scaling up 3D) pharmacophore screening to handle billions of molecules is challenging due to the computational cost of pharmacophore alignment (Warr et al., 2022). Despite efforts to optimize these algorithms (Wolber et al., 2008; Permann al., 2021) and the development of various pre-filtering techniques (Seidel et al., 2010), the alignment step remains critical bottleneck.
*Corresponding author

---

Published as a conference paper at ICLR 2025
Phermacophore virtual screening
PhermacoMatcc
Ligand
Encoder
Embedding vectors
Conformer generation
Pharmacophore generation
Decision function
Hittisst
Encoder
Query
Ligand-pratein complex
Pharmacophore generation
Query
embedding
Figure 1: Overview of the PharmacoMatch workflow: Conformer and pharmacophore generation from ligands and query creation, for example from ligand-protein complex, precede pharma- cophore screening. The encoder model converts the screening database into embedding vectors, stored for later use. hitlist is generated by comparing the query embedding with the database embeddings.
In this work, propose perform 3D pharmacophore screening by using learned representations. Specifically, our PharmacoMatch model employs a graph neural network (GNN) encoder to map 3D pharmacophores into order embedding space (Ying al., 2020), and predicts pharmacophore matching through vector comparisons. The embedding vectors for the screening database are com- puted once and then used to quickly generate hitlist based the query embedding (Figure 1). Our key contributions are as follows.
We conceptualize pharmacophore matching representation learning problem, enabling the use of order embeddings for efficient and scalable pre-screening.
We develop GNN encoder that generates meaningful vector representations from 3D pharmacophores. The model is trained in a self-supervised manner on unlabeled data, employing contrastive loss objective to capture partial ordering relationships between queries and targets in the learned embedding space. We design augmentation strategies specifically suited for the task of pharmacophore matching.
We thoroughly analyze the learned embeddings and validate the practical utility of our ap- proach as an effective pre-screening tool through experiments on virtual screening bench- mark datasets.
We empirically demonstrate that our method is significantly faster than existing solutions, offering a practical advantage for screening billion-come compouund libraries. A promising direction for further improving efficiency is combining our approach with approximate retrieval techniques. Our work thus represents key step towards vector databases for virtual screening. Additionally, our refor- mulation of pharmacophore screening neural subgraph matching directly captures its subgraph matching nature. To the best of our knowledge, we are the first to introduce order embeddings in virtual sereening, providing novel perspective for the field. Furthermore, learning presents a promising strategy for addressing data scareity in drug discovery. We hope that the combined insights from study contribute to advancing representation learning for virtual screening applications.
of sulf-supervised
RELATED WORK
Pharmacophore alignment algorithms Alignment algorithms compute rigid-body transforma- tion, the pharmacophore alignment, to match query's pharmacophoric feature pattern to database ligands. A scoring function then evaluates the pharmacophore matching by considering both the number of matched features and their spatial proximity. The alignment is typically preceded by fast filtering methods that prune the search space based on feature types, feature point counts, and quick distance checks. Only molecules that pass these filters undergo the final, computationally
2

---

Published as a conference paper at ICLR 2025
expensive 3D alignment step, which is usually performed by minimizing the root square de- viation (RMSD) between pairs of pharmacophoric feature points (Dixon et al., 2006; Seidel et al., 2010). The algorithm by Wolber et al. (2006) creates smoothed histograms from the neighborhoods of pharmacophorie feature points for pair assignment using the Hungarian algorithm, followed by alignment with Kabsch's method (Kabsch, 1976). A recent implementation by Permann et al. (2021) improves on runtime and accuracy by using search strategy that maximizes pairs of matching phar- macophoric feature points. Alternatively, shape-matching algorithms like ROCS (Hawkins et al., 2007) and Pharao (Taminau et al., 2008) model chemical features by Gaussian volumes, optimizing for volume overlap.
Machine learning for virtual screening A common approach to using machine learning for vir- tual screening is to train models on measured bioactivity values. However, these models are con- strained by the scarcity of experimental data, which is both costly and challenging to obtain (Li et al., 2021). Unsupervised training of target-agnostic models for virtual screening avoids depen- dence on labeled data, but remains relatively unexplored. DrugClip (Gao al., 2023), which is not based on pharmacophores, approaches virtual screening as a similarity matching problem be- tween protein pockets and molecules, using multi-modal learning approach where a protein and a molecule encoder create a shared embedding space for virtual screening. Sellner et al. (2023) used the Schrödinger pharmacophore shape-screening score to train a transformer model on phar- macophore similarity, which is different objective than pharmacophore matching. PharmacoNet (Seo Kim, 2023) a pre-screening tool that uses instance segmentation for pharmacophore gen- eration in protein binding sites and employ deep learning for pharmacophore modeling, but not for the alignment nor matching.
algorithm for binding pose estimation. They
PRELIMINARIES
Pharmacophore representation In this work, we treat 3D pharmacophores as attributed point clouds (Mahé et al., 2006; Kriege & Mutzel, 2012). A pharmacophore P can be represented by set of pharmacophoric feature points P = {(r,, d:) E T x D}; with the Cartesian coordinates r; and the descriptor d, of the pharmacophoric feature point pt- The descriptor set 'D covers the following pharmacophoric feature types: hydrogen bond donors (HBD) and acceptors (HBA), halogen bond donors (XBD), positive (PI) and negative electrostatic interaction sites (NI), hydrophobic interaction sites (H), and aromatic moieties (AR). Directed feature types like HBD and HBA can be associated with vector component, but for simplicity, we will omit this information in our study. We further denote the set of pair-wise distances between feature points as R = <I|r: - ryl12 | S i,j < IPI). The pharmacophore P can be represented as a complete graph G(P) (Vp, Ep, Ap), where Vp = .o...... 1| } denotes the set of nodes with node attributes Ap(1:) 14, and Ep Vp x Vp denotes the set of edges, with the edge attribute of exj defined as the pair-wise Euclidean distance Ap(eiy)= = |r, - 122 between the positions of nodes and j. The edges are undirected, edge eij can be identified with edge Eje- The label set = U R is the union of the feature descriptor set and the set of pair-wise distances, and Ap represents a labelling function À : VU E + L that assigns a label to the corresponding vertex v or edge e. This representation is invariant to translation and rotation.
Subgraph matching Two graphs G = (V1, FI, A1) and G2 = (V2, Ez, d2) are isomorphic, denoted by G1 = G2 if there exists edge-preserving bijection f Vi + V/2 such that V(u, v) E Eg. Additionally, we require the preservation of node and edge labels, such >z(s(v)), and V(u, v) E Xi ((u, v))
Ei : (f(u.,
>2((f(u), f()))) Let Go = (Vo, Eq, Xa) be a query graph, G7 (Vr, Em,,r, Ar) a larger target graph, and GH = EV,,EH, Xh) a subgraph of GT such that VH S Vr, and Eн S Er. The objective of subgraph matching is to decide, whether GQ is subgraph isomorphic to GT, denoted by GQ S Gr. requiring the existence of a non-empty set of subgraphs = (G G5 z Go} that are isomorphic to Go-
that Vo E Vi (v)
Pharmacophore matching In its most general setting, pharmacophore matching seeks to match all pharmacophoric feature points of a query pharmacophore Pq with the corresponding feature points of larger target pharmacophore Pr. Let Pи S Pr denote a subset of the feature points of Pr. Then Po matches Pr after alignment if there exists a bijection g : Po - Pn such that di E Pq= d; dg(s) and [r, - r2 PT, where rr is the radius of tolerance sphere. It
3

---

of a given pharmacophore Published as a conference paper at ICLR 2025
is thereby sufficient that query feature points mapped into the tolerance sphere of their target counterpart. For simplicity, we assume the same tolerance radii among all pharmacophoric feature points. The ultimate goal of pharmacophore matching is to retrieve molecules from a database. A matching pharmacophore is always linked to a corresponding ligand molecule via a look-up table.
When represented graphs Go = G(P)), Gи = G(PH), and GT = G(PT), this task boils down the node-induced subgraph matching of a query pharmacophore graph Go to a target pharmacophore graph Gr. The tol- erance sphere, however, weakens the require- ment on edge label matching. An approximate matching Aq((u, v)) s (u), f()))) is sufficient if the difference between )o((u, v)) and Aи ((f(v) f(u))) is less than 2rt, where Tr represents the tolerance radius of each phar- macophoric feature point. This ensures that the query points fall within the tolerance spheres of the target points (compare Figure 2). Our problem formulation of pharmacophore match- ing relies on relative distances instead of the ab- solute positioning of pharmacophoric features
Volld Matching
Invalid Matching
abo I
val
2ny
abot In
Ir r5)) Bry
Target
Query
Fr
Illustration of the pharmacophore matching objective: The aim is to match the phar- macophoric points of a query with the correspond- ing points of a target pharmacophore such that the query points fall within the tolerance sphere of the target points, with a tolerance radius TT.
Figure 2:
and is therefore independent of prior alignment.
METHODOLOGY
Overview In the following we introduce PharmacoMatch, a novel contrastive learning framework with the aim to encode query-target relationships of 3D pharmacophores into an embedding space. We propose to train a GNN encoder model in a self-supervised fashion, as illustrated in Figure 3. Our model is trained database (Davies et al., 2015; Zdrazil et al., 2023) and learns pharmacophore matching solely from augmented examples, comparing positive and negative pairs of query and target pharmacophore graphs, while optimizing an order embedding loss to extract relevant matching patterns.
approximately 1.2 million unlabeled small molecules from the ChEMBL
Unlabeled data for contrastive training To span the pharmaceutical compound space, we down- load set of drug-like molecules sourced from the ChEMBL (2024) website in the form of Simplified Molecular Input Line Entry System (SMILES) strings (Weininger, 1988) and curate an unlabeled dataset using the open-source Chemical Data Processing Toolkit (CDPKit) (Seidel, 2024) (see Ap- pendix A. I for details). After an initial data clean-up, which includes the removal of solvents and counter ions, adjustment of protonation states to physiological pH, and climination of duplicate structures, the dataset contains approximately 1.2 million small molecules. To ensure a zero-shot setting in validation experiments, remove all molecules from the training data that also ap- pear in the test sets. Finally, we generate a low-energy 3D conformation and the corresponding pharmacophore for each ligand.
Model input We represent the node labels (X,(
graph G(P) = (Vp, Ep, Xp) as en-ooed (OHE) feature vectors h = =(h (hi..., p()). We employ a distance encoding to represent pair-wise distances, which was inspired by the SchNet ar- chitecture (Schütt et al., 2018). The edge attributes of edge Cav are derived from the edge label Ae(eu)) and represented by a radial basis function 6₁ (r, r.) = where centers μ4 were taken from uniform grid of K points between zero and the distance cutoff at 10 Å, and the smoothing factor 3 represents a hyperparameter. To this end, the pharmacophore P represented by a data point x = [h. e] which is a tuple of the feature matrix h [PP)X|P| and the distance-encodings e E (P|X
- re ll )2)
GNN encoder
The encoder input is the pharmacophore graph representation x h, el, with the feature matrix and the edge attributes e. Node feature embeddings are generated by initially passing the OHE feature matrix through single dense layer without activation function. We
4

---

Published as a conference paper at ICLR 2025
3. Model training
Augmantation
f 9
tol3 g
5 (RGy aT)
Shared weichts
Augmentation
P
Faa
Pru
Prp
F6.
Rg
61) ~ 3.
4(+) m. M
b. Embedding geometry
F
Pqa
Pq.
P
Dimension
Target region
Encoder
Rot embedtinng
Target embedding
Py: Target pharmacaphore
Augmantation
((--) Loss fundtion
Pg: Query pharmacophore
quary-aanel
Off-target
Pg Negative query phamacophhree
Negadive quary-largel pair
P4: Negative target pharmacophore
Figure 3: (a) The encoder model learns order embedding space by comparing augmented phar- macophores. (b) Illustration of the embedding space, where pharmacophores matching a query are positioned to the upper right. (c) Augmentation strategies for model training involve generating positive and negative query-target pairs on-the-fly by combining node deletion with varying degrees of node displacement. Negative pairs are also created by shuffling the batch, mapping query phar- macophores to random target pharmacophores.
then update the node representations through message passing using the edge-conditioned convo- lution operator (NNConv) by Gilmer al. (2017); Simonovsky & Komodakis (2017), which was originally designed for representation learning on point clouds and 3D molecules, to aggregate dis- tance information into the learned node representations (see Appendix A.3 for details). We connect successive convolutional layers using NenseNet-style skip connections (Huang et al., 2017). Graph- level read-out is achieved by additive pooling of the updated feature matrix h E R |P into graph representation q E Rmh. which then projected to the final output embedding z E RR. by multi-layer perceptron. The employed loss function requires to map the final representation to the non-negative real number space. accomplish this by using the absolute values the learnable weights for the last linear transformation immediately after the final ReLU unit (see Appendix A.4 for details).
Loss function In order to encode query-target relationships of pharmacophores into the embed- ding space, we employ the loss function by Ying et al. (2020). The key insight that subgraph relationships can be effectively encoded in the geometry of order embedding space through partial ordering of the corresponding vector embeddings. Let Zq the embedding of graph Go, 2r R. a GNN encoder to map pharmacophore graphs g to embedding vectors E R.D The partial ordering 2o N Zr reflects, whether GQ is subgraph isomorphic GT: z S zr [i], Vi e {1.... D} iff GQ Z Gr. The following max max-margin objective can be used to train the GNN encoder fe on this relation:
the embedding of graph Gr, and fe
L(zq, ZT)
Σ (2q.2))EPos
E(zq, ZT)
Σ (eQ)2))ENeg
max(0, Elzq, 2T)}
XP + R4 reflects violation of the partial ordering the embedding vector pair: E(zq. Zr) | Zq zT 11.. Pos is the set of positive pairs per batch, these are pairs of query zq and target graph embedding Zy with a subgraph-supergraan is the set of negative examples, these pairs of query and target embedding vectors that violate this relationship. The positive and negative pairs are generated on-the-fly via augmentation during training.
The penalty function E R.
relationship. and Neg
5

---

Published as a conference paper at ICLR 2025
Augmentation module The PharmacoMatch model correlates the matching of query and tar- get pharmacophore with the partial ordering of their vector representations. Positive pairs represent successful matchings, while negative pairs serve as counter examples. In order to create these pairs from unlabeled training data, we define three families of augmentations T, which are composed of random point deletions and positional point displacements. For positive pairs, valid queries created by randomly deleting some nodes from a pharmacophore P. leaving at least three, and dis- placing the remaining nodes within tolerance sphere of radius rT. This augmentation, denoted as (2(-) Ti produces the positive pair (t₂ (P), P). Negative pairs highlight examples of unsuccess- ful matching, using three strategies to capture different undesired outcomes. The first strategy intro- duces positional mismatches by displacing the pharmacophorie feature points of P to the boundary of the tolerance sphere without deleting any points. Specifically, we compute the mean position of
the points, PP >!' Fí, and shift each point p; by IT along the direction r; p.. This approach avoids the unintended creation of positive pairs, which can occur with random sampling, and ensures that displacements do not cancel out. The resulting augmentation, denoted t²(-) ~ T2. generates the negative query-target pair (t2 (P), P). This strategy demonstrated better model performance than random positional displacement. Our second strategy teaches the model that every pharmacophoric feature point in the query should correspond to a point in the target. This is achieved by deleting some target nodes, using an augmentation operator ta(-) ~ Ta. where Ta involves node deletion without displacement. As result, the query in the pair target. With the third strategy, we train the model to avoid matching queries with targets that are significantly different. This approach involves randomly mapping queries (P:) the incorrect targets P,, where i # (for more details, see Appendix A.2).
(P)) only partially matches its
Model details & design choices Our GNN encoder model is implemented with three convolu- tional layers with an output dimension of 64. The MIP has a depth of three dense layers with a hidden dimension of 1024 and output dimension of 512. The final model was trained for 500 epochs using an Adam (Kingma, 2014) optimizer with a learning rate of 10 10-3 The margin of the best performing model
set to a 100. The default tolerance radius rr in CDPKit's pharma- cophore screening is set to 1.5 A, and we use the same value for the node displacement during model training ensure consistency with the alignment algorithm in subsequent evaluations. We design a curriculum learning strategy for learning on pharmacophore graphs, detailed in Appendix A.5, along with details on model training and hyperparameter optimization. Systematic ablation studies of model architecture components and the impact of different augmentation strategies, model size, and embedding dimension on model performance are provided in Appendix A.5.
Decision function for model inference We use the trained GNN encoder Se to precompute vector embeddings Zy of the database pharmacophores. These are queried with the pharmacophore embed- ding zq by verification of the partial ordering constraint through penalty function E(.,,) (Equation 1), which shall not exceed a threshold t. This leads to decision function g I x R - {0, 1):
iff Ef 2T) < t (1
g(zq, ZT)
otherwise
evaluating 1, if the partial ordering on Zo and zT reflects a pharmacophore matching, and 0 oth- erwise. In the following, we refer to equation (2) as matching decision function. In practice, we recommend a decision threshold of t = 6500, determined from our benchmark experiments.
EXPERIMENTS
We designed the embeddings to reflect the type and relative positioning of pharmacophoric feature points. Comparison of embedding vectors via the matching prediction function should emulate the matching of the underlying pharmacophores. To get a better understanding of the encoder's latent space, we investigate these properties as follows:
1. Pharmacophoric feature point perception: We investigate the learned embedding space qualitatively through dimensionality reduction.

---

Published as a conference paper at ICLR 2025
2. Positional perception: We investigate the influence of positional changes the output of the matching decision function.
3. Virtual sereening performance: The performance our model is evaluated using ten DUD-E (Mysinger et al., 2012) targets, and the produced hitlists are compared with the performance and runtime of the CDPKit (Seidel, 2024) alignment algorithm. We further evaluate the pre-screening performance of our model and compare it against the Pharma- coNet (Seo & Kim, 2023) model using the DEKOIS2.0 (Bauer et al., 2013) and LIT-PCBA (Tran-Nguyen et al., 2020) benchmark datasets.
Benchmark datasets We perform experiments on the DUD-E benchmark dataset (Mysinger et al., 2012), which is commonly used to evaluate the performance of molecular docking and structure- based screening. The complete benchmark contains 102 protein targets, each accompanied by active and decoy ligands in the form of SMILES strings (Weininger, 1988) and the PDB template (Burley et al., 2017) of the ligand-receptor complex. We randomly select ten different protein targets for proof-of-concep comparison with the CDPKit alignment algorithm. For our wee-screnning experi- ment,
the DEKOIS2.0 (Bauer et al., 2013) dataset, which contains 80 targets, each with 40 actives and 1,200 decoys, as well as the LIT-PCBA (Tran-Nguyen et al., 2020) dataset, consisting of 15 target sets with 7,761 confirmed actives and 382,674 inactive compounds. Database ligands are processed according to the data curation pipeline outlined in the Methodology section, with the exception that we sample up to 25 conformations per compound. Further information is provided in the Appendix A.6.
5.1 PHARMACOPHORIC FEATURE POINT PERCEPTION
We conduct a qualitative analysis through dimensionality reduction to gain a first intuition for the properties of the learned embedding space. The partial ordering of graph representations in the embedding space, based on the number of nodes per graph, is essential for encoding query-target relationships. This ordering property of the embedding space can be visualized using principal component analysis (PCA). Figure 4a displays the first two principal component axes of the learned representations, with the representations labeled according to the number of pharmacophoric feature points of the corresponding pharmacophores. This visualization demonstrates how the embedding vectors are systematically ordered relative to the number of nodes in each pharmacophore graph. Similarly, the Uniform Manifold Approximation and Projection (UMAP) algorithm (McInnes et al., 2020), a dimensionality reduction technique that preserves the local neighborhood structure of high- dimensional data, was employed. Figure 4b shows the UMAP representation of the embeddings, labeled by the number feature points of a specific type. This visualization suggests that pharma- cophores with similar set of points are mapped proximally within the embedding space.
5.2 POSITIONAL PERCEPTION
We define family of augmentations Trp to randomly delete nodes from a pharmacophore P and displace the remaining nodes by a radius rp. We sample augmentations (.)
Tro with increasing radius TD taken from a uniform grid of m distances between 0 and 10 Å. For given batch of pharmacophores (P,,.., Pa). we generate the query-target pairs ((rr)(P)), PI)..., (tro (Pn), p.)]. We then evaluate the decision function g(..)) (Equation 2) on the corresponding vector representations and calculate the mean of the decision function across all pairs against an increasing radius Tp. which is illustrated in Figure 4c. Without node displacement, the mean matching decision function is close to 1, indicating that the model recognizes pharma- cophores with reduced node sets as valid queries. With a displacement of approximately 1.5 A, the mean matching decision value drops to 50% demonstrating the model's 's consideration of the chosen tolerance radius. Beyond a displacement of 1.5 A, the decision function further decreases, ap- proaching a plateau at approximately 6 Å. The results show that our model integrates 3D- 3D-positional information of pharmacophoric feature points into the learned representations.
5.3
ALIGNMENT PREDICTION AND PERFORMANCE AS A PRE-SCREENING TOOL.
Hitlist ranking Each benchmark set is comprised of a pharmacophore query Po and a set of ligands L. (Tn,.., .., J. where each ligand L, is associated with set of pharmacophores

---

Published as a conference paper at ICLR 2025
Dimensionaltty
Dimensionality roduction
band perception
Hydrophszki
Douays
PC2
02-
Number poins
Figure 4: (a.) Dimensionality reduction of the ADA target's embedding space via PCA, with em- beddings labeled by pharmacophoric feature point count. (b.) Dimensionality reduction via UMAP, with embeddings labeled by pharmacophoric feature point type. (c.) Experimental validation of the model's perception of 3D point positions, showing the mean matching decision function versus the displacement radius rp of the augmentation, with a decision threshold set to / 6500.
(PI..., Pk. h and label yt. which indicates whether the ligand active or decoy. The task to rank the database ligands w.t.. the query, based scoring function F x P 7 R.. The ranking score vi ligand L; calculated through aggregation of the pharmacophore scores ((FF(PQ, P.)..., F(Pq. ht)). where ① is aggregation operator. PharmacoMatch H + RP and evaluates the penalty function RF x ? + R.. A low penalty corresponds to high ranking. The ranking score of database ligand Lę is calculated as 8/s = mmin((F(zz)) 21)..., E(zq, )aa)
transforms the query G(Pq) + zo and the set of pharmacophores (G.P)).... {zi....,Zyy via encoder model fe +
Ground truth for alignment prediction We evaluate the alignment prediction performance of PhatmacoMatch by relative comparison with the alignment algorithm implemented in the open- source software CDPKit (Seidel, 2024), which utilizes clique-detection followed by Kabsch align- ment (Kabsch, 1976). The ligand-receptor complex of the respective DUD-E target used gen- erate structure- brased pharmacophore query with the CDPKit. The CDPKit alignment algorithm only returns exact matches, meaning that queries with an excessive number points may yield no results. Consequently, pharmacophore modeling often requires user interaction to reduce the num- ber of points the query. For our comparison, we refined the initial query pharmacophores to a subset of 5-7 pharmacophoric feature points, which is common range in pharmacophore mod- eling. These points were selected to ensure that the query yields meaningful enrichment for the CDPKit algorithm. Only after this refinement did we proceed to compare PharmacoMatch against this ground truth. The alignment of a query Po and target Pr is evaluated with an aligument score S P P + R,, which takes into account the number of matched features and their geometric fit (further details are provided in the Appendix A.6). The ligand ranking score is calculated as d) = max ((S((((SPP))) P1). S(Pq. PP.) Fi), the highest alignment score represents the score for the database ligand. Analogous to equation (2), we can also define matching decision function ø based on the alignment score, where = \Pol:
Pr) = f PT)
(3)
otherwise
Evaluation metries Both algorithms rank database ligands to produce a hitlist. We assess the performance of PharmacoMatch the benchmark using two approaches. First, we demonstrate that the PharmacoMatch penalty .(..) correlates with the matching decision function o(- of the alignment algorithm. We evaluate both functions against all pharmacophores in dataset w.r.f query Pq- The outputs are compared by generating the corresponding receiver operating character- istic (ROC) curves, and the performance is quantified using the area under the ROC curve (AUROC) metric. We argue that this relative performance is the key metric for evaluating the alignment pre- diction performance of our model. Second, we compare the absolute screening performance of
8

---

Published as a conference paper at ICLR 2025
Table 1: Method comparison and screening performance of the PharmacoMatch algorithm and the CDPKit alignment algorithm ten different DUD-E protein targets (see Appendix A.6 for details). BEDROC values are calculated with a = 20, as recommended by Truchon & Bayly (2007), AUROC and BEDROC are reported in percent. Confidence intervals are calculated us- ing bootstrapping (Efron, 1979), with standard deviations reported based on 100 resampled datasets.
Protein
Relative performance
Absolute screening performance PharmacoMatch
Absolute screening performance
target
CDPKit
AUROC
AUROC BEDROC
1F
EF5%
EF10%
AUROC BEDROC
EF5%
EF EF10%
ACES
96.0
= 0.2
58
2
18
844
1.4
3.
0.3
0.2
55 +
16 2
5.5
1.3
3.0
2.1
2.1000 0.2
ADA
98.3
1 0.2
83
1
44
16.7
4.1 9.5
1.0
5.7
0.4
9411
82
13
53.6
4.3 15.9
1..9
8.4
0.4
ANDR
98.8
0.1 0.1
76
11
2
15.8
43
71
26
22
12.6=21
4.4
3.7 0.3
EGFR
91.1
63
=
11:1
3.1
2.3
1.6
70
26
122
4.6
3.2
FAI0
84.3
111
0.2
0.1
0.1
55
611
0.0
0.7
1.2
KIT
81.5
56
1=1
0.0
0.2
0.7
2
2
1.1008
1.2
1.3
PIK1
79.5 2
62
9
1.5
0.7
1.8
75
3
5.7
10.2
5.5
SRC
96.6 - 0.2
79
27
6.0
5.3
=
4.6
80
28
28
11.1
11.111
5.3
0.4
4.3
THRB
89.9
= 0.5
70
=1
22
6.9
1.0
4.8
0.4
3.3
#
0.2
79
35
22
11.8
=
1.5
7.2
+
0.4
4.5
0.2
UROK
83.8
1 0.2
60
1
2
4
0.6
0.7
0.5
0.2
0.4
1
0.2
91
1
55
1
24.5
2.8
10.4
0.9
8.2
0.4
model and the alignment algorithm using the ligand ranking score th. The primary objective of vir- tual screening is to find active compounds among decoys. The AUROC metric is used to evaluate the overall classification performance w.r.t. activity label y/. A drawback of this metric that does not reflect the early enrichment of active compounds in the hitlist, which is of significant in- terest in virtual screening. Early enrichment is assessed using the enrichment factor (EF,,) and the Bonzzmann-enhanced discrimination of ROC (BEDROC) metric (Truchon & Bayly, 2007), which assigns higher weights to better-ranked samples (definitions in Appendix A.6). Note that these per- formance metries are entirely dependent on the chosen query. Therefore, our primary performance metric in this experiment is the relative performance AUROC. For completeness, we also report absolute virtual screening metries. We conduct this proof-of-concept comparison using randomly selected subset of 10 DUD-E targets, serving as focused case study.
Alignment prediction performance Our results, comparing PharmacoMatch with the CDPKit alignment algorithm across the selected targets, are summarized in Table (ROC plots are provided in the Appendix A.6). We observe robust correlation between the hitlists generated by the two algorithms, demonstrating the effectiveness of our approach. This correlation varies by target, re- flecting the sensitivity of virtual screening the chosen query. Although the alignment algorithm achieves generally higher AUROC scores and early enrichment, our method consistently produces hitlists with competitive performance across several targets. In terms of runtime, PharmacoMatch significantly outperforms the alignment algorithm. We compare the time required for alignment, em- bedding, and vector matching per pharmacophore. Aligmment performed in parallel on an AMD EPYC 7713 64-Core Processor with 128 threads, while pharmacophore embedding and matching are run on an NVIDIA GeForce RTX 3090, with both devices having comparable purchase prices and release dates. Creating vector embeddings from pharmacophore graphs takes 67 + 7 jus per pharmacophore. approximately as long as aligning a query to a target with 66 j.s. However, the embedding process only needs to be performed once. Subsequently, the preprocessed vector data can be used for vector matching, which takes 0.3 + 0.1 jas, being approximately two orders of magnitude faster than the alignment. Additionally, vector comparison is independent of the query size, an advantage not shared by the alignment algorithm. Although executed on different hardware, this comparison highlights the speed-gain of our algorithm.
Applicability as a pre-screening tool We demonstrated how PharmacoMatch predicts pharma- cophore matching and compared its performance with the CDPKit alignment algorithm. A key aspect of this comparison is the role of user interaction in query design, since the alignment algo- rithm only identifies match when all pharmacophoric feature points align, and while it is possible to define feature points as optional, this flexibility significantly increases runtime due to the com- binatorial explosion of possible query patterns. Consequently, handling larger queries with ten or

---

Published as a conference paper at ICLR 2025
Table 2: Comparison of PharmacoMatch with PharmacoNet (Seo & Kim, 2023) as a pre-screening tool. The reported performance metrics were averaged over all targets of the DEKOIS2.0 and LIT-PCBA datasets. AUROC and BEDROC are given in percent, BEDROC values are calculated with a 80.5, as reported in Seo & Kim (2023). PharmacoNets runtime per ligand depends on the number of atoms, here we report their measurement for ligands with 70 heavy atoms. The runtime for PharmacoMatcc was calculated as average over all ligands in the benchmark.
DEKOIS2.0
LIT-PCBA
Runtime
AUROC
BEDROC
FFo.5%
EF1%
EFs%
AUROC
BEDROC
EFo
EF1%
EFs%
per ligand (s)
PharmacoNet
62.5
12.3
4.4
4.2
2.9
3.1
5.2
PharmacoMatch (ours)
60.9
15.1
5.5
3.2
57.4
5.0
3,5
2.2
3.3
more points becomes computationally prohibitive. PharmacoMatch operates independently of the query's size. This characteristic makes PharmacoMatch particularly suitable as a pre-screening tool. By using structure-based queries without refinement step. it can filter extensive datasets ef- fectively, reducing the computational burden before engaging more resource-intensive methods or requiring user interaction. To evaluate this case, we compare our method with the recent pre- screening tool PharmacoNet (Seo & Kim, 2023), which employs image segmentation to generate pharmacophore queries and a parameterized analytical function for hidlst creation. To replicate an automated pre-screening scenario, we exclude user interaction from our workflow. Specifically, we generate an interaction pharmacophore using CDPKit and use it directly as query for Pharma- coMatch. We assess the pre-screening performance of PharmacoMatch and PharmacoNet on the DEKOIS2.0 (Bauer al., 2013) and LIT-PCBA (Tran- Fraan-Ngyenn et al., 2020) benchmark datasets. Performance is evaluated using the average AUROC, BEDROC, and EF metrics across all targets (see Appendix A.6 for more details). The results of this comparison are summarized in Table 2. Our evaluation demonstrates that PharmacoMatch outperforms PharmacoNet with respect to early en- richment. Importantly, PharmacoMMatc'ss runtime is three orders of magnitude faster. Given that the primary goal of pre-screening is rapid and cost-effective filtering before using computationally expensive methods, we argue that this substantial improvement in runtime makes PharmacoMatch an effective pre-screening tool.
Practical considerations & limitations There are two options for integrating our model into virtual screening pipeline. First, the PharmacoMatch model can be used in place of the alignment algorithm to generate hitlist of potential active compounds, which is suitable for quickly producing a compound list for experimental testing. Second, our method can serve as an efficient geesscrennnng tool for very large databases, reducing the number of molecules from billions to millions. after which the slower alignment algorithm can be applied to this filtered subset. Note that alignment will still be necessary if visual inspection of aligned pharmacophores and corresponding ligands is desired. An important limitation is that the embedding process is not lossless and might lead to reduced 3D geometric precision, when compared to an alignment algorithm. The employed B(3)-invariant encoder cannot distinguish pharmacophore from its mirror image, potentially increasing the false positive rate. Addressing these limitations will be part of future work.
6 CONCLUSION
We have presented PharmacoMatch, contrastive learning framework that creates meaningful phar- macophore representations for virtual screening. The proposed method tackles the matching of 3D pharmacophores through vector comparison an order embedding space, thereby offering a valu- able method for significant speed-up of virtual screening campaigns. PharmacoMatch is the first machine-learning based solution that approaches pharmacophore virtual screening via an approxi- mate neural subgraph matching algorithm. We confident that method will help improve on existing virtual screening workflows and contribute to the assistance of medicinal chemist the complex task of drug discovery.
10

---

Published as a conference paper at ICLR 2025
REPRODUCIBILITY STATEMENT
The source code of this project can be found under the following link: https://githut com/ mo melinfo - viienaa/Pha http///aaaaaaatteeeeeeoooo
Please follow the instructions in the repository for re- production of our results. Training and data can be downloaded here: https: /tti. org/ 10. 60.4604/m00 m9. wwt....................... figshare. 27061081
ACKNOWLEDGEMENT
We thank Roxane Jacob, Nils Morten Kriege, and Christian Permann from the University of Vi- enna, Klaus-Jürgen Schleifer from BASF SE, and Andreas Bergner from Sestainnge--ngglheeee RCV GmbH & Co KG for fruitful discussions and proof-reading of the manuscript. We also thank the anonymous reviewers for their suggestions and comments. Financial support received for the Christian Doppler Laboratory for Molecular Informatics in the Biosciences by the Austrian Federal Ministry of Labour and Economy, the National Foundation for Research, Technology and Devel- opment, the Christian Doppler Research Association, Boehringer-Ingelleinn RCV GmbH & Co KG and BASF SE is gratefully acknowledged.
AUTHOR CONTRIBUTIONS
Conceptualization: DR, OW. Methodology: DR, OW. Data curation & analysis: DR. Code im- plementation model training: DR. CDPKit software & support: TS. Investigation: DR. Writing (original draft): DR. Writing (review and editing): DR, OW, TS, TL. Rebuttal: DR. Funding ac- quisition: TL. Resources: TL. Supervision: OW, TL. All authors have given approval to the final version of the manuscript.
REFERENCES
Matthias R. Bauer, Tamer M. Ibrahim, Simon M. Vogel, and Frank M. Boeckler. Evaluation and optimization of virtual screening workflows with dekois 2.0 public library of challenging docking benchmark sets. J. Chem. Inf. Model., 53(6):447-11462 2013.
Shaked Brody, Uri Alon, and Eran Yahav. How attentive are graph attention networks? arxiv preprint arxiv:2 /4491, 2022.
Stephen K. Burley, Helen M. Berman, Gerard J. Kleywegt, John L. Markley, Haruki Nakamura,
and Sameer Velankar. Protein Data Bank (PDB): The Single Global Macromolecular Structure Archive, pp. 627-641. Springer New York, New York, NY, 2017. ISBN 978- 1990--9-99770
ChEMBL. Chembl database embl-ebi. https: htop:::/www.... ebi ac uk/ htwp/chemb//..mm//// 2024. Ac- cessed: 2024-01-06.
Mark Davies, Michat Nowotka, George Papadatos, Nathan Dedman, Anna Gaulton, Francis Atkin- son, Louisa Bellis, and John P. Overington, ChEMBL web services: streamlining access to drug discovery data and utilities. Nucleic Acids Res., W3621):W612-W620 2015.
Steven L. Dixon, Alexander M. Smondyrev, Eric H. Knoll, Shashidhar N. Rao, David E. Shaw, and Richard Friesner. Phase: new engine for pharmacophore perception, 3d qsar model development, and 3d database screening: 1. methodology and preliminary results. J. Comput. Aided Mol. Des., 20:647-671, 2006.
B. Efron. Bootstrap Methods: Another Look at the Jackknife. Ann. Stat., 7(1):1 26, 1979.
William Falcon and The Torch Lightning team. PyTorch Lightning, March 2019. URL https: //github. com Lig.tiittgngggggggmgmgggggggg.nt AI hal//iggtning..
Matthias Fey and Jan Eric Lenssen. Fast graph representation learning with pytorch geometric. arxiv
preprint arxiv: wwwxiv1903... 0244822288.,,.. 2019.
Bowen Gao, Bo Qiang, Haichuan Tan, Yinjun Jia, Minsi Ren, Minsi Lu, Jingjing Liu, Wei- Wei-Ying Ma, and Yanyan Lan. Drugelip: Contrastive protein molecule representation learning for virtual screening. In Advances in Neural Information Processing Systems, volume 36, pp. 44595 44595-44614 44614. Curran Associates, Inc., 2023.
11

---

Published as a conference paper at ICLR 2025
Justin Gilmer, Samuel S. Schoenholz, Patrick F. Riley. Oriol Vinyals, and George E. Dahl. Neural message passing for quantum chemistry. In International conference on machine learning, pp- 1263-1272 PMLR, 2017.
Paul C. D. Hawkins, A. Geoffrey Skillman, and Anthony Nicholls. Comparison of shape-matching and docking as virtual sereening tools. J. Med. Chem., 50(1):74-82, 2007.
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Deep residual learning for image recog- nition. arxiv preprint 1rxiv:1512.03338. 2015.
Weihua Hu, Bowen Liu, Joseph Gomes, Marinka Zitnik, Percy Liang, Vijay Pande, and Jure Leskovec. Strategies for pre-training graph neural networks. arxiv preprint arxiv:1905.12260 2020.
Gao Huang, Zhuang Liu, Laurens Van Der Maaten, and Kilian Q. Weinberger. Densely connected convolutional networks. Proceedings of the IEEE conference on computer vision and pattern recognition, pp. 4700-4708, 2017.
Wolfgang Kabsch. A solution for the best rotation to relate sets of vectors. Acta Crystallogr A, 32(5):922-925 222,, 1976.
Diederik Kingma. Adam: method for stochastic optimization. arXiv preprint arXiv: 2014.
Nils Kriege and Petra Mutzel. Subgraph matching kernels for attributed graphs. In Proceedings of the 29th International Conference on Machine Learning, pp. 291-298, 2012.
Hongjaan Li, Kam-Heung Sze, Gang Lu, and Pedro Ballester. Machine-learning scoring functions for structure-basee based virtual screening. Wiley Interdiscip. Rev. Comput. Mol. Sci., 11(():e1478, 2021.
Christopher A. Lipinski, Franco Lombardo, Beryl W. Dominy, and Paul J. Feeney. Experimen- tal and computational approaches to estimate solubility and permeability in drug discovery and development settings. Adv. Drug Deliv. Rev., 23(1):3-25, 1997.
Eugenio J. Llanos, Wilmer Leal, Duc Luu, Jürgen Jost, Peter F. Stadler, and Guillermo Restrepo. Exploration of the chemical space and its three historical regimes. Proceedings of the National Academy of Sciences, 16662)):
2019.
Pierre Mahé, Liva Ralaivola, Véronique Stoven, and Jean-Philippe Vert. The pharmacophore kernel for virtual screening with support vector machines. J. Chem. Inf. Model, 46(5):2003-2014 2006.
Leland McInnes, John Healy, and James Melville. Umap: Uniform manifold approximation and projection for dimension reduction. arxiv preprint 1rxiv:180..03426. 2020.
Michael M. Mysinger, Michael Carchia, John Irwin, and Brian K. Shoichet. Directory of useful decoys, enhanced (dud-e): better ligands and decoys for better benchmarking. J. Med. Chem., 55 (14):6522-659) 2012.
Christian Permann, Thomas Seidel, and Thierry Langer. Greedy 3-point search (g3ps) a novel algorithm for pharmacophore alignment. Molecules, 26(23):7201, 2021.
Kristof T. Schütt, Huziel E. Sauceda, Pieter-Jan Kindermans, Alexandre Tkatchenko, and Klaus- Robert Müller. Schnet-a deep learning architecture for molecules and materials. Chem. Phys., 148(24), 2018.
Thomas Seidel. Chemical data processing toolkit source code repository. https://github. com/molinfo-ooeooooooooocooooooe www.nnfovvinnaa/aaaaaaaae CCPKit, 2024. Accessed: 2024-01-06.
Thomas Seidel, Gökhan Ibis, Fabian Bendix, and Gerhard Wolber. Strategics for 3d pharmacophore- based virtual screening. Drug Discov. Today Technol., (14))=2211e222) 2010.
Thomas Seidel, Christian Permann, Oliver Wieder, Stefan M. Kohlbacher, and Thierry Langer. High-quality conformer generation with conforge: algorithm and performance assessment. J. Chem. Inf. Model, 63(17):5549-5577 2023.
12

---

Published as a conference paper at ICLR 2025
Manuel S. Sellner, Amr H. Mahmoud, and Markus A. Lill. Enhancing ligand-based virtual screening with 3d shape similarity via a distance-aware transformer model. bioRxiv preprint 2023.
bioorri::2023.1111..775355
Seonghwan Seo and Woo Youn Kim. Pharmaconet: Accelerating large-scale virtual screening by deep pharmacophore modeling, arxiv preprint arxiv:2310.0068) 2023.
Chao Shen, Xujun Zhang, Chang- Yu Hsieh, Yafeng Deng, Dong Wang, Lei Xu, Jian Wu, Dan Li, Yu Kang, Tingjun Hou, al. A generalized protein potein- scoring framework with balanced scoring, docking, ranking and screening powers. Chem. Sci., 14(30):812 1298811292222 8146, 2023.
Alexander N. Shivanyuk, Sergey V. Ryabukhin, A. Tolmachev, A. V. Bogolyubsky, D. M. Mykytenko, A. A. Chupryna, W. Heilman, and A. N. Kostyuk. Enamine real database: Mak- ing chemical diversity real. Chemistry today, 25(6):58 59,88999000 2007.
Martin Simonovsky and Nikos Komodakis. Dynamic coneicondd filters in convolutional neu- ral networks on graphs. In Proceedings of the IEEE conference on computer vision and pattern recognition, pp. 3693-3702, 2017.
Gregory Sliwoski, Sandeepkumar Kothiwale, Jens Meiler, and Edward W. Lowe. Computational methods in drug discovery. Pharmacol, Rev., (6(1):334-36) -35,--95,, 2014.
Jonatan Taminau, Gert Thijs, and Hans De Winter. Pharao: pharmacophore alignment and optimiza- tion. Mol. Graph. Model., 2712):161-169, 2008
Viet-Khoa Tran-Nguyen, Célien Jacquemard, and Didier Rognan. Lit-poba: an unbiased data set for machine learning and virtual screening. J. Chem. Inf. Model., 60(9):4263-427) 42733 2020.
Jean-François Truchon and Christopher I. Bayly. Evaluating virtual screening methods: good and bad metrics for the "early recognition" problem. J. Chem. Inf. Model, 47(2):488-508, 2007.
Aaron M. Virshup, Julia Contreras Garcia, Peter Wipf, Weitao Yang, and David N. Beratan. Stochas- tic voyages into uncharted chemical space produce a representative library of all possible drug-like like compounds. J. Am. Chem. Soc., 135(19):7266-7306 2013.
Wendy A. Warr, Marc C. Nicklaus, Christos A. Nicolaou, and Matthias Rarey. Exploration of ultralarge compound collections for drug discovery. Chem. Inf. Model, 62(9)::021-2033 2022.
David Weininger. Smiles, a chemical language and information system. 1. introduction to method- ology and encoding rules. J. Chem. Inf. Comput., 28(1):31-36, 1988.
Camille-Georges Wermuth, C. R. Ganellin, Per Lindberg, and L. A. Mitscher. Glossary of terms used in medicinal chemistry (iupac recommendations 1998). Pure Appl. Chem., 70(5): 1129-
1143, 1998.
Gerhard Wolber and Thierry Langer. Ligandscout: 3-d pharmacophores derived from protein-bound ligands and their use virtual screening filters. J. Chem. Inf. Model, 45(1):160-16)) 2005.
Gerhard Wolber, Alois A. Dornhofer, and Thierry Langer. Efficient overlay of small organic molecules using 3d pharmacophores. Comput. Aided Mol. Des., 20(12))773-788, 2006.
Gerhard Wolber, Thomas Seidel, Fabian Bendix, and Thierry Langer. Molecule-pharmacophore superpositioning and pattern matching in computational drug design. Drug Discov. Today, 13 (1-2):23-29, 2008.
Rex Ying, Zhaoyu Lou, Jiaxuan You, Chengtao Wen, Arquimedes Canedo, and Jure Leskovec. Neural subgraph matching. arxiv preprint arxiv:2007.00 .3092, 2020.
Barbara Zdrazil, Eloy Felix, Fiona Hunter, Emma J. Manners, James Blackshaw, Sybilla Corbett, Marleen de Veij, Harris loannidis, David Mendez Lopez, Juan F. Mosquera, Maria Paula Magari- nos, Nicolas Bosc, Ricardo Arcila, Tevfik Kizilören, Anna Gaulton, A. Patrícia Bento, Melissa F. Adasme, Peter Monecke, Gregory A. Landrum, and Andrew R. Leach. The ChEMBL Database in 2023: a drug discovery platform spanning multiple bioactivity data types and time periods. Nucleic Acids Res., 52(D1):D1 180-D1192, 2023.
13

---

Published as a conference paper at ICLR 2025
A APPENDIX
A.1 DATASET CURATION & STATISTICS
Unlabeled training data was downloaded from the ChEMBL database represent small molecules with drug-like like properties. At the time of data download, the ChEMBL database contained 2,399, 7,73 unique compounds. We constrained the compound category to "small molecules" and enforced ad- herence to the Lipinsky rule of five (Lipinski et al.. 1997). specifically setting violations to "0," resulting in a refined set of 1,348, 115 compounds available for download. The molecules were ac- quired in the form of Simplified Molecular Input Line Entry System (SMILES) (Weininger, 1988) strings. Subsequent to data retrieval, we conducted preprocessing using the database cleaning func- tionalities of the Chemical Data Processing Toolkit (CDPKit) (Scidel, 2024). This process involved the removal solvents and counter ions, adjustment of protonation states to a physiological pH value, and elimination of duplicate structures, where compounds differing only in their stereo figuration were regarded as duplicates. To prevent data leakage, we carefully removed all structures from the training data that would occur in one of the test sets used for our benchmark experi- ments. The final set was comprised of 1 (,211,088 compounds. For each compound within the dataset, 3D conformation was generated using the CONFORGE (Seidel et al., 2023) conformer generator from the CDPKit, which was successful for 1,220, 104 compounds. To enhance batch diversity, we generated only one conformation per compound for contrastive training. Subsequently, 3D pharma- cophores were computed for each conformation, with removal of pharmacophores containing less than four pharmacophoric feature points. The ultimate dataset comprised 1,217,361 distinct pharma- cophores. Figure 5 shows the frequency of pharmacophores with specific pharmacophoric feature point count in the training data. On average, a pharmacophore consists of 12.4 pharmacophoric fea- ture points, with the largest pharmacophore in the dataset containing 32 points. Hydrophobic feature points and hydrogen bond acceptors are the most prominent, while hydrogen bond donors and aro- matics occur less frequently. Ionizable feature points and halogen bond donors comparatively rare.
1.0
0.8
Frequency
0.6
0.2
# Feature points
# Hydrophobic points
Aromatic points
Negative lonizable points
0.8
Frequency
0.6
0.4
0.2
0.0
10
0
10
10
10
# Positive ionizable points
# Hydragen bond donor points
Hydragen bond acceptor points
Ha Halogen bond donor points
Figure 5: Pharmacophoric feature point statistics of the training data. The respective histograms display the total number of pharmacophoric feature points and the number of points of specific types per pharmacophore in the training data. The complete training dataset contains 1,217,361 distinct pharmacophores.
14

---

Published as a conference paper at ICLR 2025
A.2 AUGMENTATION MODULE
The augmentation module receives the initial pharmacophore Xo [bo, rol, with the initial OHE feature matrix hg and the Cartesian coordinates ro- Edge attributes of the complete graph were calculated from the pair-wise distances between nodes after modifying the input according to the augmentation strategy, which combines random node deletion and random node displacement. The module outputs the modified tuple
[h, el with the feature matrix h and the edge attributes
Node deletion Random node deletion involved removing at least node, with the upper bound determined by the cardinality of the set nodes V; of graph G,- To ensure the output graph retained least three nodes, the maximum number of deletable nodes [V]] The number of nodes to delete was drawn uniformly at random.
Node displacement There are two modes for the displacement of pharmacophoric feature points, displacement within the tolerance sphere, and displacement onto the surface of the tolerance sphere. For simplicity, assumed the same tolerance sphere radius TT across different pharmacophoric feature types. For displacement within the tolerance sphere, we calculated coordinate displace- ment (Ax, Ay, A2) from spherical coordinates & ((0, 2TT) and cos U(-1, 1), which were drawn at random from a uniform distribution:
A Arsin co Ay Ar Osin Φ, A: Ar cos
where Ar = Fr V and u ~ ((00 1). Displacement the nodes onto the tolerance sphere surface was achieved by caleulating the mean of the positions of the pharmacophoric feature points, / 1. Σ ra and displacing each point p; by rr the direction r; μ.. The displacement away [T
from the center ensures that displacement directions do not cancel each other.
A.3 MESSAGE PASSING NEURAL NETWORK
Convolution irregular domains like graphs is formulated as message passing, which generally be described as:
11(h) (1)
(k-1)
(h (h! (1(()1) .,
((t111), eij))
(5)
JEN(S)
where h,"h E FR denotes the node features of node at layer k, h (R-1) E RR denotes the node features of node at layer 1, eij E R.. the edge features of the edge from node to node j. (8) and g)(k) are parameterized, differentiable functions, and + is an aggregation operator like, g., the summation operator (Fey & Lenssen, 2019). In encoder architecture, employed the following condiiiodd convolution operator, which was proposed both by Gilmer et al. (2017) and Simonovsky & Komodakis (2017):
h (k)
Oh
(k-1) +
h (k-1) (ei))) JEN()
Σ
(6)
where θ IRF denotes learnable weights and (+) R - - FXF" denotes neural network, in our case an MI.P with one hidden layer. These transformations map node features h into a latent representation that combines pharmacophoric feature types with distance encodings.
A.4 ENCODER IMPLEMENTATION
The encoder was implemented as a GNN fe : 9 + RD that maps a given graph G to the abstract representation vector E R.. The architecture is comprised of an initial embedding block, three subsequent convolution blocks, followed by pooling layer, and projection block.
15

---

Published as a conference paper at ICLR 2025
Embedding block The embedding block receives the pharmacophore graph G; the tuple x, [h, e,}, with the OHE feature matrix h, and the edge attributes e,, Initial node feature embeddings created from the OHE features with fully-connected (FC) dense layer with learnable weights W and bias b:
h C Wh; +
Convolution block The convolution block consists graph convolution layer, which is imple- mented as edge-conditioned convolution operator (NNConv), the update rule is described in Section A.3. The network further consists of batch normalization layers (BN), GELU activation functions, and dropout layers. The hidden representation h of graph G; is updated at block / as follows:
[b e:] - (NNConv + BN GELU + concat h/) - dropout) h
where hy represents the latent representation after activation. Updating the feature matrix / times yields the final node representations of the pharmacophoric feature points.
Pooling layer We employed additive pooling for graph-level read-out r,, which aggregates the set of [V] node representations (hı,-.., hivl): of Graph G3 by element-wise summation:
M
qi
h&
Projection block The projection block maps the graph-level read-out to the positive real number space and is implemented as a multi-layer perceptron MLP F - R.P where is the dimension of the vector representation before and 1) the dimension after the projection. The block consists of sequential layers of FC layers, BN, ReLU activation, and dropout:
q + [FC + BN - ReLU + Dropout] + q. k+1
(10)
The final layer is FC layer without bias and with positive weights, only:
2, - abs W
(11)
Matrix multiplication of the positive learnable weights W and the output the last ReLU activation function produces the final representation z, E R
Compleee grash
Graphh cansodatin
Gruphiieend nentioui
Projection
84
1024
NOCO
OHE gharmanahhooi
Node iestuee umbbddiinge
Nodess
fealuros
1000
HerCow Batch nem Astimion
Desped
Carlesion mooddiaaiess
Eoge pall-waa dilance encodings
Global add tive pooling
M.P
Figure 6: Architecture of the GNN encoder model.
16

---

Published as a conference paper at ICLR 2025
A.5 MODEL IMPLEMENTATION AND TRAINING
Implementation dependencies The GNN was implemented in Python 3.10 with PyTorch (v2.0.1) and the PyTTrch Geometric library (v2.3.1) (Fey & Lenssen, 2019). Both, model and dataset, were implemented within the Pa TyTorch Lightning (Falcon & The PyTorch Lightnnng team, 2019) frame- work (v2.1.0). Model training was monitored with Tensorboard (v2.13.0). CDPKit (v1.1.1) was employed for chemical data processing. Software was installed and executed on a Rocky Linux (v9.4) system with x86-64 architecture.
Model training Training was performed on a single NVIDIA GeForce 3090 RTX graphics unit with 24 GB GDDR6X. Training runs were performed for a maximum of 500 epochs with batch size 256 pharmacophore graphs. Curriculum learning was applied by gradual enrichment of the dataset with increasingly larger pharmacophore graphs. At training start, only pharmacophore graphs with 4 nodes were considered. After 10 subsequent epochs without considerable minimization of the loss function, pharmacophore graphs with one additional node were added to the training data. This approach allows the model to start with very simple examples, gradually increasing the difficulty of the matching task. The loss function was minimized with the Adam (Kingma, 2014) optimizer, we further applied gradient clipping. A training on the full dataset took approximately 48 hours with the above hardware specifications.
Hyperparameter optimization & model selection Hyperparameters were optimized through random parameter selection, the tested ranges are summarized in Table 4. Unlabeled data was split into training and validation data with a 98:2 ratio. Training runs were compared using the AUROC value on the validation data. This was calculated by treating the positive and negative pairs as bi- nary labels, and the predictions were based their respective order embedding penalty, which was calculated with function E(--) ) in Equation (I). Hyperparameter optimization was performed on a reduced dataset with 100,000 graphs, which took approximately hours per run. The best perform-
ing models were retrained on the full dataset. The hyperparameters of the final encoder model are summarized in Table 3. After model selection, the final model performance was tested on virtual screening datasets.
Table 3: Hyperparameters of the best performing encoder model
Hyperparameter
batch size
256
dropout convolution block
0.2
dropout projection block
0.2
max. epochs
500
hidden dimension convolution block
64
hidden dimension projection block
1024
output dimension convolution block
1024
output dimension projection block
512
learning rate optimizer
0.001
margin for negative pairs
100.0
number of convolution blocks depth of the projector MI.P edge attributes dimension sampling sphere radius positive pairs
3
5
1.5
sampling surface radius negative pairs
1.5
Ablation studies To evaluate the importance of various model parameters, we conduct a series of ablation studies using the best-performing model. these experiments, we systematically alter one parameter at time and assess its impact classification performance using the validation hold-out set. Our findings reveal several key insights. The embedding dimension of the learned representa- tions can be reduced to 128 without loss in performance. The encoder requires least 32 dimensions to remain effective. Skip-connections are critical for model performance, with DenseNet-styke con- nections slightly outperforming ResNet-style (He et al., 2015) alternatives. Interestingly, the choice
17

---

Published as a conference paper at ICLR 2025
Table 4: Tested hyperparameter ranges for model training.
Hyperparameter
Parameter range
dropout
[0.2, 0.3, 0.4, 0.5] [0.1. 0.5, 1, 2, 5, 10, 100, 1000] [64, 128, 256, 512, 1024] [0.25, 0.5, 1.0, 1.5]
margin for negative pairs output dimension projection block
displacement sphere radius vy of positive pairs
Table 5: Tested hyperparameter ranges for ablation studies.
Parameter
Parameter range
Validation AUROC
Tolerance radius
[0.0, 0.5, 1.0, 1.5, 2.0]
[0.91, 0.93, 0.94, 0.94, 0.93]
Encoder dimension
[8. 16, 32, 64, 96]
[0.92, 0,93, 0.94, 0,94, 0.94]
Embedding dimension
[32, 64, 128, 256, 512, 1024]
[0.91, 0.93, 0.94, 0.94, 0.94, 0.93]
Skip-connection
[dense, res, none]
[0.94, 0.93, 0.74]
GNN Layer
INNConv, GINE, CFConv, GAT]
[0.94, 0.94, 0.94, 0.93]
Projector layers
[1, 2, 31
[0.94, 0,94, 0.94]
Convolution layers
[1, 31
[0.94, 0.94, 0.94]
Margin
[0.01, 0.1, 2, 10, 100, 1000] [0.88, 0.90, 0.92. 0.92, 0.93. 0.93, 0.94, 0.94]
of message- puessageppasing layer, whether NNConv, graph isomorphism operator (GINE) (Hu et al., 2020), graph attention operator (GAT) (Brody et al., 2022), or continuous-filte convolutional layers (CF- Conv) (Schütt et al., 2018), has minimal impact on performance. The depth of the projector and encoder also does not significantly affect results. In contrast, the margin value plays significant role in model performance. While larger values enhance performance, excessively high margins can lead to training instability. margin 100 provides an effective balance between these factors. The displacement radius for augmentations in creating positive pairs is most effective at 1.5 Å, while removing node displacement degrades model performance.
A.6 VIRTUAL SCREENING
DUD-E dataset details General information about the DUD-E targets is summarized in Table 6. For each target we downloaded the receptor structure from the PDB and created the corresponding interaction pharmacophore with the CDPKit. Vector features were converted into undirected phar- macophoric feature points. The resulting pharmacophore queries (Figure 7) were used in our virtual screening experiments.
Table 6: DUD-E targets that were selected for bechmarking experiments this study.
Active Ligands
Active Conformations
Decoy Ligands
Decoy Conformations Points
Query
Target PDB code Ligand ID)
ACES
le66
HUX
451
10048
26198
567122
6
ADA
2elw
FR6
2166
5448
125035
ANDR
2am9
TES
269
3039
14333
211968
6
EGFR
2rgp
HYZ
541
12468
35001
755017
FA10
3k16
443
537
13343
28149
638831
5
KIT
3g0e
B49
3703
10438
224364
5
PLKI
2owb
626
107
2531
6794
152999
6
SRC
3c18
PD5
523
11868
34407
737864
6
THRB
1ype
UIP
461
11494
26894
626722
7
UROK
1sqt
UI3
162
3450
9837
199204
6
18

---

Published as a conference paper at ICLR 2025
AGES
ADA
ANDR
EGFR
FA10
KI
FLK1
SRC
THRB
UROK
Figure 7: Structure-basse based pharmacophore queries of targets of the DUD-E benchmark dataset.
DEKOIS2.0 dataset details General information about the DEKOIS2.0 targets (Bauer et al., 2013). Each target is associated with PDB four-letter code and corresponding ligand ID. For each target, we downloaded the receptor structure and its respective ligand from the PDB and generated the interaction pharmacophore using the CDPKit. These queries were used without further refine- ment in our pre-screening experiments. For targets containing the small molecule ligand in multiple binding pockets, we randomly selected one pocket for pharmacophore generation. The SIRT2 target was excluded from evaluation because its structure does not contain a small molecule ligand. Actives and decoys were processed analogous to the ligands in the DUD-E benchmark.
LIT-PCBA dataset details For each target, the downloaded files include SMILES file contain- ing the active and inactive ligands, along with receptor structures provided as PDB files. To ensure a fair comparison with PharmacoNet, created pharmacophore queries using the same receptor PDB files as in their study. These PDB files were selected based on a methodology described by (Shen al., 2023) and are listed in Table S4 of the respective supporting information. Interaction phar- macophores were generated using CDPKit and used without further refinement in our pre- gre-screening experiments. Active and inactive ligands were processed following the same protocol as for the DUD-F benchmark.
CDPKit alignment scoring function The CDPKit implements alignment as a clique-detection algorithm and computes rigid-body transformation via Kabsch's 's algorithm to align the pharma- cophore query Pa to the pharmacophore target Pr. The goodness of evaluated with a geometric scoring function S: P R
S(2,, Pr) SMFP(PQ, Pr) + Sceom (Pq, Pr)
(12)
where SMFP PxP - N counts the number of matched feature pairs and SGeom PxP evaluates their geometric fit.
(0,
Runtime measurement We measured alignment runtimes using the psdscreen tool from the CDP- Kit with 128 threads on an AMD EPYC 7713 64-Core Processor, while embedding and matching runtimes with PharmacoMatch were recorded using an NVIDIA GeForce RTX 3090 GPU with 24 GB GDDR6X. We additionally used 8 CPU worker nodes for the embedding step. Runtime per pharmacophore was estimated by dividing the total runtime by the number of pharmacophores in each dataset, with the final estimate taken as the mean of ten runs. The results report the mean and standard deviation of these estimates across all ten datasets.
Performance metrics
The enrichment factor (EF) is formally defined
14 n
F.
(13)
with n,, the number of true active compunds in the top a% of the hitlist ranking and the total number of active compounds
19

---

Published as a conference paper at ICLR 2025
The Boltzmann Enhanneed Discrimination of Receiver Operating Characteristic (BEDROC) metric (Truchon & Bayly. 2007) is similar to the AUROC metric but assigns greater weight to samples ranked higher in the hitlist. The degree of this reweighting is controlled by the early recognition parameter a. The formal definition is as follows:
Σ Ei= a
Rs sinh (a/2) R. (1 cosh (a/22) cosh (a/22 aRR.)
BEDROC,
(14)
ea(1 R.)
where n the number of active compounds in the dataset, N is the total number of compounds, R n//V is the ratio of active compounds in the dataset, and r; rank,-1 the normalized ranking of the i-th active compound.
The area under receiver operating characteristic (AUROC) performance metrics our virtual screening experiments the selected DUD-E targets are derived from the ROC curves presented in Figures 8 and Figure 9. The performance metrics for the pre-screening experiments are calculated as follows: for each target, we calculated the AUROC, BEDROC. and enrichment factors (EF) for the top 0.5%. 1%, and 5% of the hit list, reporting the average values for these metrics. Metrics for PharmacoNet were reported as presented in their original manuscript and were recalculated.
A.7 EMBEDDING SPACE VISUALIZATION
UMAP visualization UMAP embeddings for visualization plots were calculated with the UMAP Python library. The 'metric' parameter was set to Manhattan distance, all other parameters are the default settings of the implementation. We tested range of hyperparameters ensure that the visualization results are not sensitive to parameter selection.
20

---

Published as a conference paper at ICLR 2025
ACCSTALLEA +0.0
ADA 2OK.OLOOK -U00
AYCA EMOROC 10M1
ES6 LMAROC- EAD
19.00 JARON -0.05
They pwsti
DNUROO C3R)
PUC SAUICC CLTB)
LROK: JUSIOC -U.0
03
24
05
24
84
02
30
02
044
64 03 13
Figure 8: Performance comparison of PharmacoMatch and the alignment algorithm. The ROC- curves display the agreement of the hitlist ranking of the two algorithms for ten targets of the DUD-E benchmark dataset.
Kilo m11
Falle
DAMNK 2011
Figure 9: Absolute screening performance of PharmacoMatch and the alignment algorithm perfor- mance for ten targets of the DUD-E benchmark dataset. The pharmacophore queries were generated from the respective PDB ligand-receptor structures.
21

---

Published as a conference paper at ICLR 2025
Actives/ Decoys
Ildrrophobic
Aromatic
Hydrogen bond denor
Decays
0
22
0
0
Hydrogen bond acceptor
Positive ionizable
Negative lonizable
Hallgeen bond donor
12
Figure 10: UMAP visualization of the vector embeddings of the ACES target.
Actives Decoys
Hydrophobic
Aromatic
Hydrogen bond donor
Decoys Actives
Hydrogen bond acceptor
Positive ionizable
Negative ionizable
Hallgenn bond donor
Figure 11: UMAP visualization of the vector embeddings of the ANDR target.
Decoys
Hydrophobic
Hydragen
Decoys Actives
0
0
0
Hydrogen bond acceptor
Positive ionizable
Negative ionizable
Hallogen bond donor
Figure 12: UMAP visualization of the vector embeddings of the EGFR target.
22

---

Published as a conference paper at ICLR 2025
Actives/ Decoys
Ilyrrophobic
Aromatic
Hydrogen bond denor
Decays Actives
0
12
0
0
13
Hydrogen bond acceptor
Positive ionizable
Negative lonizable
Hallgen bond donor
18
Figure 13: UMAP visualization of the vector embeddings of the FA target.
Actives Decoys
Hydrophobic
Aromatic
Hydrogen bond donor
Deccos Actives
Hydrogen bond acceptor
Positive ionizable
Negative ionizable
Hallgenn bond donor
Figure 14: UMAP visualization of the vector embeddings of the KIT target.
Decoys
Hydrophobic
Hydragen
Decoys
Actives
0
10
0
0
Hydrogen bond acceptor
Positive ionizable
Negative lonizable
Hallogen bond donor
Figure 15: UMAP visualization the vector embeddings of the PLKI target.
23

---

Published as a conference paper at ICLR 2025
Actives/ Decoys
Ilyrrophobic
Aromatic
Hydrogen bond denor
Decays Actives
0
16
0
0
Hydrogen bond acceptor
Positive ionizable
Negative lonizable
Halogen bond
16
Figure 16: UMAP visualization of the vector embeddings of the SRC target.
Actives Decoys
Hydrophobic
Aromatic
Hydrogen bond donor
Decoys Actives
Hydrogen bond acceptor
Positive ionizable
Negative ionizable
Hallgenn bond donor
16
Figure 17: UMAP visualization of the vector embeddings of the THRB target.
Decoys
Hydrophobic
Hydragen
Derays
Actives
0
11
0
0
10
Hydrogen bond acceptor
Positive ionizable
Negative ionizable
Hallogen bond donor
Figure 18: UMAP visualization of the vector embeddings of the UROK target.
24


Reviews:
[{"review_id": "y1yq3iZ1lI", "reviewer": "Reviewer_GwAp", "summary": "This paper introduces a novel contrastive learning approach based on neural subgraph matching, i.e., PharmacoMatch, and the authors claim that it reinterprets pharmacophore screening as an approximate subgraph matching problem and enables efficient querying of conformational databases by encoding query-target relationships in the embedding space.", "questions": "See comments in the Weaknesses part.", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "+ This paper is well-written and nicely organized. \n+ The proposed framework is novel and is nicely motivated.", "weaknesses": "- I suggested the authors consider conducting more experiments over other datasets instead of only DUD-E.\n- I wonder if the proposed contrastive learning approach can be applied to other domain datasets?\n- This paper does not provide any unique and novel insights about why the proposed architecture that works well for the current dataset.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "h9ql1bnT26", "reviewer": "Reviewer_p5tV", "summary": "The authors introduce PharmacoMatch, a deep learning approach that reframes 3D pharmacophore screening as a neural subgraph matching problem, using a graph neural network trained through contrastive learning to encode pharmacophore structures into an embedding space. Their method achieves comparable screening performance to traditional alignment-based approaches while being approximately two orders of magnitude faster in matching speed, making it particularly valuable for screening extremely large molecular databases. The model is trained in a self-supervised manner on over 1.2 million unlabeled molecules from ChEMBL, learns to encode both structural and spatial relationships of pharmacophoric points, and demonstrates robust performance across multiple DUD-E benchmark datasets in a zero-shot setting.", "questions": "1. How sensitive is the model to the choice of tolerance radius (r_T = 1.5Å)? Could you provide an analysis?\n2. How did you select the 10 DUD-E targets? Could you demonstrate the method's robustness on more targets?\n3. What is the impact of different augmentation strategies? For example, how much does performance degrade if you remove node deletion or displacement?\n4. Could you compare PharmacoMatch with recent methods like DrugClip and PharmacoNet on the same benchmarks?\n5. The speed advantage is clear, but can the author better justify the conceptual novelty? (Contrastive learning is not new, GNN for molecules is well-established and order embeddings have been used before)\n6. How much 3D geometric precision is actually lost during embedding? Could you quantify the tradeoff between speed and accuracy?\nI wonder if there are specific types of 3D arrangements where the method may fail.", "limitations": "", "rating": 5, "confidence": 4, "soundness": 2, "presentation": 2, "contribution": 2, "strengths": "The acceleration is impressive, with a thorough evaluation of embeddings and runtime analysis. The model provides a practical impact for screening billion-compound libraries. Besides, the reformulation of pharmacophore screening as neural subgraph matching is creative combining self-supervised training approach using augmentation strategies.", "weaknesses": "Despite mentioning recent works like DrugClip and PharmacoNet in the related work section, there are no direct comparisons with these methods. It seems the paper only compares against the traditional CDPKit alignment algorithm, missing comparisons with other more current learning approaches. Besides, there is no comparison or detailed discussion with simpler baseline models (e.g., basic GNN architectures without contrastive learning)\n\nThe discussion on the model details is not sufficient. Lacks systematic ablation studies of model architecture components (e.g., the impact of different GNN layers, the importance of skip connections). Missing analysis of the impact of different augmentation strategies on model performance and investigation of how embedding dimension and model size affect performance, as well as contrastive loss function impacts.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "vTRGiHXAnl", "reviewer": "Reviewer_MX9K", "summary": "The authors propose a contrastive learning approach that emphasizes augmentation strategies by incorporating the concept of pharmacophore in neural subgraph matching. Furthermore, they apply this concept to ligand-based virtual screening, demonstrating that the results are well-aligned with CDPKit’s alignment algorithm. Although the learned representations effectively capture the proposed pharmacophore concepts, the performance in virtual screening—a primary objective of the model—does not appear sufficiently strong. This is primarily due to the lack of benchmarking against other models and the use of additional evaluation metrics.", "questions": "## **Methods:**\n**Q1. Clarification of Pharmacophore Representation:** In the section detailing pharmacophore representation, it is mentioned that $\\mathcal{L}$ comprises only pharmacophoric descriptors. However, distances are subsequently incorporated into the representation. The notation and description should be revised for consistency and clarity.\n\n**Q2. Model Input Specification:** Within the model input section, node labels are currently denoted as $V_p$, which is a set of node. It may be more appropriate to represent these as node element to enhance the clarity.\n\n**Q3. Negative Data Augmentation Strategy:**\nThe current approach limits displacement directions to a single direction to avoid cancellation effects. Allowing for all displacement directions except for the case of cancellation can highly increase the diversity of negative training data, might lead to better model performance. Are there specific reasons why the authors chose to use only a single direction for negative displacements?\n\n## **Results:**\n\n**Q4. Evaluation on Recent Datasets:** To better assess the method’s applicability to real-world scenarios, evaluation on more recent datasets like LIT-PCBA[1] is recommended. Additionally, utilizing metrics beyond AUROC, such as enrichment factors (EF) (similar to BEDROC for early recognition of hit candidates), could provide a better understanding of the model’s performance in virtual screening tasks.\n\n**Q5. Benchmarking with other methods:** Since one of the main contributions of this paper is \"*fast virtual screening in the embedding space and evaluate the performance of our method through experiments on virtual screening benchmark datasets*\", the authors should benchmark their virtual screening performance with other models. It seems that there are no significant differences in the objectives or methods of the compared models relative to the current work. Are there reasons why the authors did not benchmarked their model with the previous works such as DrugClip or PharmacoNet?\n\n**Q5-1. Ambiguous criteria for showing virtual screening performance:**  The authors compared their performance with CDPKit's alignment algorithm. It only make sense that the propose method does well on virtual screening only if CDPKit alighment algorithm's performance is already good enough. Authors may provide the performance of CDPKit alignment algorithm for virtual screening explicitly in their manuscript.\n\n**Q6: Limitation in Protein-Specific Virtual Screening Capability:**\nThe methodology seems similar to ligand-based approaches, where ligand structures are pre-generated and graph matching determines activity. This may limit the model's applicability, since considering only ligand structure of protein-ligand complex cannot fully consider protein pocket. For example, if protein has large binding site that various ligands with different pharmacophoric sites can interact with, considering only a single ligand might give a bias during virtual screening on the protein target. Why did authors adopted ligand-based approach instead of protein-based one, such as using protein binding site's pharmacophore instead its binding ligand?\n\n## **References:**\n[1] Tran-Nguyen, Viet-Khoa, Célien Jacquemard, and Didier Rognan. \"LIT-PCBA: an unbiased data set for machine learning and virtual screening.\" Journal of chemical information and modeling 60.9 (2020): 4263-4273.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 2, "presentation": 2, "contribution": 1, "strengths": "**S1. Alignment of Contribution and Results:**\n\tThe study presents a coherent alignment between its contributions and the resulting outcomes. The objectives set forth by the authors are consistently addressed throughout the work.\n\n**S2. Effective Representation Learning via contrastive learning:**\n\tThe approach to representation learning through a contrastive learning with data augmentation appears to function as intended. The learned representations for pharmacophores are well-clustered in embedding space.", "weaknesses": "**W1. Limited Methodological Novelty:**\n\tThe proposed methods do not introduce significant novel approaches. As mentioned in the manuscript, the concepts of Neural Subgraph Matching or neural network architectures are already proposed by other previous works. While the formulation applied to pharmacophores—particularly the augmentation strategies for contrastive learning—is noteworthy, it may not sufficiently advance the general methodologies handled in the ICLR community. If there's any other novelty compared to previous works, the points should be clearly explained in the manuscript. Authors may propose novel input processing schemes or model architectures better suited to the pharmacophore graph, or improve neural subgraph matching techniques.\n\n**W2. Insufficient Experimental Results:**\n\n- **W2.1 Lack of Comprehensive Benchmarks:**\n\tFor a study proposing algorithms intended for virtual screening, it is essential to benchmark against established methods such as DrugClip or PharmacoNet to demonstrate comparative advantages. The absence of such comparisons, without a compelling justification, weakens the evaluation of the proposed method’s efficacy. Additionally, reliance on the outdated DUD-E benchmark and the arbitrary selection of only 10 targets out of 102 weakens the robustness of the experimental validation. (see Q4, 5)\n\t\n- **W2.2 Ambiguous goal of benchmark experiment:**\n\tThe goal of \"*to achieve comparable values between our model and the alignment algorithm*\" would be meaningful only if the alignment method inherently guarantees superior results compared to the previous methods, which is not adequately demonstrated in the current manuscript. (see Q5-1)\n\t\n- **W2.3 Suitability for Virtual Screening:**\n\tThe method appears to focus on ligand interactions with only parts of the protein pharmacophore, potentially neglecting the comprehensive information of the global protein binding site. This partial consideration might introduce bias, and it remains unclear whether the method can outperform existing techniques that utilize complete protein-ligand information. Empirical evidence showing superior performance in this regard would strengthen the study. (see Q6)", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "HgicfIPakk", "reviewer": "Reviewer_2qzD", "summary": "In this paper, the authors propose a contrastive learning approach for pharmacophore screening, based on subgraph matching. The key idea is to employ approximate subgraph matching for querying conformational database, a main step in pharmacophore screening. The subgraph matching is done through a contrastive learning approach by encoding query-target relationships in the embedding space. Their model has been validated based on benchmark dataset including DUD-E.", "questions": "1) From Methodology point of view, is there any major novelty in the current paper?\n2) Other than the efficiency in terms of the runtime, any other clear advantage of the current model?", "limitations": "", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "It is interesting to reinterpret pharmacophore screening problem as an approximate subgraph matching problem.", "weaknesses": "From the methodology point of view, the paper lacks novelty. The contrastive learning framework and augmentation module are all rather standard approach in GNN models. The contribution is not significant. Further, the performance is not very impressive as shown in Table 1. Even though the authors have emphasized that \"our goal is to achieve comparable values between our model and the alignment algorithm\", the only advantage of the current model seems to be the runtime efficiency.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (Poster)