3D-QSAR Studies, Molecular Docking, Molecular Dynamic Simulation, and ADMET Proprieties of Novel Pteridinone Derivatives as PLK1 Inhibitors for the Treatment of Prostate Cancer

Overexpression of polo-like kinase 1 (PLK1) has been found in many different types of cancers. With its essential role in cell proliferation, PLK1 has been determined to be a broad-spectrum anti-cancer target. In this study, 3D-QSAR, molecular docking, and molecular dynamics (MD) simulations were applied on a series of novel pteridinone derivatives as PLK1 inhibitors to discover anti-cancer drug candidates. In this work, three models—CoMFA (Q² = 0.67, R² = 0.992), CoMSIA/SHE (Q² = 0.69, R² = 0.974), and CoMSIA/SEAH (Q² = 0.66, R² = 0.975)—of pteridinone derivatives were established. The three models that were established gave Rpred2 = 0.683, Rpred 2= 0.758, and Rpred 2= 0.767, respectively. Thus, the predictive abilities of the three proposed models were successfully evaluated. The relations between the different champs and activities were well-demonstrated by the contour chart of the CoMFA and CoMSIA/SEAH models. The results of molecular docking indicated that residues R136, R57, Y133, L69, L82, and Y139 were the active sites of the PLK1 protein (PDB code: 2RKU), in which the more active ligands can inhibit the enzyme of PLK1. The results of the molecular dynamic MD simulation diagram were obtained to reinforce the previous molecular docking results, which showed that both inhibitors remained stable in the active sites of the PLK1 protein (PDB code: 2RKU) for 50 ns. Finally, a check of the ADME-Tox properties of the two most active molecules showed that molecular N° 28 could represent a good drug candidate for the therapy of prostate cancer diseases.


Introduction
Among the most infectious and dangerous diseases in recent years is cancer disease, which is manifested by the abnormal regulation of different pathways; it remains the second most prevalent disease and one of the major health problems in the world [1]. The research into anti-cancer drugs targeting key factors is essential to regulate this major problem [2]. Among the methods used to discover drug candidates are computational methods to reduce the cost of drugs and increase their effectiveness [3]. Polo-like kinases (PLKs), a serine-threonine kinase, have five family members (PLK1-5) that play a key role in mitosis and have been proven to be necessary for centrosome maturation and bipolar spindle formation [4]. PLK1 is the most investigated of all PLK family members. PLK1 overexpression has been found in many types of different cancers (lung cancer, prostate

Database and Biological Activity
To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model.

Database and Biological Activity
To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC 50 ), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model.

Database and Biological Activity
To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model.

Database and Biological Activity
To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model.

Database and Biological Activity
To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model.

Database and Biological Activity
To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model. To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model. To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model. To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model. To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model. To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model. To build the different CoMFA and CoMSIA models, we based them on a set of experimental data (28 derivatives, Table 1) synthesized by Zhiwei Li et al. [6]. To evaluate their anti-cancer biological activity (IC50), we divided the experimental data set (80% training set (22 derivatives) and 20% test set (6 derivatives)) [7] and the training set to construct a model and used the test set to evaluate the performance of the built model.

Molecular Alignment and Generation of the Models
First of all, it was necessary to align all the molecules because molecular alignment is one of the most important steps for the generation of the CoMFA and CoMSIA models [14]. In the database where we grouped all the molecules of the new pteridinone derivatives, we performed a rigid distill alignment using SYBYL-X 2.1 software [15]. The database of all molecules was minimized by using the standardized Tripos force field and employing the Gasteiger-Huckel atomic partial charges. Furthermore, 0.005 kcal/mol Å

Molecular Alignment and Generation of the Models
First of all, it was necessary to align all the molecules because molecular alignment is one of the most important steps for the generation of the CoMFA and CoMSIA models [14]. In the database where we grouped all the molecules of the new pteridinone derivatives, we performed a rigid distill alignment using SYBYL-X 2.1 software [15]. The database of all molecules was minimized by using the standardized Tripos force field and employing the Gasteiger-Huckel atomic partial charges. Furthermore, 0.005 kcal/mol Å

