We used the protein with Protein Data Bank (PDB) ID 4HG7 (chain A) (Berman et al., 2000) as the receptor. This crystal structure has a 1.6 Å resolution and contains a co-crystallized scaffold, a small molecule non-covalently bound to the MDM2 binding site of p53 (Zhang et al., 2022). We used Chimera (Pettersen et al., 2004) and its Dock Prep tool to prepare the receptor structure and the co-crystallized ligand. Solvent molecules and the co-crystallized small molecule sulfate ion were removed, missing residues were completed, and partial charges for histidines in the receptor were assigned using the AM1-BCC method (Jakalian et al., 2002).

For the molecules in the dataset, we used Open Babel (O’Boyle et al., 2011) to generate 3D conformations using the MMFF94 force field. Once the structures were prepared, we performed molecular docking using Smina (Morris et al., 2009; Oyedele et al., 2022). The search grid was defined using the co-crystallized ligand as a reference. The grid size was set to 30 Å in each axis, providing enough space for all ligands to rotate freely. We selected one pose for each molecule, meaning the one with the lowest (most negative) Smina docking score. Smina uses a stochastic sampling method to create different poses, but each docking run repeats the optimization eight times by default. This means that the best pose for each molecule is usually stable, even if the docking is repeated (Torres et al., 2019; Masters et al., 2020). Therefore, we refer to this pose as the most representative one for each molecule in the dataset.

To confirm that the docking poses were chemically and physically consistent, and not influenced by artifacts from Smina, we applied orthogonal validation using a reliable control set, following the method of Bender et al (Nat. Protoc. 2021) (Bender et al., 2021). First, we grouped the actives from our dataset into clusters by 2D similarity, which was calculated using ECFP4 fingerprints with a Tanimoto cutoff of 0.35. From each cluster, we kept only the most potent compound. Finally, this process gave us a control set of 30 active ligands that were structurally diverse and showed the highest experimental potencies.

These ligands were analyzed with Protein–Ligand Interaction Profiler (PLIP) to detect hydrophobic and aromatic contacts (Schake et al., 2025). First, we grouped the actives from our dataset into clusters by 2D similarity, which was calculated using ECFP4 fingerprints with a Tanimoto cutoff of 0.35. From each cluster, we kept only the most potent compound. This yielded a structurally diverse control set of 30 active ligands with the highest experimental potencies. Each ligand–MDM2 complex was then analyzed with the Protein–Ligand Interaction Profiler (PLIP) to identify hydrophobic and aromatic contacts. PLIP reported the interacting atoms, residues, and distances, confirming the expected key contacts across representative ligands. Furthermore, to evaluate whether Smina produced physically realistic docking poses, we used PoseBusters (Buttenschoen et al., 2024) which applies 18 criteria to assess the physical plausibility of docking poses, including intra- and intermolecular consistency. PoseBusters outputs a boolean descriptor (PB_valid) indicating the overall pose validity.

Protein–ligand extended connectivity (PLEC) fingerprints were generated using the PLEC function from the ODDT Python package (Wójcikowski et al., 2015), following the method described in previous studies (Tran-Nguyen et al., 2023). The fingerprint size was set to 4092, with ligand and protein fragment depths of 1 and 5, respectively. Additionally, 3D Grid-based (Grid) features were computed using the RDKitGridFeaturizer function from the DeepChem Python package (Feinberg et al., 2018), based on the approach described in previous studies (Gómez-Sacristán et al., 2025). The extracted features included extended connectivity fingerprints (ECFP) (Rogers and Hahn, 2010), spatial ligand-based interaction fingerprints (SPLIF) (Da and Kireev, 2014), hydrogen bonds, and salt bridges. The selected parameters were ecfp_power = 9, splif_power = 9, and voxel_width = 16.0, resulting in a total of 2052 features. To improve efficiency, feature extraction was performed in parallel using 20 cores, and the resulting features were saved as.npy files for further analysis.

