A General Transformer-Based Multi-Task Learning Framework for Predicting Interaction Types between Enzyme and Small Molecule

Predicting enzyme–small molecule interactions is critical for drug discovery and generally understanding the biochemical processes of life. While recent deep learning approaches have shown promising results, several challenges remain: the lack of a comprehensive training dataset, architectures lacking communication between representations of enzymes and small molecules, the tendency to sim- plify the problem as enzyme–substrate vs. enzyme—non-interacting, and thereby, misclassify enzyme–inhibitor pairs as substrates, and the negligence of the true impact of data leakage on the model’s performance. To address these issues, we present EMMA (Enzyme–small Molecule interaction Multi-head Attention), a transformer-based multi-task learning framework designed to learn pairwise interaction signals between enzymes and small molecules, that can well general- ize to out-of-distribution data. EMMA operates directly on the SMILES string representation of small molecules and enzyme sequences, with two classifica- tion heads that distinguish enzyme—non-interacting, enzyme—substrate, and enzyme—inhibitor pairs. By evaluating EMMA under five distinct data-splitting regimes that control for different types of data leakage, we demonstrate that EMMA achieves a strong and robust performance, particularly for previously unseen combinations of enzymes and small molecules. Further, a deeper anal- ysis highlights that the topological properties of the enzyme—small molecule interaction network are crucial for the model performance and its ability to gen- eralize, yet again stressing the decisive role of well-designed training datasets for successful model training. 1 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint 1 Main Enzymes are proteins that catalyze essential biological reactions, playing key roles in metabolism, energy transfer, and signal transduction [1]. Understanding enzyme— small molecule interactions is fundamental in biology and drug design. These interac- tions can be categorized into three types: (i) enzyme—substrate, when a small molecule binds to the enzyme’s active site and is converted into a product, (ii) enzyme—non- interacting, where a small molecule with very low affinity or incompatibility for the enzyme cannot bind to the active site, and (iii) enzyme—regulator pairs, where small molecules such as inhibitors or activators modulate enzymatic activity [1]. The functions of many enzymes are still poorly understood, as only about 3% of human proteins have high-quality functional annotations supported by experimental evidence [2]. Computational methods, particularly machine learning, have revolution- ized biological research by enabling large-scale data analysis and uncovering insights that were previously unattainable through manual methods [3]. However, predicting enzyme—small molecule interactions remains a challenge due to their larger complex- ity compared even to drug—target interactions. This complexity stems from the high specificity of enzymes, dynamic and transient interactions with ligands, and conforma- tional flexibility during catalysis [4–6]. Additionally, enzymes are often multifunctional and promiscuous binders [7, 8]. The evolutionary diversity of enzymes and the lack of high-quality data exacerbate these challenges [5, 6], making the design of efficient computational methods to predict enzyme functionality even more challenging. Computational methods for predicting interactions of enzymes and small molecules can be docking-based [9–12] and non-docking [13–17] (ligand-based or data-driven), each with its own trade-offs in accuracy, interpretability, and scalability. Since docking- based methods are prohibitively slow for large datasets [18], this work addresses non-docking prediction, which is typically done by means of machine learning , or more precisely, nowadays deep learning. In such non-docking methods, training a model to predict all interaction types (e.g., substrates, regulators, and non-interactors) remains challenging due to the lack of comprehensive and curated enzyme—small molecule datasets. Consequently, many studies simplify the task to a binary classification problem, distinguishing known enzyme—substrate pairs from either experimentally validated [13] or synthetically generated [14–17,19] non-substrates. However, the term non-substrate can refer to both small molecules that fail to bind an enzyme due to low affinity or structural incompatibility, as well as to tightly binding regulating ligands (inhibitors or activators). We use the term non-interacting small molecule to avoid this confusion. Since experimental data on non-interacting ligands is scarce, synthetically gen- erated enzyme—non-interacting pairs are often used to train computational models, which may introduce biases and data shifts caused by the generation process. Fur- thermore, simplification to a binary classification task introduces a critical limitation: the model is trained to differentiate high-affinity enzyme—substrate interactions from presumed enzyme—non-interacting pairs, and thus is prone to misclassifying other high-affinity binders such as inhibitors or cofactors as substrates, leading to incorrect predictions. Furthermore, many datasets are skewed towards interactions with energy- transfer small molecules (e.g., the dataset from [14] contains 1,379 small molecules 2 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint in total and ∼28% of interactions involve an energy-transfer small molecule). This imbalance makes these interactions easily predictable, often inflating model perfor- mance metrics, but they are less informative for truly learning interaction patterns with diverse substrate molecules. Another cause of artificially overinflated performance in models for enzyme— small molecule interaction predictions that has to be accounted for is data leakage, in particular inter-sample similarity leakage that occurs when similar data points appear in both training and test sets [20]. To minimize the impact of data leakage in biological datasets and ensure robust model evaluation, various tools have been developed, including DeepChem [3], LoHi-Splitter [21], GraphPart [22], and DataSAIL [23]. Among these tools, DataSAIL is the only one capable of minimizing inter-sample similarity leakage in two-dimensional datasets, such as the enzyme—small molecule interactions. In this study, we present a comprehensive framework for enzyme—small molecule interaction prediction that addresses these challenges through three key contribu- tions: (i) we constructed the novel EMMI (Enzyme—sMall Molecule Interaction) dataset encompassing experimentally confirmed enzyme—substrate, enzyme—non- interacting, and enzyme—inhibitor pairs, significantly expanding beyond previous binary formulations; (ii) we implemented rigorous data splitting protocols using Data- SAIL to control inter-sample similarity leakage between training and test sets; and (iii) we developed a transformer-based multi-task learning architecture that processes both enzyme and small molecule features. The model employs two classification heads: an interaction head, which distinguishes between low- and high-affinity interactions, and a subclass head, which specifies the interaction type: either inhibitor or substrate. Given the lack of comparable models, we benchmarked EMMA against random forests and demonstrated its effectiveness in distinguishing between enzyme—substrate, enzyme— non-interacting, and enzyme—inhibitor pairs, while maintaining robustness against data leakage. 2 Results 2.1 The EMMI dataset We created the EMMI (Enzyme-sMall Molecule Interactions) dataset comprising a total of 147,224 samples, divided into three types: enzyme—non-interacting, enzyme— substrate, and enzyme—inhibitor pairs. The dataset comprises 15,652 unique enzymes and 28,977 unique small molecules, resulting in an overall enzyme-to-small molecule ratio of 0.54 (Table 1). The dataset is well-balanced in size across all enzyme—small molecule relation types, while still capturing diversity in both enzyme and small molecule entries. We intentionally kept the number of enzyme—non-interacting pairs (79,668 pairs) higher than the combined number of enzyme—substrate and enzyme— inhibitor pairs (67,556 pairs) to better reflect the biological setting, where most enzymes do not interact with the majority of small molecules. The ratio of enzyme— substrate to enzyme—inhibitor samples is ∼1. The higher enzyme-to-small molecule ratio in the enzyme—substrate subset highlights its enzyme-centric nature, whereas the enzyme—non-interacting and —inhibitor subsets are more small molecule-rich. 3 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint T able 1 EMMI dataset statistics and composition after and (before) clustering-based down-sampling. Enz/mol: enzyme /molecule; enz: enzyme; NI: non-interacting; sub: substrate; inh: inhibitor Dataset samples Unique entries Enz/mol ratio Enzyme ID molecule ID Enz–NI 79,668 (319,413) 4,351 (4,351) 14,283 (164,317) 0.305 (0.027) Enz–sub 34,221 (34,221) 11,343 (11,343) 6,569 (6,569) 1.730 (1.730) Enz–inh 33,335 (807,173) 5,997 (5,997) 10,436 (564,613) 0.575 (0.011) Subclass 67,556 (841,394) 17,340 (17,340) 17,005 (571,182) 1.020 (0.030) Interaction 147,224 (1,161,365) 15,652 (15,652) 28,977 (691,959) 0.540 (0.023) A single enzyme and a single small molecule can be a member in the pairs with different labels (substrate, inhibitor, or non-interacting), but the distribution of these labels among enzymes and small molecules is markedly different in the EMMI dataset (Figure 1). The enzyme distribution reveals substantial overlap between substrate-binding and inhibitor-binding enzymes: while 8,104 enzymes are uniquely substrate-binding, 1,630 are uniquely inhibitor-binding, and 1,145 are uniquely non- interacting, there are 4,773 multi-label enzymes with 1,266 of them having three labels, reflecting the multifunctional nature of enzymes [24, 25]. In contrast, the distribution of small molecules shows far less overlap: with 1,966 multi-label small molecules and only 349 small molecules being simultaneously labeled as substrates, inhibitors, and non-interacting. Most small molecules are strongly partitioned into either substrates (5,586), inhibitors (8,852), or non-interacting (12,518), consistent with their specialized roles towards enzymes. Collectively, these results highlight an important asymmetry; enzymes are often multifunctional and promiscuous binders [7], while small molecules tend to exhibit more exclusive binding roles due to their intrinsic chemical specificity and the selective optimization processes of drug discovery, which deliberately favor narrow, well-defined interactions to maximize efficacy and minimize off-target effects [26–28]. (Figure 1). 2.2 Types of inter-sample similarity leakage We propose that inter-sample similarity leakage can be further partitioned into two subtypes, each affecting the model’s performance differently. Before defining these subtypes, it is helpful to represent the enzyme–small molecule interaction network as a bipartite graph, where enzymes and small molecules form two disjoint node (vertex) sets, and edges between node sets represent observed interactions. When this network is split into training and test sets based on either node type (enzyme or small molecule), inter-sample similarity within the vertex sets can cause information leakage in two forms. The first form is single-label inter-sample similarity leakage (SLSL). This happens when similar data points—e.g., similar enzymes (or their close variants) interact with multiple small molecules—appear in both the training and test sets with the same label (Figure 2 a). This form of leakage can artificially boost a model’s performance: after a few training iterations, the model may effectively memorize the features of the bridging enzymes, leveraging their repeated presence across multiple substrates rather 4 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint 1145 8104 406 1630 1534 1567 1266 Non-interacting enzyme Substrate-binding enzyme Inhibitor-binding enzyme MBC: 3.0 MD: 2.0 MBC: 6.7 MD: 2.0 MBC: 11.4 MD: 6.0 MBC: 1.4 MD: 2.0 MBC: 12.1 MD: 14.0 MBC: 11.3 MD: 5.0 MBC: 13.0 MD: 17.0 a 12518 5586 398 8852 983 236 349 Non-interacting small molecule Substrate Inhibitor MBC: 0.0 MD: 1.0 MBC: 0.0 MD: 1.0 MBC: 11.7 MD: 7.0 MBC: 0.0 MD: 2.0 MBC: 10.7 MD: 7.5 MBC: 11.2 MD: 5.0 MBC: 13.7 MD: 26.0 b Fig. 1 Distribution of class labels in the EMMI dataset for enzymes (a) and small molecules (b). MBC: median betweenness centrality. MD: median degree. Note: Both MBC (log+1) and MD are log-scaled. than learning genuine interaction patterns—in other words, relying on correlation through repetition rather than true causality. Theoretically, single-label nodes with high connectivity—in terms of degree or betweenness centrality (BC)—that bridge the training and test sets are more likely to mislead the model toward memorization (Figure 2 a vs. c), due to the dominant availability of their features during training (Figure 2 e vs. g). The impact of SLSL can be more severe in tree-based models, such as random forests and XGBoost, where feature importance is often ranked based on how frequently features are used to create decision splits. The second form is multi-label inter-sample similarity leakage (MLSL). This occurs when the same node or very similar nodes are present in both the training and test sets but are associated with different labels (Figure 2 b and d). This form of leakage can negatively affect a model’s performance, particularly if the model tends to memorize rather than learn interaction signals. However, these multi-label nodes are more likely to lead the model toward learning interaction signals if they have low to moderate connectivity. In contrast, when these nodes are highly connected, the model is more likely to memorize the shared node across many interactions, because it is present in many sample pairs during the training process. If this node (or its close variants) is labeled differently in the test set, the model is likely to fail. Unlike single-label nodes, multi-label nodes tend to have a higher median BC and degree (Figure 1). Overall, we believe that the tendency of the model to memorize or learn the inter- action signals is not only related to the magnitude and the form of the inter-sample similarity leakage, but also the topology of the enzyme—small molecule interaction network in terms of connectivity of nodes in the training and test sets (Figure 2 b vs. d and f vs. h); single-label nodes with high connectivity are more prone to being mem- orized by the model due to the dominant availability of their features during training, whereas multi-label nodes with moderate connectivity can trigger the model to learn interaction signals. 5 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint 8 Substrate Non-interactingInhibitorEnzyme 1 2Memorization enzyme 1 1 2 a e 4 5 3 1 3 1 4 1 5 Train Test 8 1 2 1 2 b f 4 5 3 1 3 1 4 1 5 Train Test Learning interaction signals 8 1 2 1 2 d h 4 5 3 1 3 Train Test Learning interaction signals 8 1 2Memorization enzyme 1 1 2 c g 4 5 3 1 3 Train Test 1 5 1 5 Fig. 2 Different forms of inter-sample similarity leakage across training and test sets. (a, c)Single- label inter-sample similarity leakage (SLSL) where single-label enzyme with high and low connectivity (degree and betweenness centrality) shared between training and test sets. (b, d) Multi-label inter- sample similarity leakage (MLSL) where multi-label enzyme with high and low connectivity shared between training and test sets. (e, f, g, h) Frequent occurrence of enzyme 1 features under high and low connectivity and its impact on model learning ability. T able 2 Overview of different splits with corresponding sample and label ratios. Enz–Inh: Enzyme–inhibitor. Enz–Sub: Enzyme–substrate. Enz–NI: Enzyme–non-interacting. Num. Samples Enz-Inh/Enz-Sub Enz-Sub/Enz-NI Enz-Inh/Enz-NI Split Training Test Training Test Training Test Training Test Random 117,779 29,445 0.974 0.974 0.430 0.430 0.418 0.418 Enzyme-based 117,771 29,453 0.985 0.932 0.424 0.454 0.417 0.423 Small molecule-based 117,863 29,361 1.033 0.786 0.403 0.547 0.416 0.430 Label-based 117,594 29,630 0.973 0.977 0.433 0.417 0.421 0.408 Two-dimensional 64,881 16,083 1.193 0.763 0.324 0.666 0.387 0.508 2.3 Data split regimes To assess model performance, we employed four one-dimensional split regimes (random, enzyme-based, small molecule-based, and label-based) as well as a two- dimensional split regime (Section 4.1.5, Table 2), ensuring an 80:20 % ratio between training and test sets in each case. For all regimes, we kept the ratio of enzyme— substrate to enzyme—inhibitor pairs close to 1 and the ratio of enzyme—interacting (enzyme—substrate or enzyme—inhibitor) to enzyme—non-interacting pairs ∼0.43. Each one-dimensional split produces training and test sets of ∼118K and ∼29K, respectively, retaining all the samples from the EMMI dataset, and the two- dimensional split yields a smaller total set due to its stricter partitioning criteria: ∼65K in the training and ∼16K in the test sets. 6 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint T able 3 Number interactions with a shared molecule/enzyme between training and test sets across split methods. MLSL: Multi-label inter-sample similarity leakage. SLSL: Single-label inter-sample similarity leakage. With shared enzymes With shared molecules Split Method MLSL SLSL MLSL SLSL Total interactions Random 93,097 18,013 28,680 72,878 212,668 Enzyme-based 0 0 26,943 66,927 93,870 Small molecule-based 76,922 13,985 0 19 90,926 Label-based 3,424 29,128 21,998 56,572 111,122 Two-dimensional 0 0 0 0 0 Each split method exhibits distinct characteristics with respect to inter-sample similarity leakage (see Table 3 and Section 4.1.5 for details). In the random split, both single- and multi-label nodes are distributed arbitrarily across train and test sets, without explicit control over MLSL or SLSL. Consequently, the random split yields the highest number of interactions involving either shared enzymes or shared small molecules between the training and test sets. If shared substructures (e.g., a benzene ring within a family of molecules) were also considered, this number increases even further; however, we report only the number of exact shared interaction partners. In the enzyme-based split, the size of SLSL interactions (66,927) is twice that of MLSL interactions (26,943) for small molecules. In contrast, in the small molecule-based split, the size of SLSL interactions (13,985) is approximately seven times smaller than MLSL interactions (76,922) for enzymes. The label-based split is designed to minimize the similarity leakage within each label set (see Section 4.1.5 for details). In this split, the MLSL is smaller than SLSL for both enzymes and small molecules; however, the number of interactions involving shared enzymes or small molecules is reduced by half compared to the random split. In contrast, in the two-dimensional split, there is no MLSL and SLSL at either level between the training and test sets. These observations on inter-sample similarity leakage in different splits can be explained by an analysis of betweenness centrality (BC) in the corresponding enzyme- small molecule networks (Figure 3). BC in the training sets shows a roughly similar distribution across all split regimes. In contrast, BC has different distributions across all split regimes in the test sets (Figure 3). Most small molecules (at least 50%) in the training and test sets have a low BC and thus are structurally unimportant in terms of network connectivity, as the first quartile is equal to the median. The distribution of BC in the upper quartile changes with the split method for both enzymes and small molecules in the test sets. Enzymes generally have higher BC values than small molecules in all training sets, and the same is true for the test sets; however, the distributions are narrower, suggesting fewer extreme central nodes. Across the test sets, the median BC of small molecules remains constant, while the median BC of enzymes varies for different split methods (see additional analysis using node degree for each split method in the supplementary material figure 1). 7 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint 0 5 10 15 20 0.0 0.1 0.2 0.3 0.4 Density a Random - Train Type Enzyme Small molecule 0 5 10 15 0.0 0.1 0.2 0.3 0.4 b Random - Test Type Enzyme Small molecule Enzyme Small molecule 0 5 10 15 BC score (log + 1 scaled) c Random - Train Enzyme Small molecule 0 5 10 15 d Random - Test 0 5 10 15 20 0.0 0.1 0.2 0.3 0.4 Density e Enzyme-Based - Train 0 5 10 15 20 0.0 0.1 0.2 0.3 0.4 f Enzyme-Based - Test Enzyme Small molecule 0 5 10 15 BC score (log + 1 scaled) g Enzyme-Based - Train Enzyme Small molecule 0 5 10 15 h Enzyme-Based - Test 0 5 10 15 20 0.0 0.1 0.2 0.3 0.4 Density i Small Molecule-Based - Train 0 5 10 15 0.0 0.1 0.2 0.3 0.4 j Small Molecule-Based - Test Enzyme Small molecule 0 5 10 15 BC score (log + 1 scaled) k Small Molecule-Based - Train Enzyme Small molecule 0 5 10 15 l Small Molecule-Based - Test 0 5 10 15 20 0.0 0.1 0.2 0.3 0.4 Density m Label-Based - Train 0 5 10 15 0.0 0.1 0.2 0.3 0.4 n Label-Based - Test Enzyme Small molecule 0 5 10 15 BC score (log + 1 scaled) o Label-Based - Train Enzyme Small molecule 0 5 10 15 p Label-Based - Test 0 5 10 15 20 BC score (log + 1 scaled) 0.0 0.1 0.2 0.3 0.4 Density q Two-Dimensional - Train 0 5 10 15 BC score (log + 1 scaled) 0.0 0.1 0.2 0.3 0.4 r Two-Dimensional - Test Enzyme Small molecule 0 5 10 15 BC score (log + 1 scaled) s Two-Dimensional - Train Enzyme Small molecule 0 5 10 15 t Two-Dimensional - Test Fig. 3 Betweenness centrality distributions for each split method. The red line indicates the median, and the whiskers represent the range of non-outlier values in the box plot 2.4 Model performance The EMMA model is a dual-stream transformer with separate streams for enzyme and small molecule representations. It respects the classical transformer block structure 8 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint (attention → residual connection and normalization → feed-forward → residual con- nection and normalization), but replaces the standard self-attention with a composite of self- and cross-attention mechanisms within each stream. The two streams are con- nected through a bidirectional cross-attention module, where each stream attends to the other (enzymes attend to molecules and molecules attend to enzymes) to capture mutual interaction signals. This dual-stream transformer is followed by two task- specific prediction heads: one for predicting the existence of an interaction (interaction task) and one for identifying the type of interaction (substrate or inhibitor; subclass task), both trained with focal loss functions (see Figure 7 and Section 4.3 for further details). 2.4.1 Benchmarking EMMA against a random forest baseline We trained EMMA and a random forest (RF) baseline on the five splits described above. The performance of the models varies depending on whether similarity leakage is controlled at the enzyme level, the small-molecule level, or both, or not at all, as in the random split (Table 4). This is in accordance with our earlier work, where we demonstrated that the amount of data leakage is directly correlated with performance metrics such as AUROC and accuracy ([23]). Interestingly, the RF model performs better than EMMA for a random split for both the interaction and subclass tasks, with an area under the receiver operating characteristic curve (AUROC) scores above 0.97 for the interaction task (IT) and subclass task (ST) (Table 4). However, random splits are known to inflate performance estimates by allowing highly similar molecules or enzymes to appear in both training and test sets, causing data leakage [23]. Consistent with this, we observed the highest degree of similarity leakage under the random split (Figure 4). In the enzyme-based split, the RF retains a modest advantage for interaction pre- diction, but EMMA surpasses it in the subclass task, achieving the highest AUROC (0.991). This indicates that EMMA captures more generalizable patterns when unseen enzymes are introduced. The gap becomes more pronounced in the small molecule- based split, where the RF performance degrades substantially (AUROC 0.734 for IT, AUROC 0.889 for ST), while EMMA maintains high scores (AUROC: 0.854 for IT, AUROC: 0.947 for ST). This demonstrates that EMMA generalizes more efficiently when challenged with unseen chemical scaffolds. The label-based and two-dimensional splits represent the most stringent evalu- ation settings. In the label-based split, the RF performance further drops in both tasks (AUROC: 0.836 IT, AUROC: 0.849 for ST), whereas EMMA remains robust (AUROC: 0.845 for IT, AUROC: 0.950 for ST). A similar pattern is observed for the two-dimensional split, where EMMA consistently outperforms RF across both tasks. These results highlight EMMA’s ability to capture transferable representations that support generalization beyond the training distribution. These results further support our assumption that for predicting enzyme—ligand interactions, tree-based models are more sensitive to sample-similarity leakage than neural network-based models. This is due to the inherent imbalance in network node connectivity and the tendency of tree-based models to prioritize features associated with highly connected nodes. 9 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint T able 4 Comparison of EMMA and the Random Forest (RF) across different data split methods. We report performance on the interaction task (IT) and subclass task (ST). Results are averaged over five runs with different random seeds, with standard deviation shown. The best-performing model for each split and task is highlighted in bold. Note: We trained two separate random forest classifiers, one for each task. AUROC: Area under the receiver operating characteristic curve Split Method Model Task AUROC (↑) Accuracy (↑) F1 score ( ↑) Random RF IT 0.976 ± 0.000 0.915 ± 0.001 0.910 ± 0.001 ST 0.989 ± 0.000 0.957 ± 0.001 0.958 ± 0.001 EMMA IT 0.944 ± 0.001 0.859 ± 0.001 0.845 ± 0.001 ST 0.984 ± 0.000 0.947 ± 0.001 0.948 ± 0.002 Enzyme-based RF IT 0.954 ± 0.000 0.878 ± 0.002 0.872 ± 0.002 ST 0.980 ± 0.000 0.937 ± 0.001 0.940 ± 0.001 EMMA IT 0.917 ± 0.002 0.821 ± 0.003 0.810 ± 0.004 ST 0.991 ± 0.001 0.934 ± 0.005 0.940 ± 0.005 Small molecule-based RF IT 0.734 ± 0.017 0.633 ± 0.016 0.681 ± 0.012 ST 0.889 ± 0.003 0.796 ± 0.011 0.814 ± 0.011 EMMA IT 0.854 ± 0.005 0.755 ± 0.010 0.772 ± 0.004 ST 0.947 ± 0.006 0.846 ± 0.011 0.878 ± 0.007 Label-based RF IT 0.836 ± 0.007 0.749 ± 0.005 0.714 ± 0.004 ST 0.849 ± 0.005 0.591 ± 0.006 0.363 ± 0.015 EMMA IT 0.845 ± 0.006 0.753 ± 0.006 0.726 ± 0.008 ST 0.950 ± 0.009 0.881 ± 0.040 0.876 ± 0.047 Two-dimensional RF IT 0.727 ± 0.012 0.648 ± 0.014 0.704 ± 0.009 ST 0.908 ± 0.003 0.769 ± 0.005 0.765 ± 0.005 EMMA IT 0.847 ± 0.008 0.760 ± 0.008 0.768 ± 0.017 ST 0.954 ± 0.004 0.884 ± 0.011 0.899 ± 0.012 Such features frequently appear in decision splits that maximize the reduction of the splitting criterion, leading to higher feature importance. 2.4.2 Effect of split regime characteristics on EMMAs performance Our results highlight the importance of different forms of inter-sample similarity leakage in enzyme—small molecule interaction prediction. It seems that learning trans- ferable representations for unseen small molecules is considerably more challenging than for unseen enzymes. However, measuring data leakage (Equation (3)) reveals that the total similarity leakage (TSL) in the small molecule-based split is higher than for the enzyme-based split (0.32 vs. 0.479, Figure 4). An explanation for these observations can be found in the ratio of enzymes to small molecules in the EMMI dataset (Table 1), which is 0.54, creating inherent challenges for splitting the data. As mentioned above, the enzyme—small molecules interaction datasets can be represented as a bipartite graph (Section 4.1.6). When applying small molecule-based splitting to this structure, the high connectivity of enzymes across small molecule partners frequently creates bridges between training and test partitions, as many enzymes interact with multi- ple small molecules. Their inclusion in training and test sets through different small 10 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint Random Small molecule-based Label-based Enzyme-based Two-dimensional Split Method 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 Leakage Score0.320 0.279 0.309 0.318 0.270 0.239 0.200 0.147 0.011 0.014 0.559 0.479 0.456 0.329 0.284 a MSL ESL TSL 0.0 0.1 0.2 0.3 0.4 0.5 Leakage Score 0.84 0.86 0.88 0.90 0.92 0.94 AUROC MSL r = 0.746 ESL r = 0.209 TSL r = 0.348 b 0.0 0.1 0.2 0.3 0.4 0.5 Leakage Score 0.95 0.96 0.97 0.98 0.99 MSL r = 0.743 ESL r = -0.144 TSL r = 0.018 c Fig. 4 (a) Train-test similarity-based data leakage for each splitting method. ( b) Pearson correla- tion between the AUROC score and data leakage for the interaction prediction head. ( c) Pearson correlation between the AUROC score and data leakage for the subclass head. MSL: small molecule- similarity leakage, ESL: enzyme-similarity leakage. molecule partners introduces data leakage, which cannot be avoided when the split- ting procedure is based on small molecules. This architectural constraint explains why small molecule-based splits consistently show higher data leakage than enzyme-based splits in our experiments. Despite the high leakage observed in the small molecule–based split, why does the model trained on the enzyme-based split achieve higher performance? The key