Molecular Alignment and Generation of the Models
First of all, it was necessary to align all the molecules because molecular alignment is one of the most important steps for the generation of the CoMFA and CoMSIA models [14]. In the database where we grouped all the molecules of the new pteridinone derivatives, we performed a rigid distill alignment using SYBYL-X 2.1 software [15]. The database of all molecules was minimized by using the standardized Tripos force field and employing the Gasteiger-Huckel atomic partial charges. Furthermore, 0.005 kcal/mol Å

Molecular Alignment and Generation of the Models
First of all, it was necessary to align all the molecules because molecular alignment is one of the most important steps for the generation of the CoMFA and CoMSIA models [14]. In the database where we grouped all the molecules of the new pteridinone derivatives, we performed a rigid distill alignment using SYBYL-X 2.1 software [15]. The database of all molecules was minimized by using the standardized Tripos force field and employing the Gasteiger-Huckel atomic partial charges. Furthermore, 0.005 kcal/mol Å 16

Molecular Alignment and Generation of the Models
First of all, it was necessary to align all the molecules because molecular alignment is one of the most important steps for the generation of the CoMFA and CoMSIA models [14]. In the database where we grouped all the molecules of the new pteridinone derivatives, we performed a rigid distill alignment using SYBYL-X 2.1 software [15]. The database of all molecules was minimized by using the standardized Tripos force field and employing the Gasteiger-Huckel atomic partial charges. Furthermore, 0.005 kcal/mol Å

Molecular Alignment and Generation of the Models
First of all, it was necessary to align all the molecules because molecular alignment is one of the most important steps for the generation of the CoMFA and CoMSIA models [14]. In the database where we grouped all the molecules of the new pteridinone derivatives, we performed a rigid distill alignment using SYBYL-X 2.1 software [15]. The database of all molecules was minimized by using the standardized Tripos force field and employing the Gasteiger-Huckel atomic partial charges. Furthermore, 0.005 kcal/mol Å 7. 18 4.787 * Test set, pIC 50 = 6-log 10 (IC 50 ).

Molecular Alignment and Generation of the Models
First of all, it was necessary to align all the molecules because molecular alignment is one of the most important steps for the generation of the CoMFA and CoMSIA models [14].
In the database where we grouped all the molecules of the new pteridinone derivatives, we performed a rigid distill alignment using SYBYL-X 2.1 software [15]. The database of all molecules was minimized by using the standardized Tripos force field and employing the Gasteiger-Huckel atomic partial charges. Furthermore, 0.005 kcal/mol Å was set as the standard convergence parameter for the Powell gradient algorithm, and 1000 rounds of iterations were performed to produce a more stable configuration of the molecule [16]. Secondly, field studies of the CoMFA and CoMSIA descriptors were computed for each lattice, with a grid distance of 1 Å and extending to 4 Å points, in all three coordinates within the defined region [17]. Steric (S), electrostatic (E), acceptor hydrogen bonding (A), and hydrophobic (H) fields were calculated using the standard Tripos force field and a Van Der Waals potential, as well as the Coulombic terms [18]. A sp 3 hybridized carbon (C) atom with a charge of +1e was utilized as the probe atom, and steric and electrostatic field energy values were truncated at 30 kcal/mol [19]. The attenuation factor and filter factor values for the column were 0.3 and 2 kcal/mol for the steric and electrostatic fields [20].
Finally, the partial least-squares (PLS) method allowed for the correlation between the different CoMFA and CoMSIA fields, which included the value of the biological activity of pteridinone derivatives [21]. This method is usually used in 3D-QSAR studies to determine the various statistical values, such as the optimal number of components (NOC), the crossvalidated coefficients (Q 2 ), the conventional coefficient (R 2 ), F-statistic values (F), and the standard error of estimation (SEE). An accurate model means higher R 2 and Q 2 statistical parameters and the smallest possible SEE [22]. To estimate the predictive capabilities of the 3D-QSAR models, the biological activity of the external test set (5 molecules) was predicted using the resulting PLS models, and the predictive correlation coefficient (R 2 pred ), based on the molecules in the test set, must be greater than 0.6 [23]. The QSAR model is considered good when the Q 2 value is greater than 0.5 [24]. During this analysis, SAMPLES remained disabled, and column filtering was set to 2.0 kcal mol −1 to speed up the analysis [25]. The leave-one-out (LOO) method involved removing one of the individual compounds from the training set and then predicting the activity of each removed compound to verify that the results of the CoMFA and CoMSIA models were predictive for the compounds that were not in the training set [26].