The use of ML scoring functions has significantly improved the performance of SBVS for specific targets (Scantlebury et al., 2023; Zhang et al., 2023; Junaid et al., 2024; Gómez-Sacristán et al., 2025). This is because ML models can extract and learn target-specific 3D molecular features and correlate them to the experimental activity of these molecules (Gómez-Sacristán et al., 2025; Tran-Nguyen et al., 2023). As a result, ML-based scoring functions for virtual screening can achieve higher accuracy compared to general scoring functions, which rely on recognizing generic target features and applying them to evaluate molecular interactions. Furthermore, previous studies have demonstrated that developing ML models using a dataset enriched with inactive compounds enhances performance (Gómez-Sacristán et al., 2025; Kurczab et al., 2014; Réau et al., 2018). Therefore, to develop an MDM2-specific ML model, and considering the aforementioned insights, we constructed and evaluated 5 ML algorithms using a dataset enriched with inactive compounds and 3D interaction features. The results are summarized in Supplementary Material: Supplementary Table S4 and Figure 4.

The performance of ML models was evaluated using the PRC and AUC (Keilwagen et al., 2014). These metrics help measure how well the model can separate active and inactive compounds, which is important for handling imbalanced datasets (Tran-Nguyen et al., 2023; Richardson et al., 2024). A good model has an AUC-PRC close to 1, meaning it can distinguish well between the two classes. In drug discovery, AUC-PRC values between 0.6 and 0.8 are considered acceptable (Gómez-Sacristán et al., 2025; Caba et al., 2024; Zhu et al., 2025; Tahir ul Qamar et al., 2022). For example, Zhou et al., in 2024 showed that ML models for virtual screening achieved a similar level of performance (Zhou et al., 2024). In summary, we evaluated the performance of ML models for discovering MDM2 inhibitors targeting the p53 binding site and compared their results to those of previous virtual screening studies to examine whether they achieved similar or better results. A more detailed comparison with traditional methods is presented in the section of Comparison with Generic Virtual Screening Models.

To achieve this, we tested five supervised learning algorithms using a classification approach to explore their predictive capabilities. The models were trained and tested to classify each compound as either active or inactive. Although the selected algorithms such as SVM, RF, XGB, ANN, and DNN can also be used for regression tasks, in this study we only applied them in classification mode. These models were chosen because they have been successfully applied in various drug discovery studies (Patel et al., 2020; Zhou et al., 2024; Caba et al., 2024; Ballester and Mitchell, 2010) and are also described in the protocol we followed, which is based on the approach proposed by Tran-Nguyen et al (2023) in their guide on ML-based scoring functions for virtual screening (Tran-Nguyen et al., 2023).

To train these ML models with representative data, a numerical representation of molecular structures is required. Therefore, we encoded the 3D interaction information of the MDM2-L protein-ligand complex at the p53 binding site using the chemical descriptors Grid and PLEC. These descriptors were applied to both active and inactive compounds from the previously constructed database. Grid and PLEC numerically encode the 3D interaction information between the protein and ligands derived from docking poses, and have been successfully used in prior virtual screening studies (Caba et al., 2024; Gómez-Sacristán et al., 2025). Notably, previous studies have reported that ML model performance improves when descriptors incorporate 3D protein-ligand interaction information.

To evaluate the effectiveness of these descriptors in our models, we examined the performance of ML algorithms using Grid-ML and PLEC-ML representations. Figure 4 shows the combined prediction results of the Grid-ML and PLEC-ML models across 10 test runs (Keilwagen et al., 2014). For the Grid dataset (Figure 4A), DNN achieved the highest AUC-PRC (0.353), followed by SVM (0.309). However, the results show that all AUC-PRC values were relatively low, indicating that Grid descriptors do not provide good performance for machine learning models.

In contrast, the performance improved significantly when using the PLEC dataset, as shown in Figure 4B. The SVM achieved the highest AUC-PRC (0.767), followed closely by RF with 0.766 and XGB with 0.756. These results suggest that PLEC descriptors provide more useful features for classification. However, in this case, the DNN had the lowest performance, with an AUC-PRC of 0.525. This indicates that more complex models are not always the best choice, as their effectiveness depends on the type of descriptor and specific target. Therefore, according to these results, simpler or more traditional models such as PLEC-SVM, PLEC-RF, or PLEC-XGB may achieve higher performance in correctly classifying true actives (inhibitors), demonstrating greater discriminatory power between the active and inactive compounds for MDM2, as evidenced by their superior performance in the PRC, which is a well-established metric for evaluating classifier performance in imbalanced datasets (Davis and Goadrich, 2006; Saito and Rehmsmeier, 2015).