of the small molecule–based split is that it allows inter-sample similarity leakage at the enzyme level. Additionally, because the number of multi-label enzymes is larger than the number of multi-label small molecules (Figure 1), the size of the MLSL is greater than the SLSL for enzymes in the small molecule–based split compared to the enzyme-based split. As discussed in section 2.2, MLSL can negatively affect model performance when the nodes are highly connected (e.g., have a high betweenness centrality (BC)). The model performance decreases for nodes with very low (bins 1–3) or very high (bins 7–10) BC (as shown Figure 5), particularly for those splits designed to minimize similarity leakage, while the highest performance is observed for nodes with moderate BC (bins 4–6) across all splits. Our results confirm that as the 11 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint 1 2 3 4 5 6 7 8 9 10 Binned betweenness 0.4 0.5 0.6 0.7 0.8 0.9 1.0AUROC a) Interaction head - Enzyme Betweenness Random Enzyme-based Small molecule-based Label-based Two-dimensional r 0.344 -0.172 -0.435 -0.020 -0.235 1 2 3 4 5 6 7 8 9 10 Binned betweenness 0.4 0.5 0.6 0.7 0.8 0.9 1.0AUROC b) Interaction head - Molecule Betweenness Random Enzyme-based Small molecule-based Label-based Two-dimensional r 0.838 0.599 0.148 0.789 0.143 1 2 3 4 5 6 7 8 9 10 Binned betweenness 0.86 0.88 0.90 0.92 0.94 0.96 0.98 1.00 1.02AUROC c) Subclass head - Enzyme Betweenness Random Enzyme-based Small molecule-based Label-based Two-dimensional r 0.382 0.669 -0.526 0.474 -0.317 1 2 3 4 5 6 7 8 9 10 Binned betweenness 0.86 0.88 0.90 0.92 0.94 0.96 0.98 1.00 1.02AUROC d) Subclass head - Molecule Betweenness Random Enzyme-based Small molecule-based Label-based Two-dimensional r -0.241 -0.167 -0.479 0.198 -0.708 Random Enzyme-based Small molecule-based Label-based Two-dimensional Fig. 5 Model performance (AUROC) across different levels of betweenness centrality. Nodes were grouped into bins according to their log-transformed betweenness centrality, and AUROC was com- puted within each bin for both enzyme- and molecule-level betweenness. Results are shown separately for the interaction and subclass heads across all data splits. Tables beneath each subplot report the Pearson correlation coefficients (r) between AUROC and betweenness bin indices for each split. BC of enzymes increases, model performance decreases. Since most highly connected enzymes tend to have multiple labels (Figure 1a), we conclude that the larger MLSL for enzymes in the small molecule–based split is responsible for its lower performance compared to the enzyme-based split. Interestingly, in the label-based split, where multi-label enzymes and small molecules can appear in both the training and test sets while minimizing similar- ity leakage within each label set, the model achieved lower performance compared to the enzyme- and small molecule-based splits, despite showing an intermediate TSL of 0.456. The size of MLSL is smaller than SLSL at both the enzyme and small molecule levels for this split (Table 3). Here, the enzyme–substrate set was divided using a small molecule-based split because its enzyme-to-small molecule ratio is greater than 1, indicating that small molecules are less connected on average; thus, split- ting by small molecules better minimizes inter-sample similarity leakage. Conversely, the enzyme–inhibitor and enzyme–non-interacting sets were divided using an enzyme- based split, as their ratios are less than 1 and enzymes are more closely connected, making enzyme-based splitting more effective at reducing leakage. By splitting the data into these three smaller problems rather than splitting the entire dataset at once, the leakage within each subset was effectively minimized. However, leakage across label sets still occurs in the form of MLSL at both enzyme and small molecule levels. 12 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint In the two-dimensional split, the most challenging evaluation scenario, which min- imizes similarity leakage across both enzymes and small molecules simultaneously, the model achieves a comparably low performance, reflecting the inherent difficulty of predicting interactions for entirely novel enzyme—small molecule combinations. In this split, the distributions of BC for both the training and test sets are narrower, suggesting fewer extreme central enzymes and small molecules compared to the one- dimensional splits. This stems from the fact that the two-dimensional split removes interactions that bridge between training molecules and test molecules, to minimize inter-sample similarity leakage. In our analysis, the number of multi-label enzymes decreased by 2% in the combined train–test data of the two-dimensional split com- pared to the EMMI dataset. In contrast, the number of multi-label small molecules showed only a slight decrease. Importantly, the ST consistently achieves better performance metrics than the IT across all splitting strategies. This difference may be rooted in how the dataset was constructed: interacting pairs (substrates and inhibitors) are supported by strong experimental evidence, either through direct biochemical evidence or experimentally validated functional GO annotations (see Section 4.1). By contrast, enzyme—non- interacting pairs are defined indirectly, via taking pairs with low binding affinity in experimental assays. As a result, the likelihood of mislabeling is higher for enzyme—non-interacting pairs than for interacting ones. Consequently, subclass-level classification benefits from more reliable ground-truth labels, whereas binary inter- action tasks suffer from noise introduced by uncertain non-interacting assignments. As systematic data on enzyme—non-interacting pairs is lacking, we believe that our strategy for generating negative data—incorporating experimentally validated enzyme—non-interacting ligand pairs—has borne fruit across different split regimes. This emphasizes the value of integrating experimental evidence and reducing reliance on synthetic negatives to ensure biologically meaningful predictions. Furthermore, the higher number of multi-label nodes in the interaction dataset, compared to the sub- class dataset (Table 1, makes the interaction task inherently more difficult than the subclass task. 3 Discussion In this study, we introduce EMMA, a Transformer-based multi-task learning frame- work that leverages a novel task formulation to jointly predict two distinct yet related binary tasks within a shared latent space using two classification heads: (i) an inter- action head that distinguishes low-affinity pairs from high-affinity ones, and (ii) a subclass head that further classifies high-affinity interactions as either inhibitory or catalytic. To train and validate the model, we compiled a comprehensive dataset that incorporates experimentally validated enzyme—small molecule pairs not only for sub- strate and inhibitor pairs, but also for non-interacting ones. This biologically grounded approach yielded exciting results under different split regimes, with EMMA achieving an average AUROC of 0.866 for the interaction task and 0.960 for the subclass task, demonstrating a high ability at rank prediction. 13 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint We specifically studied the effect of inter-sample similarity leakage on model per- formance and showed that higher leakage does not always improve performance. We further divided the leakage into single-label and multi-label inter-sample similarity leakage (SLSL and MLSL, respectively). We demonstrated how SLSL can artifi- cially enhance model performance, whereas MLSL is harmful for it. Furthermore, we demonstrated that the topology of the enzyme—small molecule interaction net- work, represented as a bipartite graph, and how it is separated into the training and test splits, could amplify or reduce the effect of the MLSL and SLSL. In particular, nodes with high connectivity are prone to be memorized, as when similarity leakage is minimized, performance on nodes with high connectivity—often corresponding to multi-label nodes—declines markedly. We believe that the model’s performance on the single- and multi-label nodes with different degrees of BC can provide insight into its ability to learn and gen- eralize. To the best of our knowledge, there is currently no established metric to evaluate this aspect of model capability. We believe our findings offer a foundation for developing such a metric. Our analysis emphasizes the complexity of predicting enzyme—substrate and enzyme—inhibitor interactions, where special attention needs to be paid to training regimes, interaction network topology, and data splitting to ensure models’ generalizability to unseen enzymes and small molecules. 4 Methods 4.1 Creating the EMMI dataset First, we have compiled a dedicated dataset of Enzyme-sMall Molecule Interactions (EMMI). The EMMI dataset contains experimentally confirmed enzyme—substrate, enzyme—non-interacting, and enzyme—inhibitor pairs. To create this dataset, we extracted data from the following sources: BindingDB [29], Brenda [30], ChEMBL [31, 32], Gene Ontology [33, 34], Iupher-PBS [35], PubChem [36], Rhea [37], Sabio-RK [38], and UniProt [39] (Figure 6). In the data-processing pipeline, first, all enzyme–small molecule pairs were extracted and each entity was mapped to a unique identifier. For enzymes, UniProt IDs were provided by the source databases. For small molecules, however, most of the sources provided only an IUPAC name or a common name. To map small molecule names to standardized identifiers, if the source database did not provide an identifier, we first attempted to retrieve PubChem Compound IDs (CIDs) using the PubChem PUG REST API [40]. For small molecule names that lacked a CID, we created a man- ually curated mapping. This dictionary integrates data from several databases (ChEBI ontology [41], HMDB [42], SMPDB [43], and KEGG [44]), and contains 610,608 unique small molecule names mapped to IDs. Our mapping strategy prioritizes a PubChem ID; if none is found, we then search for a ChEBI ID, followed by a KEGG ID. These identifiers were then used to associate enzymes with sequences and small molecules with SMILES strings. Only pairs with experimentally validated annotations were retained (see Section 4.1.1 for details). Pairs with invalid sequences or SMILES strings were excluded. Since the maximum position embeddings (also called token size limit) in ESM2 t30 [45] and MolFormer XL [46] are 1,024 and 202, respectively, we 14 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint EMMA DatasetExtractEnzyme-Substrate ExtractEnzyme-Inhibitor Extract Enzyme-Non-Interacting Mapping to Identifier M apping tosequence and SMILES ExcludeCofacotors and EnergyMolecules Experimentally Validation BindingDB Brenda ChEMBL UniProt PubChem Rhea Iupher-PBS Gen Ontology Sabio-RK Fig. 6 Abstract view of data processing steps to create EMMI dataset from various sources. only retained enzyme sequences of up to 1,500 amino acids (∼47% longer than the embedding limit of the ESM2 t30), with only ∼11% of enzymes in the EMMI dataset have a length between 1,025 and 1,500 amino acids. For SMILES strings, although the majority already fall well below the 202-token limit, we permitted SMILES strings of up to 512 characters; as a result, only less than 1% of small molecules exceed the token limit in the EMMI dataset. These thresholds were chosen to minimize the number of truncated enzymes or small molecules while maximizing the number of enzymes, thereby achieving a more balanced enzyme-to-molecule ratio within each label set (e.g., enzyme—substrate, enzyme—inhibitor, and enzyme—non-interacting) and in the EMMI dataset. SMILES strings were subsequently canonicalized using RdKit [47], and duplicates were removed at the level of enzyme–small molecule pairs. Finally, we filtered out all interactions involving known coenzymes, cofactors, and energy-transfer small molecules (Supplementary table 1). Details of the data mining and processing procedures for each database are provided in the Supplementary material Section S1. 4.1.1 Experimentally confirmed data points We distinguish between two types of enzyme—small molecule pairs. The first category includes enzyme—small molecule pairs supported by quantitative binding or functional assay measurements, referred to as binding affinity-related assay (BAA) pairs. For those interactions, one or multiple of the metrics in Table 5 have been measured. These assays provide direct evidence of molecular interaction. We refer to these pairs as BAA-confirmed pairs. The second category consists of enzyme–small molecule pairs inferred from exper- imentally validated Gene Ontology annotations (GOA). In this case, the presence of 15 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint T able 5 Descriptions of the experimental measurements for enzyme-small molecules interactions. Name Description Ki The inhibition constant is a direct measure of inhibitor binding affinity, independent of assay conditions. Kd The dissociation constant indicates ligand binding strength (e.g., from surface plasmon resonance (SPR) or isothermal titration calorimetry (ITC)). Ka The association constant is the inverse of K d, quantifying binding affinity. Km The Michaelis constant reflects enzyme-substrate affinity in catalytic reactions. IC50 The half-maximal inhibitory concentration is a functional measure of inhibition, depen- dent on assay setup (e.g., substrate concentration, incubation time). ECX The X-Maximal effective concentration is a functional measure of the concentration of a compound that produces X% of its maximal biological effect. Like IC50, it is dependent on assay conditions and reflects apparent potency rather than an intrinsic binding constant. Inhibition Reports the extent of enzyme activity reduction either as a percentage at a given concen- tration, or as the compound concentration required to achieve a defined level of inhibition. MIC The Minimum Inhibitory Concentration defines the lowest concentration of a compound required to prevent visible growth of a microorganism. a functional GO ID with associated experimental evidence code (EEC) for an enzyme implies functional relevance of the corresponding interacting small molecule. To map functional GO IDs to EECs, we processed the ∼1.2 billion GOAs and selected only those supported by EECs. This curated dataset was then used as a filter to retain only enzyme–small molecule pairs for which both the UniProt ID and the functional GO ID could be jointly mapped to an EEC-supported term. We refer to them as EEC-confirmed pairs. This classification framework allows us to capture both direct biochemical evidence and biologically supported associations. 4.1.2 Conversion of experimentally measured values The conversion of IC50 to Ki is necessary because Ki is a thermodynamic constant that reflects the true binding affinity of an inhibitor, independent of assay conditions such as substrate concentration or reaction time. In contrast, IC 50 values are assay-specific and cannot be directly compared across different studies without normalization. According to Kalliokoski et al. [48], although IC 50 values are more abundant in public datasets like ChEMBL, they are only moderately more variable than K i values (by approximately 25%). The authors showed that IC 50 values from mixed assays can be used reliably in large-scale analysis if converted appropriately, and that the average ratio between IC 50 and K i values across many protein–ligand pairs is approximately 2, which corresponds to assay conditions where the substarte concentration [ S] = Km in the Cheng-Prusoff equation [49]: Ki = IC50 1 + [S] K m (1) Thus, converting IC 50 to K i improves comparability of inhibitor-binding data in our integrated dataset. This heuristic approximation also avoids explicit dependence on substrate concentration or Km especially for working with large-scale heterogeneous 16 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint databases (e.g., PubChem, ChEMBL), where assay-specific details such as [S] and Km are often missing or unstandardized. While it is an approximation, it preserves the relative ranking of binding affinities across different inhibition levels and enables effective comparison of IC 50 and Ki values within a unified framework. We also included direct conversions of K a to Kd using: Kd = 1 Ka (2) For ECX assays, we did not find a standard or reliable method for converting values to Ki or Kd. Therefore, we retained the originally reported values. Furthermore, EC X values can reflect the effective concentration required to elicit a biological response, either from an inhibitor or a substrate. For large-scale analysis, it is not feasible to determine from the assay descriptions whether the reported EC 50 value is for an inhibitory or catalytic interaction. Therefore, we used them only to identify non- interacting pairs and considered only EC X values where 1 ≤ X ≤ 50. 4.1.3 Assigning active and inactive interaction For BAA-confirmed pairs, activity comments, when available, are used to label pairs as active or inactive based on the assay values. In PubChem, a general threshold of 10 µM is applied for BAAs (except for K m) where pairs with values above 10 µM are considered as inactive, while values below this threshold indicate interaction [50]. However, exceptions exist; in some cases, some pairs from functional assays with assay values above 10 µM are annotated as active, while others with values below 10 µM are labeled as inactive. These exceptions are assay-specific and comparable only under certain conditions. These exceptional cases make it risky to assign labels based on a universal thresh- old, especially for data points without activity comments. However, despite these

and the assay-specific nature of functional measurements, lower assay val- ues generally indicate a stronger binding affinity and the higher values indicate a weaker binding affinity for BAAs. Hence, we defined a number of binding strength categories as follows: strong (≤ 1 µM), moderate to weak (1 µM < x < 10 µM), and very weak ≥ 10 µM, corresponding to BAA, except K m. We only selected pairs with strong and very weak binding affinity for the EMMI dataset: pairs with assay values ≥ 10 µM are labeled as non-interacting, while pairs with values ≤ 1 µM are labeled as interacting. The conversion of IC 50 to Ki values could further effectively minimize the risk for IC50 assay above 10 µM which are annotated as active (Section 4.1.2). Furthermore, for functional assays with values below 10 µM that are labeled as inactive, the 1 µM cutoff is sufficiently low to minimize the number of such assays that could be incorrectly classified as interacting. The Km values of enzymes vary widely. For most enzymes, Km lies between 0.1 µM and 100,000 µM. The median K m is approximately 100 µM, with about 60% of all Km values falling in the range of 10–1000 µM [51, 52]. To define high-affinity and moderate interactions, we selected pairs with K m ≤ 10 µM. To mitigate potential 17 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint T able 6 Database filtering criteria with assay type specifications. Note: All IC50 values were converted to Ki. To indicate the source of each measurement, the original assay names from the respective databases were retained. Inh: inhibitor; Sub: substrate; NI: non-interacting; BAA: binding affinity-related assay; EEC: experimental evidence code; Act. Comt: activity comment; CAS: classified as a substrate; CAI: classified as an inhibitor; n.a.: not available. Database Binding Assay Type Threshold Evidence Act. Comt BindingDB Inh IC 50, Ki ≤ 1 µM BAA n.a. NI IC 50, Ki, Kd, EC50 ≥ 10 µM BAA n.a. ChEMBL Sub K m ≤ 10 µM BAA Active/n.a. Inh IC 50, Ki, MIC, Inhibition ≤ 1 µM BAA Active/n.a. Inhibition ≥ 70% BAA Active/n.a. NI IC 50, EC1-50, Ki, Kd ≥ 10 µM BAA Inctive/n.a. Inhibition ≤ 10% BAA Inctive/n.a. Km ≥ 100, 000 µM BAA Inctive/n.a. PubChem Inh IC 50, Ki ≤ 1 µM BAA Active NI IC 50, Ki, Kd, EC50 ≥ 10 µM BAA Inactive Sabio-RK Sub K m ≤ 1 µM BAA n.a. Inh IC 50, Ki ≤ 1 µM BAA n.a. NI IC 50, Ki ≥ 10 µM BAA n.a. Iuphar-BPS Inh IC 50, Ki ≤ 1 µM BAA n.a. Brenda Sub Km ≤ 10 µM BAA n.a. Inh IC 50, Ki ≤ 1 µM BAA n.a. NI IC 50, Ki ≥ 10 µM BAA n.a. Km ≥ 100, 000 µM BAA n.a. Sub - - EEC CAS Inh - - EEC CAI UniProt Sub - - EEC CAS Inh - - EEC CAI RHEA Sub - - EEC CAS GO.obo Sub - - EEC CAS noise, we define pairs with Km ≥ 100, 000 µM as non-interacting pairs. The thresholds for the MIC and inhibition assays were derived through an empirical assessment of the ChEMBL database, based on assays for which activity comments are available. Furthermore, within the interacting pairs, for those pairs with a K i, IC50, MIC or inhibition value or classified as inhibitor (CAI), we sub-labeled as enzyme—inhibitor pairs. Similarly, enzyme—small molecule pairs with a K m value or classified as a substrate (CAS) are sub-labeled as enzyme—substrate. 4.1.4 Clustering-based down-sampling After data labeling, to ensure balanced learning, we implemented clustering-based down-sampling for enzymes with a high number of associated inhibitors and non- interacting small molecules (Table 1). We first converted all small molecules to extended-connectivity fingerprints (ECFPs). Specifically, we computed Morgan finger- prints with a radius of 3 (ECFP6), a bit length of 2048, and chirality enabled using RdKit [47]. For each enzyme, the associated small molecules were clustered based on 18 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint their fingerprints using the K-means clustering algorithm, with a fixed random seed to ensure deterministic results. The global number of clusters K per enzyme was tuned to balance the classes (four clusters for inhibitors and six for non-interactors). From each cluster, we selected the small molecule closest to the cluster centroid. In cases where multiple small molecules were equidistant from the centroid, a deterministic tie-breaking procedure ensured reproducibility. This strategy allowed us to control label imbalance and retain structurally diverse representatives of both inhibitors and non-interacting small molecules for each enzyme. The composition of the EMMI dataset, including the contribution of each source database after cluster-based down-sampling, is summarized in Supplementary table 2. 4.1.5 Data splitting and data leakage We split the EMMI dataset into training (80%) and test (20%) sets using the fol- lowing methods implemented in DataSAIL [23]: enzyme-based, small molecule-based splits, and a two-dimensional split (Table 2). In the enzyme-based split, highly similar enzymes are confined to either the training or test set; and MLSL and SLSL can occur at the small molecule level (Table3). Likewise, in the small molecule-based split, struc- turally similar small molecules are restricted to a single set to prevent overlap; however, MLSL and SLSL can happen at the enzyme level. In the two-dimensional split, no similar enzyms and small molecules appear in more than one split, and consequently, no form of inter-sample similarity leakage (MLSL or SLSL) is permitted. We further propose a label-based split strategy, motivated by both the bipartite graph structure of the interaction datasets (Section 4.1.6) and a divide-and-conquer principle. In the bipartite graph, the enzyme-to-small molecule ratio acts as an indica- tor of which node type (enzymes or small molecules) has higher average connectivity or degree. By splitting the dataset based on the more connected partner, we minimize the number of bridges between the training and test sets. This approach effectively reduces inter-sample similarity leakage compared to splitting based on node type, which has a lower average connectivity. The divide-and-conquer principle complements this approach by treating each label set (e.g., enzyme—substrate) as a separate sub- problem. By performing splits independently within these smaller subproblems based on the enzyme-to-small molecule ratio, we can effectively control inter-sample simi- larity leakage than by splitting the entire dataset at once. In practice, in the EMMI dataset, since the enzyme-to-small molecule ratio is greater than 1 in the enzyme— substrate set, we split it using a small molecule-based split. Conversely, since this ratio is less than 1 in the enzyme—inhibitor and enzyme—non-interacting sets, we applied enzyme-based splits (Table 1). This process yielded three separate train/test sets, which we subsequently combined into one unified train set and one unified test set. Unlike the two-dimensional split, the label-based split allows MLSL and SLSL at both enzymes and small molecules (Table 3). However, by dividing the splitting problem into label-specific subproblems, this divide-and-conquer strategy minimizes similarity leakage more effectively within each label set than other one-dimensional splits applied to the entire dataset. 19 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint Finally, a random split was used as a baseline, in which interactions were assigned randomly to the training and test sets, allowing both MLSL and SLSL across enzymes and molecules. For all split regimes described above, 20% of the training data were further reserved as a validation set for hyperparameter tuning. We quantify the leakage between any two splits s and s′, e.g., training and test, of a dataset D as described in Joeres et al. [23] as Leakages,s′ = P u∈s,v∈s′ w(u, v)P {u,v}∈( D 2) w(u, v) . (3) u and v are data points in D, which can be enzymes, substrates, or interactions thereof. w : D × D → R is a similarity function between data points in D. The numer- ator measures the similarity between distinct splits, while the denominator normalizes by total similarity in the dataset. Higher leakage values indicate greater similarity between the corresponding splits. 4.1.6 Analysis of the bipartite interaction graph The enzyme—molecule interaction datasets can be formally represented as a bipartite graph where enzymes and molecules form two disjoint vertex sets, with edges repre- senting observed interactions. We computed the exact betweenness centrality (BC) for all nodes using the Brandes algorithm and applied a log(x + 1) transformation to stabilize variance and enable comparison [53]. It allows us to quantify the extent to which individual enzymes or molecules act as bridges within the interaction landscape and how deeply the network is interconnected (Figure 3). To further examine how nodes’ position in the interaction network of test sets affects the model performance, we conducted a binned betweenness analysis across all test sets. Nodes were sorted by their log-transformed BC and divided into equal- width bins, where higher bin indices correspond to higher average BC values. For each bin, we evaluated model performance (AUROC) separately for enzyme- and molecule- level betweenness for both the interaction and subclass tasks. This analysis enables us to assess whether predictive performance systematically varies with the structural importance of nodes within the bipartite network (Figure 5). 4.2 Implementation EMMA is implemented in Python 3.12.0 [54] using PyTorch 2.6.0 [55] as the primary deep learning framework. Data processing and numerical computations rely on Pan- das 2.2.3 [56] and NumPy 1.26.4 [57]. For domain-specific tasks, we employ BioPython 1.85 [58], RDKit 2025.03.2 [47], ChemDataExtractor 1.3.0 [59], ChEMBL Webresource Client 0.10.9 [60], and LibChEBIpy 1.0.10 [61]. Machine learning utilities are provided by Scikit-learn 1.6.1 [62] and SciPy 1.15.2 [63]. Visualization is performed with Mat- plotlib 3.10.1 [64] and Seaborn 0.13.2 [65]. Additional tools include TQDM 4.67.1 [66] for progress monitoring, DataSAIL 1.1.1 [23] for dataset splitting, Diamond 2.1.11 [67] for sequence alignment, Transformers 4.48.1 [68] for pretrained models, Fair-ESM 2.0.0 [69] for protein embeddings, and Py-cd-hit 1.1.4 [70] for sequence clustering. 20 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint Experiment tracking and logging are performed with Weights & Biases 0.19.8 [71], while Requests 2.32.3 [72] facilitates API queries, and Colorama 0.4.6 [73] enhances command-line output. 4.3 EMMA architecture We developed a transformer-based multi-task learning architecture that we called the Enzyme–small Molecule interaction Multi-head Attention (EMMA) model, to jointly predict the strength of enzyme–small molecule interactions and classify their interac- tion types (Figure 7). The LXMERT encoder inspired our encoder design [74], which employs dual-stream transformer mechanisms for vision and language modalities. Unlike LXMERT, EMMA omits intermediate feed-forward blocks inside the single- modality encoders and self-attention sub-layers inside the cross-modality encoder (see the original paper [74] for more information). We found that for the enzyme—small molecule interaction dataset, these additional layers introduced unnecessary complex- ity without improving performance. Instead, EMMA directly integrates the outputs of the self-modality encoders into the cross-modality encoder, simplifying the archi- tecture while maintaining its predictive power. Our model comprises three primary components: modality-specific encoders, dual-stream transformer mechanisms, and task-specific prediction heads. Each enzyme is represented using a single, mean-pooled, 640-dimensional embed- ding vector derived from the ESM2 t30 model [45], and each small molecule is encoded with a single, mean-pooled, 768-dimensional embedding vector from MolFormer XL [46] (Figure 7 a). To enable cross-attention between these two modalities, enzyme embeddings are projected into the same 768-dimensional latent space as the molecule embeddings using a learnable linear transformation, while small molecule embeddings are linearly transformed within the same dimensionality (Figure 7 b); both projections are followed by layer normalization. To model intra-modal dependencies, self-attention layers are applied independently to enzyme and small molecule embeddings. Each stream is processed by a multi- head attention mechanism with 32 heads, followed by residual connections and layer normalization. This is followed by cross-attention, where enzymes attend to the fea- tures of small molecules and vice versa, enabling the model to learn modality-specific interaction signals. Residual connections and layer normalization also follow each cross-attention. The resulting enzyme and small molecule representations are further processed through feed-forward neural networks with GELU activation and resid- ual connections and normalization (Figure 7 c). The refined representations are then concatenated to form a 1536-dimensional interaction feature vector (Figure 7 d). This shared representation is passed through two independent single-layer percep- trons: one for predicting binary interactions (interacting or non-interacting) and the other for subclassifying interacting pairs into enzyme-substrate and enzyme-inhibitor pairs. Each head consists of a linear layer (1024 units) followed by batch normalization, ReLU activation, dropout regularization (Figure 7 e). 21 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint + CrossSelf FF+ + + CrossSelf FF+ + ESM2_t30 MolFormer Interaction head Subclass head a b c d e Fig. 7 EMMA architecture. Self, Cross, and FF are abbreviations for self-attention, cross-attention, and feed-forward, respectively. Red arrows represent residual connection followed by normalization (summation sign). (a) Modality-specific encoders. ( b) Projection block. (c) Dual-stream transformer mechanisms. (d) Concatenation. ( e) Task-specific prediction heads. The model is trained using a composite loss function. A binary focal loss [75] is applied to the interaction head, while a second focal loss is used for subclass pre- dictions on interacting pairs only. Optimization is performed using AdamW[76] with learning rate scheduling via ReduceLROnPlateau, and gradient clipping is applied to ensure stability. To account for differences in data distribution and class balance across different dataset splits, we manually tuned key training and loss parameters for each split. Specifically, we adjusted the batch size, initial learning rate, weight decay, and the focal loss hyperparameters (gamma and alpha) for each split method. This allowed the model to better accommodate variations in sample size and label imbal- ance between splits, ensuring stable optimization and improved convergence for both tasks. To prevent overfitting and reduce unnecessary training time, we implemented an early stopping strategy based on the validation loss. A maximum of 100 epochs was set for all split methods. After each epoch, the model’s performance on the validation set is evaluated. If the validation loss does not improve for three consecutive epochs (patience = 3), training is terminated early. We also implemented a RF classifier as a baseline model. Two RF classifiers were trained: one for the interaction task and one for the subclass prediction task. The mod- els were trained with 100 trees, unlimited depth, the Gini criterion, and a minimum of two samples per split. Evaluation metrics (AUROC, accuracy, and F1-score) were computed separately for both tasks during training and validation. 5 Data Availability The train-test splits are publicly available at https://github.com/kalininalab/EMMA. EMMI dataset along with the raw and intermediate data generated during prepro- cessing are archived on Zenodo (https://doi.org/10.5281/zenodo.17280853). 22 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint 6 Code Availability All source code and implementation details are openly accessible at the GitHub repository https://github.com/kalininalab/EMMA. 7 Author Contribution V.A.E. conceived the study, developed the methodology, and implemented the pipeline. O.V.K. and R.J. contributed to study conception, validation, and jointly supervised the project. All authors wrote and edited the manuscript. 8 Acknowledgements R.J. was partly funded from the HelmhotzAI project XAI-Graph. O.V.K. acknowl- edges financial support from the Klaus Faber Foundation. 9 Competing Interests The authors declare no competing interests.

[1] Berg, J.M., Tymoczko, J.L., Stryer, L.: Biochemistry, pp. 524–560. W. H. Freeman, ??? (2012) [2] UniProt Consortium: Uniprot: the universal protein knowledgebase in 2021. Nucleic Acids Res. 49(D1), 480–489 (2021) https://doi.org/10.1093/nar/ gkaa1100 [3] Ramsundar, B., Eastman, P., Walters, P., Pande, V., Leswing, K., Wu, Z.: Deep Learning for the Life Sciences. O’Reilly Media, Sebastopol, CA (2019). https: //www.amazon.com/Deep-Learning-Life-Sciences-Microscopy/dp/1492039837 [4] Callender, R., Dyer, R.B.: The dynamical nature of enzymatic catalysis. Acc. Chem. Res. 48(2), 407–413 (2015) https://doi.org/10.1021/ar5002928 [5] Hammes, G.G., Benkovic, S.J., Hammes-Schiffer, S.: Flexibility, diversity, and cooperativity: pillars of enzyme catalysis. Biochemistry 50(48), 10422–10430 (2011) https://doi.org/10.1021/bi201486f [6] Riziotis, I.G., Ribeiro, A.J.M., Borkakoti, N., Thornton, J.M.: Conformational variation in enzyme catalysis: a structural study on catalytic residues. J. Mol. Biol. 434(7), 167517 (2022) https://doi.org/10.1016/j.jmb.2022.167517 [7] Tawfik, O.K., S, D.: Enzyme promiscuity: a mechanistic and evolutionary perspective. Annual review of biochemistry 79(1), 471–505 (2010) 23 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint [8] Kuo, Y.-M., Henry, R.A., Andrews, A.J.: Measuring specificity in multi- substrate/product systems as a tool to investigate selectivity in vivo. Biochim. Biophys. Acta Proteins Proteom.1864(1), 70–76 (2016) https://doi.org/10.1016/ j.bbapap.2015.08.011 [9] Hermann, J.C., Ghanem, E., Li, Y., Raushel, F.M., Irwin, J.J., Shoichet, B.K.: Predicting substrates by docking high-energy intermediates to enzyme structures. Journal of the American Chemical Society 128(49), 15882–15891 (2006) [10] Suh, D., Schwartz, R., Gupta, P.K., Zev, S., Major, D.T., Im, W.: Charmm-gui enzydocker for protein–ligand docking of multiple reactive states along a reaction coordinate in enzymes. Journal of Chemical Theory and Computation 21(4), 2118–2128 (2025) [11] Ram´ ırez-Palacios, C., Wijma, H.J., Thallmair, S., Marrink, S.J., Janssen, D.B.: Computational prediction of ω-transaminase specificity by a combination of dock- ing and molecular dynamics simulations. Journal of chemical information and modeling 61(11), 5569–5580 (2021) [12] Hao, J.-H., Zheng, D.-J., Ye, Y.-H., Yu, J.-T., Li, X.-Y., Xiong, M.-J., Jiang, W.-H., He, K.-P., Li, P.-Y., Lv, Y.-S.,et al.: Atomevo: a web server combining pro- tein modelling, docking, molecular dynamic simulation and mmpbsa analysis of candida antarctica lipase b (calb) fusion protein. Bioresources and Bioprocessing 9(1), 53 (2022) [13] Pertusi, D.A., Moura, M.E., Jeffryes, J.G., Prabhu, S., Biggs, B.W., Tyo, K.E.J.: Predicting novel substrates for enzymes with minimal experimental effort with active learning. Metab. Eng. 44, 171–181 (2017) https://doi.org/10.1016/j. ymben.2017.09.016 [14] Kroll, A., et al.: A general model to predict small molecule substrates of enzymes based on machine and deep learning. Nat. Commun. 14, 2787 (2023) https://doi. org/10.1038/s41467-023-38481-w [15] Du, Z., Fu, W., Guo, X., Caragea, D., Li, Y.: Fusionesp: Improved enzyme– substrate pair prediction by fusing protein and chemical knowledge. Journal of Chemical Information and Modeling 65(6), 2806–2817 (2025) [16] Kroll, A., Ranjan, S., Lercher, M.J.: A multimodal transformer network for protein-small molecule interactions enhances predictions of kinase inhibition and enzyme-substrate relationships. PLOS Computational Biology 20(5), 1012100 (2024) [17] Nie, Z., Zhang, H., Jiang, H., Liu, Y., Huang, X., Xu, F., Tian, Y., Chen, J., Zhang, W.-B.: Multi-purpose enzyme-substrate interaction prediction with progressive conditional deep learning (2024) 24 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint [18] Hassan, N.M., Alhossary, A.A., Mu, Y., Kwoh, C.-K.: Protein-ligand blind dock- ing using quickvina-w with inter-process spatio-temporal integration. Scientific reports 7(1), 15451 (2017) [19] Qian, W., Wang, X., Huang, Y., Kang, Y., Pan, P., Hsieh, C.-Y., Hou, T.: Deep learning-driven insights into enzyme–substrate interaction discovery. Journal of Chemical Information and Modeling 65(1), 187–200 (2024) [20] Bernett, J., Blumenthal, D.B., Grimm, D.G., et al. : Guiding questions to avoid data leakage in biological machine learning applications. Nat. Methods 21, 1444– 1453 (2024) https://doi.org/10.1038/s41592-024-02362-y [21] Steshin, S.: Lo-hi: practical ml drug discovery benchmark. Adv. Neural Inf. Process. Syst. 36, 64526–64554 (2023) [22] Teufel, F., G´ ıslason, M.H., Almagro Armenteros, J.J., Johansen, A.R., Winther, O., Nielsen, H.: Graphpart: homology partitioning for biological sequence analy- sis. NAR Genomics and Bioinformatics 5(4), 088 (2023) [23] Joeres, R., Blumenthal, D.B., Kalinina, O.V.: Data splitting to avoid information leakage with datasail. Nature Communications 16, 3337 (2025) https://doi.org/ 10.1038/s41467-025-58606-8 [24] Copeland, R.A.: Evaluation of Enzyme Inhibitors in Drug Discovery: A Guide for Medicinal Chemists and Pharmacologists. Wiley-Interscience, Hoboken, NJ (2005) [25] Silverman, R.B., Holladay, M.W.: The Organic Chemistry of Drug Design and Drug Action. Academic press, ??? (2014) [26] Hopkins, A.L.: Network pharmacology: the next paradigm in drug discovery. Nature chemical biology 4(11), 682–690 (2008) [27] Imming, P., Sinning, C., Meyer, A.: Drugs, their targets and the nature and number of drug targets. Nature reviews Drug discovery 5(10), 821–834 (2006) [28] Keser¨ u, G.M., Makara, G.M.: The influence of lead discovery strategies on the properties of drug candidates. nature reviews Drug Discovery 8(3), 203–212 (2009) [29] Gilson, M.K., Liu, T., Baitaluk, M., Nicola, G., Hwang, L., Chong, J.: Bindingdb in 2015: A public database for medicinal chemistry, computational chemistry and systems pharmacology. Nucleic Acids Research 44(D1), 1045–1053 (2016) https://doi.org/10.1093/nar/gkv1072 [30] Chang, A., Jeske, L., Ulbrich, S., Hofmann, J., Koblitz, J., Schomburg, I., Neumann-Schaal, M., Jahn, D., Schomburg, D.: BRENDA, the ELIXIR core data 25 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint resource in 2021: new developments and updates. Nucleic Acids Research 49, 498–508 (2021) https://doi.org/10.1093/nar/gkaa1025 [31] Zdrazil, B., Felix, E., Hunter, F., Manners, E.J., Blackshaw, J., Corbett, S., Veij, M., Ioannidis, H., Lopez, D.M., Mosquera, J.F., Magarinos, M.P., Bosc, N., Arcila, R., Kizil¨ oren, T., Gaulton, A., Bento, A.P., Adasme, M.F., Monecke, P., Landrum, G.A., Leach, A.R.: The chembl database in 2023: a drug discovery plat- form spanning multiple bioactivity data types and time periods. Nucleic Acids Research 52(D1), 1180–1192 (2024) https://doi.org/10.1093/nar/gkad1004 [32] Davies, M., Nowotka, M., Papadatos, G., Dedman, N., Gaulton, A., Atkinson, F., Bellis, L., Overington, J.P.: Chembl web services: streamlining access to drug discovery data and utilities. Nucleic Acids Research 43(W1), 612–620 (2015) https://doi.org/10.1093/nar/gkv352 [33] Ashburner, M., Ball, C.A., Blake, J.A., Botstein, D., Butler, H., Cherry, J., Davis, A., Dolinski, K., Dwight, S.S., Eppig, J., Harris, M., Hill, D., Issel-Tarver, L., Kasarskis, A., Lewis, S., Matese, J., Richardson, J., Ringwald, M., Saito, C., Schroeder, A., St. Pierre, S., Wu, G.: Gene ontology: tool for the unification of biology. Nature Genetics 25(1), 25–29 (2000) https://doi.org/10.1038/75556 [34] Consortium, T.G.O.: The gene ontology knowledgebase in 2023. Genetics 224(1), 031 (2023) https://doi.org/10.1093/genetics/iyad031 [35] Harding, S.D., Armstrong, J.F., Faccenda, E., Southan, C., Alexander, S.P.H., Davenport, A.P., Spedding, M., Davies, J.A.: The iuphar/bps guide to phar- macology in 2024. Nucleic Acids Research 52(D1), 1438–1449 (2024) https: //doi.org/10.1093/nar/gkad944 [36] Kim, S., Chen, J., Cheng, T., Gindulyte, A., He, J., He, S., Li, Q., Shoemaker, B.A., Thiessen, P.A., Yu, B., Zaslavsky, L., Zhang, J., Bolton, E.E.: Pubchem 2025 update. Nucleic Acids Research 53(D1), 1516–1525 (2025) https://doi.org/ 10.1093/nar/gkae1059 [37] Bansal, P., Morgat, A., Axelsen, K.B., Muthukrishnan, V., Coudert, E., Aimo, L., Hyka-Nouspikel, N., Gasteiger, E., Kerhornou, A., Neto, T.B., Pozzato, M., Blatter, M.-C., Ignatchenko, A., Redaschi, N., Bridge, A.: Rhea, the reaction knowledgebase in 2022. Nucleic Acids Research, 1016 (2021) https://doi.org/10. 1093/nar/gkab1016 [38] Wittig, U., Rey, M., Weidemann, A., Kania, R., M¨ uller, W.: Sabio-rk: an updated resource for manually curated biochemical reaction kinetics. Nucleic Acids Research 46(D1), 656–660 (2018) https://doi.org/10.1093/nar/gkx1065 [39] UniProt Consortium: Uniprot: the universal protein knowledgebase in 2023. Nucleic Acids Res. 51(D1), 523–531 (2023) https://doi.org/10.1093/nar/ gkac1137 26 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint [40] Kim, S., Thiessen, P.A., Cheng, T., Yu, B., Bolton, E.E.: An update on pug-rest: Restful interface for programmatic access to pubchem. Nucleic Acids Research 46(W1), 563–570 (2018) [41] Hastings, J., Owen, G., Dekker, A., Ennis, M., Kale, N., Muthukrishnan, V., Turner, S., Swainston, N., Mendes, P., Steinbeck, C.: Chebi in 2016: Improved ser- vices and an expanding collection of metabolites. Nucleic acids research 44(D1), 1214–1219 (2016) [42] Wishart, D.S., Guo, A., Oler, E., Wang, F., Anjum, A., Peters, H., Dizon, R., Sayeeda, Z., Tian, S., Lee, B.L.,et al.: Hmdb 5.0: the human metabolome database for 2022. Nucleic acids research 50(D1), 622–631 (2022) [43] Frolkis, A., Knox, C., Lim, E., Jewison, T., Law, V., Hau, D.D., Liu, P., Gautam, B., Ly, S., Guo, A.C.,et al.: Smpdb: the small molecule pathway database. Nucleic acids research 38(suppl 1), 480–487 (2010) [44] Kanehisa, M., Goto, S.: Kegg: kyoto encyclopedia of genes and genomes. Nucleic acids research 28(1), 27–30 (2000) [45] Lin, Z., Akin, H., Rao, R., Hie, B., Zhu, Z., Lu, W., Smetanin, N., Santos Costa, A., Fazel-Zarandi, M., Sercu, T., Candido, S., et al.: Language models of protein sequences at the scale of evolution enable accurate structure prediction. bioRxiv (2022) [46] Ross, J., Belgodere, B., Chenthamarakshan, V., Padhi, I., Mroueh, Y., Das, P.: Large-scale chemical language representations capture molecular structure and properties. Nature Machine Intelligence 4(12), 1256–1264 (2022) [47] Landrum, G.: RDKit: Open-source cheminformatics. https://www.rdkit.org. Accessed: 2023-01-01 (2016) [48] Kalliokoski, T., Kramer, C., Vulpetti, A., Gedeck, P.: Comparability of mixed ic50 data – a statistical analysis. PLOS ONE 8(4), 61007 (2013) https://doi.org/ 10.1371/journal.pone.0061007 [49] Yung-Chi, C., Prusoff, W.H.: Relationship between the inhibition constant (ki) and the concentration of inhibitor which causes 50 per cent inhibition (i50) of an enzymatic reaction. Biochemical pharmacology 22(23), 3099–3108 (1973) [50] PubChem BioAssay: BACE-1 Ki Assay (BACE-1 HTRF FRET Assay) from US Patent US10071998: ”Iminothiadiazine dioxides bearing an amine-linked sub- stituent as BACE inhibitors, compositions, and their use”.https://pubchem.ncbi. nlm.nih.gov/bioassay/1795532#section=Description. Accessed: 17 August 2025 (2025) [51] BioNumbers Database: KM (Michaelis constant) of some enzymes. https:// 27 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint bionumbers.hms.harvard.edu/bionumber.aspx?s=n&v=3&id=105085. Accessed: 17 August 2025 (2025) [52] BioNumbers Database: Median KM for the entire data set of enzymes and nat- ural substrates. https://bionumbers.hms.harvard.edu/bionumber.aspx?s=n&v= 1&id=111413. Accessed: 17 August 2025 (2025) [53] Brandes, U.: A faster algorithm for betweenness centrality. Journal of mathemat- ical sociology 25(2), 163–177 (2001) [54] VanRossum, G., Drake, F.L.: The Python Language Reference vol. 561. Python Software Foundation Amsterdam, The Netherlands, ??? (2010) [55] Ansel, J., Yang, E., He, H., Gimelshein, N., Jain, A., Voznesensky, M., Bao, B., Bell, P., Berard, D., Burovski, E., Chauhan, G., Chourdia, A., Constable, W., Desmaison, A., DeVito, Z., Ellison, E., Feng, W., Gong, J., Gschwind, M., Hirsh, B., Huang, S., Kalambarkar, K., Kirsch, L., Lazos, M., Lezcano, M., Liang, Y., Liang, J., Lu, Y., Luk, C., Maher, B., Pan, Y., Puhrsch, C., Reso, M., Saroufim, M., Siraichi, M.Y., Suk, H., Suo, M., Tillet, P., Wang, E., Wang, X., Wen, W., Zhang, S., Zhao, X., Zhou, K., Zou, R., Mathews, A., Chanan, G., Wu, P., Chin- tala, S.: PyTorch 2: Faster Machine Learning Through Dynamic Python Bytecode Transformation and Graph Compilation. In: 29th ACM International Conference on Architectural Support for Programming Languages and Operating Systems, Volume 2 (ASPLOS ’24). ACM, ??? (2024). https://doi.org/10.1145/3620665. 3640366 . https://docs.pytorch.org/assets/pytorch2-2.pdf [56] McKinney, W., et al. : Data structures for statistical computing in python. scipy 445(1), 51–56 (2010) [57] Harris, C.R., Millman, K.J., Van Der Walt, S.J., Gommers, R., Virtanen, P., Cournapeau, D., Wieser, E., Taylor, J., Berg, S., Smith, N.J., et al. : Array programming with numpy. nature 585(7825), 357–362 (2020) [58] Cock, P.J., Antao, T., Chang, J.T., Chapman, B.A., Cox, C.J., Dalke, A., Fried- berg, I., Hamelryck, T., Kauff, F., Wilczynski, B., et al. : Biopython: freely available python tools for computational molecular biology and bioinformatics. Bioinformatics 25(11), 1422 (2009) [59] Swain, M.C., Cole, J.M.: Chemdataextractor: a toolkit for automated extrac- tion of chemical information from the scientific literature. Journal of chemical information and modeling 56(10), 1894–1904 (2016) [60] Davies, M., Nowotka, M., Papadatos, G., Dedman, N., Gaulton, A., Atkinson, F., Bellis, L., Overington, J.P.: Chembl web services: streamlining access to drug discovery data and utilities. Nucleic acids research 43(W1), 612–620 (2015) [61] Swainston, N., Hastings, J., Dekker, A., Muthukrishnan, V., May, J., Steinbeck, 28 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint C., Mendes, P.: libchebi: an api for accessing the chebi database. Journal of Cheminformatics 8(1), 11 (2016) [62] Pedregosa, F., Varoquaux, G., Gramfort, A., Michel, V., Thirion, B., Grisel, O., Blondel, M., Prettenhofer, P., Weiss, R., Dubourg, V.,et al.: Scikit-learn: Machine learning in python. the Journal of machine Learning research 12, 2825–2830 (2011) [63] Virtanen, P., Gommers, R., Oliphant, T.E., Haberland, M., Reddy, T., Courna- peau, D., Burovski, E., Peterson, P., Weckesser, W., Bright, J., van der Walt, S.J., Brett, M., Wilson, J., Millman, K.J., Mayorov, N., Nelson, A.R.J., Jones, E., Kern, R., Larson, E., Carey, C.J., Polat, ˙I., Feng, Y., Moore, E.W., Vander- Plas, J., Laxalde, D., Perktold, J., Cimrman, R., Henriksen, I., Quintero, E.A., Harris, C.R., Archibald, A.M., Ribeiro, A.H., Pedregosa, F., van Mulbregt, P., SciPy 1.0 Contributors: SciPy 1.0: Fundamental Algorithms for Scientific Com- puting in Python. Nature Methods 17, 261–272 (2020) https://doi.org/10.1038/ s41592-019-0686-2 [64] Hunter, J.D.: Matplotlib: A 2d graphics environment. Computing in science & engineering 9(03), 90–95 (2007) [65] Waskom, M., Botvinnik, O., O’Kane, D., Hobson, P., Lukauskas, S., Gemperline, D.C., Augspurger, T., Halchenko, Y., Cole, J.B., Warmenhoven, J., Ruiter, J., Pye, C., Hoyer, S., Vanderplas, J., Villalba, S., Kunter, G., Quintero, E., Bachant, P., Martin, M., Meyer, K., Miles, A., Ram, Y., Yarkoni, T., Williams, M.L., Evans, C., Fitzgerald, C., Brian, Fonnesbeck, C., Lee, A., Qalieh, A.: Mwaskom/seaborn: V0.8.1 (September 2017). https://doi.org/10.5281/zenodo.883859 . https://doi. org/10.5281/zenodo.883859 [66] Costa-Luis, C.O.: tqdm: A fast, extensible progress meter for python and cli. Journal of Open Source Software 4(37), 1277 (2019) [67] Buchfink, B., Xie, C., Huson, D.H.: Fast and sensitive protein alignment using diamond. Nature methods 12(1), 59–60 (2015) [68] Wolf, T., Debut, L., Sanh, V., Chaumond, J., Delangue, C., Moi, A., Cistac, P., Rault, T., Louf, R., Funtowicz, M., et al. : Transformers: State-of-the-art natu- ral language processing. In: Proceedings of the 2020 Conference on Empirical

in Natural Language Processing: System Demonstrations, pp. 38–45 (2020) [69] Rives, A., Meier, J., Sercu, T., Goyal, S., Lin, Z., Liu, J., Guo, D., Ott, M., Zitnick, C.L., Ma, J., et al.: Biological structure and function emerge from scal- ing unsupervised learning to 250 million protein sequences. Proceedings of the National Academy of Sciences 118(15), 2016239118 (2021) [70] Li, W., Godzik, A.: Cd-hit: a fast program for clustering and comparing large 29 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint sets of protein or nucleotide sequences. Bioinformatics 22(13), 1658–1659 (2006) [71] Biewald, L.: Experiment Tracking with Weights and Biases. Software available from wandb.com (2020). https://www.wandb.com/ [72] Reitz, K.: Requests documentation. EEUU: Requests Documentation (2016) [73] colorama. https://pypi.org/project/colorama/. Accessed: 2024-10-17 [74] Tan, H., Bansal, M.: Lxmert: Learning cross-modality encoder representations from transformers. arXiv preprint arXiv:1908.07490 (2019) [75] Lin, T.-Y., Goyal, P., Girshick, R., He, K., Doll´ ar, P.: Focal loss for dense object detection. In: Proceedings of the IEEE International Conference on Computer Vision, pp. 2980–2988 (2017) [76] Loshchilov, I., Hutter, F.: Decoupled weight decay regularization. arXiv preprint arXiv:1711.05101 (2017) 30 .CC-BY-NC 4.0 International licenseperpetuity. It is made available under a preprint (which was not certified by peer review) is the author/funder, who has granted bioRxiv a license to display the preprint in The copyright holder for thisthis version posted October 11, 2025. ; https://doi.org/10.1101/2025.10.09.681419doi: bioRxiv preprint