Molecular Docking
The objective of molecular docking is to give a prediction of the structure of the ligand-receptor complex using computational methods. Docking can be performed in two interdependent steps: first, sampling the conformations of the ligand in the active site of the protein; then, classifying these conformations via a scoring function. The software we used to perform the molecular docking, such as that discovery in 2021 to delete molecules of water and view the ligand/protein interaction [27], was Auto Dock Tools 1.5.6 [28], and we used Vina to execute molecular docking [29]. We used the co-crystalline structure BI-2536/PLK1 (PDB code: 2RKU, with a resolution of 1.95 Å) [13] as a docking model and performed molecular simulation docking of the most active molecules with PLK1 [30,31]. Key amino acid residues that facilitated docked ligand binding to the PLK1 active site were R136, R57, Y133, L69, L80, H138, Y67, Y82, and E132. After the preparation of the ligand and the receptor, we effected molecular docking on the two most active molecules, N • 17 and 28.

Molecular Dynamic (MD)
Based on the molecular docking results, the two best-docked ligands with the highest activity were chosen for molecular dynamics (MD) simulations to identify the molecular recognition between the ligand and the protein [32]. MDs were performed for 50 nanoseconds using GROMACS 5.0 software and the GROMOS9643a1 force field [33]. The SOC water model was chosen to simulate the MD in explicit solvation [34]. The ligand loading parameters were obtained from the Dundee prodrg 2.5 servers [35]. Other input parameters were selected, such as the SOC water model, the type of triclinic box, and the type of salt (Na + , Cl − ), to be neutralized [36]. The system was equilibrated at a temperature of 300 K and a pressure of 1 bar, with canonical NVT and isobaric NPT sets, respectively [37]. The MD simulations were performed for a time of 50 ns, with the temperature and pressure stable, a time step of 2 fs, and a long-range interaction threshold of 1 nm [38].

Synthetic Accessibility and ADMET Prediction
After the stability study of the two most active molecules (molecule N • 17 and 28), it was necessary to study their pharmacokinetic and pharmacodynamics properties. First, it was necessary to verify the synthetic accessibility of these compounds, and second, it was necessary to study the pharmaceutical properties of each compound. The pkCSM [39] and SwissADME servers [40] were used to evaluate the synthetic accessibility and pharmaceutical properties of the proposed compounds [41].

Distill Rigid Alignment
The structural alignment of molecules is critical to both the predictive accuracy of the 3D-QSAR model and the reliability of the contour models. The database of pteridinonederived molecules was aligned on the common core using the distill rigid alignment technique in Sybyl X-2.1. Molecular N • 28 (most active) was considered a template. Based on Figure 2, all 3D molecular structures were superimposed (Figure 2A) on the common core ( Figure 2B). The structural alignment of molecules is critical to both the predictive accuracy of the 3D-QSAR model and the reliability of the contour models. The database of pteridinone-derived molecules was aligned on the common core using the distill rigid alignment technique in Sybyl X-2.1. Molecular N° 28 (most active) was considered a template. Based on Figure 2, all 3D molecular structures were superimposed (Figure 2A) on the common core ( Figure 2B).

Generation of the CoMFA and CoMSIA Models
The Table 2 shows the predicted and observed activity values of the different CoMFA and CoMSIA models and their residuals.

Generation of the CoMFA and CoMSIA Models
The Table 2 shows the predicted and observed activity values of the different CoMFA and CoMSIA models and their residuals. The 3D-QSAR models were suggested to quantitatively explicate and predict the different effects of the molecules on the anti-proliferative activities of a series of twentyeight pteridinones. The variables in the training set were put into PLS cross-validation analysis to identify the proper statistical metrics for every module. The Table 3 shows the results for the different statistical modules and their own parameters. The following table demonstrates the results for the different statistical modules and their own parameters. The percentage of CoMFA model contributions of the two fields (steric and electrostatic) explained 81.4% and 18.6% of the variance, respectively, the values of Q 2 cv and R 2 were 0.67 and 0.992, respectively, the optimal number of principal components (NOC) used was 9, the F-value was 27.47, and the lower value of SEE was equal to 0.035.
In the CoMSIA study, the evaluation analysis of the three selected models had various combinations of the fields, such as stereoscopic (S), hydrophobic (H), electrostatic (E), and hydrogen bond acceptor (A). According to The results showed that steric, electrostatic, hydrogen bond acceptor, and hydrophobic fields played an important role in these models. Among these four fields, the steric and hydrophobic fields were the most important interactions between the ligand and the receptor protein.