Within our Grid-ML and PLEC-ML frameworks, specific ML models such as PLEC-RF, and PLEC-SVM showed a remarkable ability to correctly classify active compounds (inhibitors), demonstrating a strong discriminatory power between active and inactive compounds for MDM2. To evaluate the predictive ability of our models in selecting the most active compounds, we used the EF1% metric calculated on the test set, which measures a model’s ability to retrieve active compounds within the top 1% of its predictions (Truchon and Bayly, 2007). Specifically, in SBVS, a high EF1% indicates that the model efficiently identifies the most promising compounds and rank them at the top of the list based on affinity or activity probability (Zhou et al., 2024; Li H. et al., 2021; Mahdizadeh and Eriksson, 2025). This is crucial for prioritizing candidates for experimental testing.

The results of our models are shown in Figure 5. For the PLEC-ML models, the mean EF1% values ranged between 2.742 and 3.450, respectively, for PLEC-RF (3.450 ± 0.00), PLEC-XGB (3.290 ± 0.03), PLEC-SVM (3.450 ± 0.00), PLEC-ANN (3.290 ± 0.05), and PLEC-DNN (2.742 ± 0.18). These results indicate that the PLEC-ML models identified between 2.742 and 3.450 times more active compounds than random selection. On the other hand, the mean EF1% values obtained for the Grid-ML models ranged from 0.471 to 2.030, respectively, for Grid-RF (1.420 ± 0.00), Grid-XGB (0.810 ± 0.03), Grid-SVM (2.030 ± 1.40), Grid-ANN (0.471 ± 0.06), and Grid-DNN (0.772 ± 0.16). Thus, the Grid-ML models only identified between 0.772 and 2.030 times more active compounds than random selection. Therefore, it is remarkable that all PLEC-ML models outperformed all Grid-ML models. This is consistent with the previously presented results.

Specifically, for PLEC, the PLEC-RF and PLEC-SVM models identified 3.450 times more active compounds, respectively, compared to random selection. This confirms that PLEC descriptors capture relevant structural information for classifying MDM2 inhibitors. Overall, these results support the use of classical ML approaches, such as RF and SVM, combined with PLEC descriptors to improve the selection of MDM2 inhibitors in SBVS processes.

Although EF1% indicates the ability to retrieve compounds in the top 1% of VS, this metric is not normalized. This means it does not confirm if a model is 100% efficient in classifying only actives in the top 1%. To address this, we measured the Enrichment Factor Normalized to 1% (NEF1%), which is calculated as the observed EF1% divided by the maximum possible EF1% for a given test set (Liu et al., 2018). NEF1% indicates how well a model classifies actives in the top 1%. A value of 1 or close to it means the model selects only actives in the top 1%. The results of this analysis are shown in Figure 6, with Supplementary Material: Supplementary Table S5.

In Figure 6, we can observe the normalized enrichment factor at 1% mean (NEF1%) for different models using Grid and PLEC descriptors. The results indicate that PLEC-RF (1.000 ± 0.00), PLEC-XGB (0.952 ± 0.01), PLEC-SVM (1.000 ± 0.00), PLEC-ANN (0.952 ± 0.02), and PLEC-DNN (0.792 ± 0.05) achieved the highest NEF1% values. This means that they efficiently classified actives within the top 1% of the ranked compounds. These models demonstrated strong predictive capabilities in prioritizing MDM2 inhibitors.

However, Grid-RF (0.410 ± 0.00) and Grid-SVM (0.590 ± 0.00) obtained moderate NEF1% values, while Grid-XGB (0.24 ± 0.02), Grid-ANN (0.136 ± 0.02), and Grid-DNN (0.227 ± 0.05) showed a low NEF1% value. This suggests that, although they retrieved active compounds, they were not effective in ranking them exclusively at the top, highlighting their inability to effectively classify actives at the top 1% of the predictions.