External Validation
From the results in Table 4, it can be seen that the four models, CoMFA, CoMSIA/SEH, and CoMSIA/SEAH, had better external prediction coefficients than the three models (R 2 pred ), which were operated to validate the external predictive abilities of the three models. The R 2 pred values of the three models were 0.683, 0.758, and 0.767, respectively. Therefore, the three models had Q 2 values greater than 0.5 [23]. Regressions of the predicted pIC 50 versus the observed pIC 50 or the predicted pIC 50 versus the observed pIC 50 through the origin should be characterized by the fact that K or K' (slopes of corresponding regression lines) is close to 1 [42]. According to Table 4, the two models chosen for external validation were CoMSIA/SEH and CoMSIA/SEAH, which obtained values of K and K' of 0.923 and 1.012, respectively, so they were close to 1. Therefore, these two models (3D-QSAR) that we chose are acceptable. Then, we used them to predict the activity of a test molecule.

Analyzation of the CoMFA and CoMSIA Contour Charts
To show the details contained in the two 3D-QSAR models, we selected the more active molecule in the series (molecule N • 17) to visualize the different fields on the molecule. The different fields of the CoMSIA and CoMFA contour charts are presented in Figures 3 and 4.

Analyzation of the CoMFA and CoMSIA Contour Charts
To show the details contained in the two 3D-QSAR models, we selected the more active molecule in the series (molecule N° 17) to visualize the different fields on the molecule. The different fields of the CoMSIA and CoMFA contour charts are presented in Figures 3 and 4.

Analyzation of the CoMFA and CoMSIA Contour Charts
To show the details contained in the two 3D-QSAR models, we selected the more active molecule in the series (molecule N° 17) to visualize the different fields on the molecule. The different fields of the CoMSIA and CoMFA contour charts are presented in Figures 3 and 4.

CoMFA Contour Chart
In the CoMFA model of the two fields, stereoscopic and electrostatic are presented in Figure 3.
The two figures represent the two electrostatic and stereoscopic contour plots from the CoMFA study that can help provide information about the regions that can decrease or increase the biological activity of pteridinone derivatives. At the steric field level, only a few green outline parts are located near the meta-position of the benzene ring, which means that the para-position increases the inhibitory activity of the more active compound. On the other hand, a large yellow contour portion that decreases the inhibitory activity of the more active compound originates from the pteridinone molecule substituents. In the electrostatic field, the blue contours indicate that it is positively charged and favors the inhibitory activity of this ligand, while the red contours indicate that this region is negatively charged and disadvantages the inhibitory activity of this ligand. After studying the figure, we found that both fields were farther from the ligand because of the ligand that has no charge.

CoMFA Contour Chart
In the CoMFA model of the two fields, stereoscopic and electrostatic are presented in Figure 3.
The two figures represent the two electrostatic and stereoscopic contour plots from the CoMFA study that can help provide information about the regions that can decrease or increase the biological activity of pteridinone derivatives. At the steric field level, only a few green outline parts are located near the meta-position of the benzene ring, which means that the para-position increases the inhibitory activity of the more active compound. On the other hand, a large yellow contour portion that decreases the inhibitory activity of the more active compound originates from the pteridinone molecule substituents. In the electrostatic field, the blue contours indicate that it is positively charged and favors the inhibitory activity of this ligand, while the red contours indicate that this region is negatively charged and disadvantages the inhibitory activity of this ligand. After studying Life 2023, 13, 127 9 of 15 the figure, we found that both fields were farther from the ligand because of the ligand that has no charge.

CoMSIA/SEA Contour Chart
In the CoMSIA/SEA model, four different fields, stereoscopic, electrostatic, hydrogen bond acceptor, and hydrophobic fields, were checked. The four contour charts are presented in Figure 4.
The CoMSIA steric contour chart (Figure 4a) of the most active molecule presented similar results to the CoMFA contour chart. However, at the electrostatic contour chart levels (Figure 4b), we observed a large blue contour covering the R substitutions, indicating that the selection of large substitution groups (CF3) is necessary for this region to increase the inhibitory activity. We also observed a small red contour. At the acceptor hydrogen bonding contour chart level (Figure 4c), we observed a large purple contour covering the meta-position of R substitutions, indicating that the selection of bulky substitution groups at the meta-position is necessary for this region to increase inhibitory activity. We also observed a small orange contour away from the R substituent. At the hydrophobic binding contour map level (Figure 4d), we observed a small cyan contour covering the meta-position of R substitutions, indicating that the selection of bulky substitution groups at the meta-position is necessary for this region to increase inhibitory activity. We observed almost no magenta contour around the R substituent.
Therefore, it can be concluded that the meta-position of substituent R allows for an increase in the inhibitory activity of the chosen molecule if this position has an attractive and bulky group.

Molecular Docking
To explain and understand the interaction between the most active molecule and its protein, we performed a molecular docking of the two most active molecules. The results found by molecular docking for the two selected compounds (molecule N • 17 and 28) are presented in Figure 5.
The first visual surveillance of the two results shows that there are three hydrogen bonds with residues R136, R57, and Y133, with a distance match to 2.20, 4.94, and 2.01, respectively. There is one halogen bond with residue L69, with a distance match to 3.70, and three alkyl bonds, L80, H183, and Y67, with a distance match to 4.25, 4.95, and 4.29, respectively, in molecule N • 17. In molecule N • 28, we also observed three hydrogen bonds with residues R136, R57, and Y82, with an equal distance of 2.79, 2.82, and 2.87, respectively, and two halogen bonds with residues L69 and Y139, with match distances 3.65 and 3.10, respectively. There are two alkyl bonds with residues E132 and Y67 and a single Pi-Pi bond with residue H183. Thus, these different interactions between the molecules and their protein mean that the two molecules have a greater inhibitory effect. The overall interactions are similar to the interactions observed between a co-crystallized ligand and protein PLK1.The results suggest that the docking result is reasonable and can be used for further simulations and analysis, such as the MD simulation.
Life 2023, 13, x FOR PEER REVIEW 10 of 16 found by molecular docking for the two selected compounds (molecule N° 17 and 28) are presented in Figure 5. The first visual surveillance of the two results shows that there are three hydrogen bonds with residues R136, R57, and Y133, with a distance match to 2.20, 4.94, and 2.01, respectively. There is one halogen bond with residue L69, with a distance match to 3.70, and three alkyl bonds, L80, H183, and Y67, with a distance match to 4.25, 4.95, and 4.29, respectively, in molecule N° 17. In molecule N°28, we also observed three hydrogen bonds with residues R136, R57, and Y82, with an equal distance of 2.79, 2.82, and 2.87, respectively, and two halogen bonds with residues L69 and Y139, with match distances 3.65 and 3.10, respectively. There are two alkyl bonds with residues E132 and Y67 and a single Pi-Pi bond with residue H183. Thus, these different interactions between the molecules and their protein mean that the two molecules have a greater inhibitory effect. The overall interactions are similar to the interactions observed between a co-crystallized ligand and protein PLK1.The results suggest that the docking result is reasonable and can be used for further simulations and analysis, such as the MD simulation.

Molecular Dynamics Simulation
The two most active ligands (molecules N° 17 and 28) were chosen for dynamic molecular simulation for 50 ns to examine their stability, with respect to enzyme PLK1. The conformational changes of the two ligands are shown in the Figure 6.