These results reinforce the effectiveness of classical ML models, such as RF and SVM, particularly when combined with PLEC descriptors. This combination proves highly suitable for identifying MDM2 inhibitors in SBVS. The higher NEF1% values for PLEC-based models confirm their effectiveness. These descriptors provide meaningful structural information, which helps prioritize active compounds better.

To assess the models’ consistent performance without overfitting and their generalization capability, we performed cross-validation on the training set and evaluated ROC-AUC scores on the independent test set (see Methods section). ROC-AUC shows the balance between sensitivity and specificity. In other words, it measures model’s ability to identify active compounds and distinguish them from inactive ones (Robinson et al., 2020; LeDell et al., 2015). Besides, it helps to estimate the performance across different datasets. Ideal ROC-AUC values are near to one because this indicates that can perfectly distinguish between active and inactive classes (LeDell et al., 2015; Robinson et al., 2020). Nonetheless, it is an ideal model. Some studies on drug discovery using ML models have showed that ROC-AUC >0.8 yields good results for predicting new or repurposed drugs (Fan et al., 2019; Monsia and Bhattacharyya, 2025; Zhou et al., 2024; Budennyy et al., 2023). In addition, this metric allows us to determine whether a model has overfitting. In brief, good generalization capability is demonstrated with ROC-AUC> 0.8 in drug discovery.

As a result, we calculated and compared the ROC-AUC on the training set using five-fold cross-validation, reporting the mean and standard deviation (SD) to assess consistent performance without overfitting. Additionally, ROC-AUC was calculated on the independent test set to evaluate the generalization capability and the quality of predictions for different MDM2 inhibitors. ROC-AUC mean and SD scores for the training set models are summarized in Supplementary Figure S2, while ROC-AUC mean and SEM scores for the test set models are reported in Figure 7.

The results in Supplementary Figure S2 show the ROC-AUC from cross-validation on the training set. Values remained high across folds for PLEC-RF (0.810 ± 0.01), PLEC-XGB (0.812 ± 0.02), PLEC-SVM (0.820 ± 0.006), and PLEC-ANN (0.798 ± 0.01) indicating low variability and minimal overfitting. In contrast, PLEC-DNN achieved a ROC-AUC of 0.714 ± 0.07, which was lower than the other models and exhibited higher variability. Figure 7 presents the mean ROC-AUC scores on the independent test set illustrating each model’s generalization ability to unseen MDM2 inhibitors. The models achieved the following mean ROC-AUC scores: PLEC-RF (0.866 ± 0.0007), PLEC-XGB (0.862 ± 0.0001), PLEC-SVM (0.871 ± 3.8 × 10⁻⁶), and PLEC-ANN (0.747 ± 0.061). These results are consistent with the training set performance, confirming that the models generalize well to unseen data. This indicates that the best-performing PLEC models not only achieve high classification ability but also maintain robustness. In contrast, PLEC-DNN remained at chance level (0.500 ± 0.0001), indicating poor generalization in this dataset.

Compared to the previous models, Grid-based models showed lower ROC-AUC values than PLEC-based models on both the test and training sets. The results in the test set (Figure 7) were: Grid-RF (0.478 ± 1.00 × 10⁻⁵), Grid-XGB (0.477 ± 2.00 × 10⁻⁵), Grid-SVM (0.513 ± 2.00 × 10⁻⁶), Grid-ANN (0.505 ± 1.00 × 10⁻²), and Grid-DNN (0.484 ± 5.00 × 10⁻⁴). Furthermore, in Supplementary Figure S2, ROC-AUC values for the training set were near 0.5. This indicates that Grid-based models have poor performance in both the training and test sets and do not have a predictive advantage over random classification because they have low generalization capacity.

In summary, these results highlight the importance of selecting an appropriate feature representation for the specific MDM2-inhibitor system. In this case, the PLEC strategy effectively captures structural and interaction information relevant for distinguishing active and inactive MDM2 inhibitors, whereas the Grid strategy does not. Although previous SBVS performance suggested that both PLEC-SVM and PLEC-RF were strong contenders, PLEC-SVM achieved a slightly higher mean ROC-AUC during cross-validation (0.820 ± 0.01) and lower standard deviation compared to PLEC-RF (0.800 ± 0.01). However, a paired t-test (p = 0.075) showed that this difference was not statistically significant. Thus, both models exhibit strong predictive performance, with PLEC-SVM showing a slight, though not statistically confirmed, advantage. In conclusion, these results confirm that PLEC-SVM and PLEC-RF are the most reliable models, combining strong training performance with good generalization to independent test data.

Some studies have shown that reclassifying docking poses with ML improves the identification of active molecules compared to generic docking methods for specific targets (Gómez-Sacristán et al., 2025; Caba et al., 2024; Li et al., 2014; Hany et al., 2025). To evaluate this in our case, we focused on MDM2 inhibitors at the p53 binding site and measured the mean PRC and AUC scores, EF1% and NEF1% between our models and classical generic methods like Smina (Morris et al., 2009), SCORCH-generic (McGibbon et al., 2023), and CNN-generic (Ragoza et al., 2017). These methods evaluate docking poses by generating scores that approximate or correlate with the binding free energy (ΔG _binding) between the protein and ligand, rather than computing it explicitly. However, each method employs different scoring functions and algorithms to estimate binding affinity.

Smina, a derivative of AutoDock Vina 1.1.2, calculates ΔG_binding using a linear regression model (Morris et al., 2009). This approach assumes that weak interactions between proteins and ligands can be modeled as linear relationships. Smina was trained with the Community Structure-Activity Resource (CSAR)-NRC HiQ 2010 dataset, which is known for its high-quality protein-ligand complex data (Koes et al., 2013). As a result, several studies have shown that it improves the performance of scoring docking poses of active compounds against decoys across multiple biological targets. Its performance was benchmarked against ten other docking programs (MVP, DockIt, FlexX, Flo+, Fred, Dock4, Glide, Gold, MOEDock, and LigandFit) (Warren et al., 2006; Chaput et al., 2016; Cross et al., 2009) using the EF1% metric. In this evaluation, Smina achieved an EF1% of 11.6, significantly outperforming the best competing programs, which had EF1% values between 5 and 7. These results highlight Smina’s superior ability to prioritize active compounds in virtual screening, on average across a range of evaluated biological targets.

SCORCH is a ML model that combines gradient-boosted decision trees (XGB), feedforward neural networks (FNN), and wide and deep neural networks (DNN). It was trained on protein-ligand datasets like PDBbind, Binding MOAD, and Iridium, along with decoys created using DeepCoy (McGibbon et al., 2023). Then, it was tested on different targets and is considered highly effective because it outperforms traditional scoring functions in virtual screening (VS) and pose prediction.

Finally, CNN-Score is a deep learning model consisting of five convolutional neural networks (CNN) with up to 20 hidden layers. It was trained on protein-ligand interaction data, including decoys from the DUD-E database, which contains 22,645 positive instances and 1,407,145 negative ones (Mysinger et al., 2012). Then, it was tested on the LIT-PCBA dataset (Tran-Nguyen et al., 2020), designed for virtual screening and ML with 149 dose-response PubChem bioassays. CNN-Score outperformed traditional scoring functions like Smina and improved pose re-scoring without overestimating detection power. We employed these high-performing models for a comparative analysis on our specific target (MDM2 at the p53 binding site), and the results are presented in Figure 8 and Supplementary Table S6.

In Figure 8, the mean PRC is shown for the best ML models that were previously evaluated. The figure includes the models based on PLEC (SVM, RF, XGB, ANN, and DNN) and the generic scoring functions SCORCH, CNN-Score, and Smina. The PLEC-ML models show a decreasing order of mean scores, starting with SVM (0.767), RF (0.766), XGB (0.756), ANN (0.665), and DNN (0.525) while generic models show lower values, such as CNN (0.131), SCORCH (0.023), and Smina (0.022).

With regard to the foregoing, the PLEC-ML models outperform the generic scoring functions in terms of AUC-PRC. Then, we can see that, although CNN has shown good results in other studies, (Ragoza et al., 2017), in the specific case of identifying inhibitors targeting the p53 binding site in MDM2, its AUC-PRC is lower than that of the generic models. This indicates that its ability to re-rank docked poses in this dataset is limited. Besides, SCORCH, which is based on decision trees enhanced with neural networks, also exhibited very low performance in our specific case, suggesting that its re-scoring capability is not optimal in this scenario. Finally, Smina, a scoring function based on linear regression, showed the worst performance in this specific case. This confirms that classical docking approaches may be less effective in distinguishing active from inactive docked poses in specific targets.