Molecular Dynamics Simulation
The two most active ligands (molecules N • 17 and 28) were chosen for dynamic molecular simulation for 50 ns to examine their stability, with respect to enzyme PLK1. The conformational changes of the two ligands are shown in the Figure 6. To study the dynamic of the protein-ligand interactions of two docked molecules, a 50 ns MD simulation was run using two complexes (PLK1-molecule N° 17 (violet) and PLK1-molecule N° 28 (green)) designed for the protein to ensure the predicted binding stability of the complex system. From Figure 6, it can be seen that the RMSD of the complex PLK1-molecule N° 17 fluctuated between 0.22 and 0.30 nm during the MD simulations, and the average RMSD was found to be 0.28 nm, from 0 to 30 ns. The curve of the complex increased slightly to the value of about 0.28 nm, and then an equilibrated system was obtained in the complex for the remaining time.
The RMSD analysis indicated that for the complex PLK1-molecule N° 28, the protein fluctuated between 0.22 and 0.36 nm during the MD simulations, and the average RMSD was found to be 0.3 nm. The RMSD analysis indicated that the designed molecules N° 28 and To study the dynamic of the protein-ligand interactions of two docked molecules, a 50 ns MD simulation was run using two complexes (PLK1-molecule N • 17 (violet) and PLK1-molecule N • 28 (green)) designed for the protein to ensure the predicted binding stability of the complex system. From Figure 6, it can be seen that the RMSD of the complex PLK1-molecule N • 17 fluctuated between 0.22 and 0.30 nm during the MD simulations, and the average RMSD was found to be 0.28 nm, from 0 to 30 ns. The curve of the complex Life 2023, 13, 127 11 of 15 increased slightly to the value of about 0.28 nm, and then an equilibrated system was obtained in the complex for the remaining time.
The RMSD analysis indicated that for the complex PLK1-molecule N • 28, the protein fluctuated between 0.22 and 0.36 nm during the MD simulations, and the average RMSD was found to be 0.3 nm. The RMSD analysis indicated that the designed molecules N • 28 and N • 17 formed a stable complex with the protein throughout the simulations.
To study the impact of the binding of the two designed molecules on the internal dynamics of the target protein for 50 ns, the RMSF values were also calculated. A maximum fluctuation of 0.68 nm was noticed in the loop region of residue 330 for both molecules.
The radius of gyration (Rg) (Figure 7) represents the change in the compactness of the protein structure over time. tions, and the average RMSD was found to be 0.28 nm, from 0 to 30 ns. Th complex increased slightly to the value of about 0.28 nm, and then an equilib was obtained in the complex for the remaining time.
The RMSD analysis indicated that for the complex PLK1-molecule N° 2 fluctuated between 0.22 and 0.36 nm during the MD simulations, and the a was found to be 0.3 nm. The RMSD analysis indicated that the designed molecu N° 17 formed a stable complex with the protein throughout the simulations.
To study the impact of the binding of the two designed molecules on dynamics of the target protein for 50 ns, the RMSF values were also calcul mum fluctuation of 0.68 nm was noticed in the loop region of residue 330 ecules.
The radius of gyration (Rg) (Figure 7) represents the change in the co the protein structure over time. For the complex PLK1-molecule N° 28, dur ns, the values of the Rg varied between 1.94 nm and 2.02 nm. After this time of the simulation, the values remained reasonably stable in the range of 1.9 complex PLK1-molecule N° 17 Rg values fluctuated over time between 1.95 nm. Therefore, the Rg plot shows that there is no major change in the fold ness of the target protein after the binding of molecule N° 28, but there is a in the levels of the target protein after the binding of molecule N° 17. Therefore, the Rg plot shows that there is no major change in the folding compactness of the target protein after the binding of molecule N • 28, but there is a small change in the levels of the target protein after the binding of molecule N • 17.
After studying the MD simulation diagram, the results obtained reinforced the previous molecular docking results. Thus, the two chosen ligands formed dynamically stable interactions with their proteins during the simulation time of 50 ns, as there was no great variation in their characteristics.

Synthetic Accessibility and Lipinski Rules
The synthetic accessibility allowed for the evaluation of the ease of synthesis in the best molecules chosen (molecules N • 17 and 28). Then, we evaluated the different properties of the five Lipinski's rules, which allowed us to determine whether a biologically active chemical was probable to have the chemical and physical characteristics be orally bioavailable. Table 5 lists the different synthetic availability properties and their Lipinski properties.
The results obtained in Table 5 indicate that molecule N • 28 checked all the rules of Lipinski, Veber, and Egan. That means the selected molecule has good oral bioavailability. We also evaluated the synthetic accessibility, and from the result obtained, we concluded that this compound has easy-to-synthesize characteristics. Molecule N • 17 has characteristics that are almost the same as molecule N • 28; it only has a small difference at the molecular mass level that exceeds the norm of 500. In conclusion, we can say that molecule N • 28 has a small attrition rate when tested in clinical trials and has a higher chance of success in reaching the commercial phase.