In addition, Supplementary Table S6 illustrates the model’s ability to retrieve active compounds within the top 1% of the ranking generated using the NEF1% metric. The results indicate that the PLEC-ML models have scores closest to 1 whereas the generic models have scores closest to 0. Firstly, the SCORCH generic model showed an NEF1% of 0.30, meaning that it can only prioritize the classification of actives by 30% in the top 1% ranking. Therefore, this indicates that this model has a lower enrichment capability in retrieving active compounds at the top positions compared to PLEC-ML models. Secondly, CNN showed an NEF1% of 0.30, meaning that it can only prioritize the classification of active compounds by 30% in the top 1% ranking. Therefore, this suggests that this model has a lower enrichment capability in retrieving active compounds at the top positions compared to PLEC-ML models. Third, Smina exhibited an NEF1% of 0 implying that this model cannot prioritize molecules in the top-ranking search for actives, specifically inhibitors of the p53 binding site in MDM2.

In summary, the results show that ML models trained with PLEC descriptors outperform traditional and generic scoring functions in identifying inhibitors of the p53 binding site in MDM2. In particular, the PLEC-SVM and PLEC-RF models have a high capability to prioritize/classify active molecules in the top 1%. This suggests that re-scoring with PLEC-SVM and PLEC-RF significantly improves the prediction of actives in the top-ranked classification performed by SBVS.

A ML model designed to predict the re-scoring of docking poses and potential inhibitors must be capable of identifying diverse molecular structures. This ensures that the model does not favor specific compounds while overlooking others with inhibitory potential. This aspect is essential to improving its applicability in drug discovery (Nogueira and Koch, 2019; Tran-Nguyen et al., 2023).

Thus, we evaluated the chemical diversity of PLEC-SVM and PLEC-RF, both of which demonstrated a high generalization capacity and strong performance in classifying molecules with inhibitory potential for p53 binding site in MDM2. This analysis ensures that the model recognizes not only compounds similar to those in the training set but also molecular structures that differ from it.

We evaluated chemical diversity by calculating pairwise Tanimoto similarity using Morgan fingerprints (2048 bits, radius 2), applying a similarity threshold of 0.7. Compounds with a Tanimoto similarity equal to or greater than this threshold were grouped together, following the approach described by Tran-Nguyen et al. (2023). The classified molecules correspond to the true active compounds identified by PLEC-SVM and PLEC-RF in the top 1%. The results are shown in Figure 9.

The results show that for PLEC-SVM there are 11 different groups of molecular structures, which were clustered based on Tanimoto similarity coefficients higher than 0.7. This indicates that the PLEC-SVM model can identify distinct molecular structures as true actives, specifically inhibitors of the p53 binding site in MDM2. In contrast, PLEC-RF identified 7 different groups of molecular structures, also clustered by Tanimoto coefficients greater than 0.7. Therefore, among the groups identified by PLEC-SVM and PLEC-RF, it is evident that both models can find molecules with structural chemical diversity. This difference is expected in ML models (Caba et al., 2024), as each can capture different patterns or chemical features, which can be beneficial for exploring chemical diversity in the search for inhibitors of the p53 binding site in MDM2. Therefore, the application of virtual screening by structure-based methods (SBVS) is highly promising in this context.

Moreover, it was observed that both PLEC-SVM and PLEC-RF models coincide on a subset of molecules with high similarity (7 common: CHEMBL1834286, CHEMBL1688304, CHEMBL3927813, CHEMBL3800573, CHEMBL3260848, CHEMBL1688271, CHEMBL4451740), since both models identified these molecules as representative or relevant according to their respective criteria (SVM using kernels and RF using decision trees).

We think that the PLEC-SVM and PLEC-RF models find different chemical patterns. The SVM model looks for more diverse molecules, while the RF model focuses on molecules with stronger evidence. Because of this, using both models together may help find more different and important molecules, which can improve virtual screening for MDM2 inhibitors at the p53 binding site.