The Various ADMET Properties
To avoid failing clinical tests due to toxicity or poor pharmacokinetics, after checking the similarity of the drugs, the two most active molecules, N • 17 and 28, were subjected to ADMET prediction ( Table 6) to verify that the pharmacokinetic (absorption, distribution, metabolism, and excretion) and pharmacodynamic (efficacy and toxicity of the drug) properties of the molecules passed the study of similarity for drugs. The water solubility of a compound reflects the solubility of the molecule in water at 25 • C [43]. The solubility in water is given in log (mol/L). If the solubility value is less than zero [44], we can say that the compound is very soluble, which means that both compounds are highly soluble in water. The intestinal absorbance values were very high, above 85%, which means that both compounds had good absorbance. Caco-2 permeability is frequently used to predict the absorption of orally administered medication; a molecule considered highly permeable, Caco-2 should give predicted values at 0.9 [45], so according to the Caco-2, the values of two selected compounds can be classified as highly permeable Caco-2. The volume of distribution (VDss) is used to analyze the distribution of drugs in different tissues in vivo. If the log (VDss) value is less than −0.15, the volume of distribution is considered relatively poor, and when the log (VDss) value is superior to 0.45, the volume of distribution is classified as relatively high [46]. Both compounds had a VDss value less than −0.15, so both compounds had a low volume of distribution for drug distribution in various tissues in vivo. In terms of the central nervous system (CNS) permeability index, molecules with values (LogPS) greater than −2 are considered to be capable of penetrating the CNS, and values of LogPS less than −3 are incapable of penetrating the CNS [47]. From the table, we can summarize that all compounds of N • 17 are capable of penetrating the CNS (values less than −3), but compound N • 28 is not capable of penetrating the CNS. CYP2D6, CYP3A4, CYP1A2, CYP2C19, and CYP2C9 are isoforms of cytochrome P450, which is a crucial detoxification enzyme in the human body and is responsible for altering drug pharmacokinetics [48]. Cytochrome P450 isoenzymes play an important role in drug metabolism in the liver [49]. From the results obtained in Table 6, the two tested molecules did not become substrates of CYP2D6 or substrates of CYP3A4. At the same level, most of the compounds tested were not inhibitors of CYP1A2, CY2C9, and CYP2C19. Excretion refers to the process by which the body gets rid of waste/toxic products. The drug excretion process can be achieved by either the kidney and/or the liver, where drugs are eliminated in the forms of urine or bile, respectively [50]. The total clearance of the drug gives a general view on the half-life of the drug; the lower its value, the higher the half-life of the compound [51]. The two compounds had a low clearance value, which means that the half-life of these two compounds is high. Predicting the toxicity of compounds is an important part of the drug design development process. Computational estimates of toxicity are not only faster than determining toxic doses in animals, it can also help reduce the number of animal experiments [52]. To check the toxicity of two molecules, the Ames toxicity test was performed, and skin sensitization was checked. From Table 6, it can be concluded that both compounds were not toxic. Thus, based on the drug similarity studies and ADMET, we chose molecule N • 28 as the PLK1 enzyme inhibitor because it checks almost all of the similarity properties of a drug. However, it is necessary to make further studies in the area where this compound is a drug.

Conclusions
This study focuses on a series of twenty-eight novel pteridinone derivatives as PLK1 inhibitors. In the first step, we constructed three models and examined them using external and internal validation to identify the radical of the molecule drifts that had an influence on the biological activity of a molecule.
The models CoMFA (Q 2 = 0.67, R 2 = 0.992), CoMSIA/SHE (Q 2 = 0.69, R 2 = 0.974) and CoMSIA/SEAH (Q 2 = 0.66, R 2 = 0.975) models were used to study molecular modeling. The three models were satisfactory according to the results of the statistical validation (R 2 pred value of CoMFA, CoMSIA/SHE and CoMSIA/SEAH models is 0.683, 0.758, and 0.767 respectively). We used these models to predict the activity of the molecules in the test set, and then we can use these models to predict the activity of new molecules as PLK1 inhibitors for prostate cancer treatment.
A molecular docking study was performed to identify the type of binding between the most active ligand and the PLK1 inhibitor. The key amino acids affecting the activity of these inhibitors, such as R136, R57, and Y133, easily formed hydrogen bonds with the selected small molecules and a halogen bond with the L69 residue, so these different bonds could allow the PLK1 inhibitors to maintain stability in the binding site.
Thus, the two selected ligands formed dynamically stable interactions with their protein during the 50 ns simulation time. Finally, an ADMET prediction of two ligands showed that only compound N • 28 could become a good drug candidate for cancer drug development.