N-alkylated 5,5-diphenylhydantoin Derivatives: Synthesis, X-ray, Spectroscopic Characterization, Hirshfeld Surface Analysis, DFT, Molecular Docking, Molecular Dynamics Simulations, and Cholesterol Oxidase Inhibitory Estimation

1. Introduction

Penicillin was discovered in the 1920s. Since then, antibacterials have been considered a cure-all for bacterial infections and have contributed to the global reduction of infections [1]. However, the effectiveness of current antibiotics is seriously threatened by rising antimicrobial resistance. Heterocyclic compounds have been crucial in remarkable advances against many deadly diseases [1,2,3]. Therefore, the development of new methods for the synthesis of biologically active heterocyclic compounds remains a major priority in organic chemistry [3,4,5]. Owing to their broad pharmacological and chemical properties, the hydantoin scaffold and its derivatives have consistently attracted significant interest from both synthetic and biological chemists in the pursuit of novel drug discovery [6]. When appropriately substituted, derivatives of hydantoin exhibit various pharmacological activities, including antibacterial [7], anticonvulsant [8], antitumor [9], antinociceptive and anti-inflammatories [10].

Imidazolidines are important due to their unique structural features. Their backbones contain both donor and acceptor groups. The two NH functional groups within the ring structure can act as hydrogen bond donors, while the two carbonyl functional groups can act as hydrogen bond acceptors. This leads to increased structural stability and activity, potentially expanding the range of biological activities exhibited by imidazolidines [3,4,5,6]. Phenytoin (5,5-diphenylhydantoin) is a significant pharmaceutical drug in the hydantoin class. It has been used for over 70 years to manage epilepsy and cardiac arrhythmias [6]. Phenytoin and its derivatives are preferred pharmaceutical agents because they have a long shelf life. However, they need to be easily broken down in the body for optimal efficacy. Despite being stable under extreme chemical and thermal conditions, it can be hydrolyzed more quickly in biological systems when hydantoinase enzymes are present [11]. Given its widespread presence and importance as a structural element in many natural and synthetic compounds, developing new structures derived from phenytoin requires substantial work in synthetic methodology and structural design.

In line with our previous research on the synthesis, spectroscopic characterization, physicochemical, and biological properties of phenytoin derivatives [12,13,14], we report on the synthesis of novel N-alkylated phenytoin derivatives (5-13). These derivatives were generated by the alkylation reactions of 3-benzyl-5,5-diphenylimidazolidine-2,4-dione (3) with various alkyl halides under solid-liquid phase transfer catalysis (PTC). The starting compound, phenytoin 2, was obtained by base-catalyzed condensation from benzyl and urea. Phenytoin 2 was then alkylated with benzyl chloride to afford compounds 3 and 4 (Scheme 1). All synthesized compounds' structures were determined using a combination of spectroscopic techniques, including FTIR, UV-Vis, ¹H and ¹³C NMR, and LCMS. X-ray diffraction analysis confirmed the structures of compounds 4, 5 and 6. More detailed analysis of these molecules was carried out using DFT calculations and Hirshfeld surface analysis, providing important information on intermolecular interactions in the crystal lattice. Theoretical studies, such as molecular docking, were carried out to test the interactions of compounds 4, 5 and 6 with cholesterol oxidase. Long-md dynamic simulations were also conducted to study the three compounds for their dynamic behavior concerning stability.

2. Results and Discussion

2.1. X-ray Diffraction Analysis

The structures of 4, 5 and 6 determined crystallographically (Figure 1 and Table 1) confirm those deduced from their NMR spectra. In 4, the two benzyl groups are folded towards the normal to the plane of the imidazolium ring over the same face of that ring giving this portion of the molecule a pincer-like conformation (Figure 1).

The imidazolium ring is slightly ruffled with N2 and C1 disposed of respectively, 0.0275(6) and -0.0270(5) Å from the mean plane of the ring (rms deviation = 0.0206 Å). The mean planes of the C4···C9 and C10···C15 benzene rings are inclined to that of the imidazolium ring by 69.97(3) and 80.73(3)°, respectively while those of the C17···C22 and the C24···C29 benzene rings are inclined to the same plane by 64.34(3) and 72.00(3)°, respectively. With the sums of the angles about N1 and N2 being 360° within experimental error, their lone pairs are involved in N→C π bonding. This is shown particularly by the N1—C2 and N2—C3 distances which are, respectively, 1.3686(13) and 1.3543(12) Å. In the crystal, inversion-related pairs of C18—H18···Cg3 and of C23—H23A···Cg5 (Cg3 and Cg5 are the centroids, respectively, of the C10···C15 and the C24···C29 benzene rings) interactions form chains of molecules extending along the c-axis direction. The chains pack with normal van der Waals contacts. For 5, The imidazolium ring is somewhat ruffled as N1 and C3 are, respectively, 0.0200(6) \Å to one side of the mean plane (rms deviation = 0.0154 Å) and 0.0214(6) Å to the opposite side. The mean planes of the C4···C9, the C10···C15 and the C17···C22 benzene rings are inclined to the above plane by 79.10(4), 74.01(4) and 72.48(4)°, respectively. The sums of the interbond angles about N1 and N2 are both 360° within experimental error indicating participation of their lone pairs in N→C π bonding. This appears to occur primarily in the N1—C2 and N2—C3 bonds whose lengths are 1.3582(13) and 1.3542(13) Å, respectively. In the crystal, C23—H23C···O1 hydrogen bonds form chains of molecules extending along the b-axis direction which are connected in pairs by C16—H16A···O1 hydrogen bonds and C18—H18···Cg1 (Cg1 is the centroid of the N1/N2/C1/C2/C3 ring) interactions. These ribbons are then linked by C7—H7···O2 hydrogen bonds and C8—H8···Cg4 (Cg4 is the centroid of the C10···C15 benzene ring) interactions to form layers of molecules parallel to (10$\overline{1}$). The molecule of 6 adopts a cup-shaped conformation with the imidazolidine ring as the base and the benzyl, the ethyl and one phenyl group as the sides. Partially responsible for this conformation is the intramolecular C24—H24B···Cg3 (Cg3 is the centroid of the C10···C15 benzene ring) interaction (Figure 1). The C17···C22 benzene ring it tilted away from the normal to the plane of the imidazolidine ring (dihedral angle between the mean planes of the two rings is 68.21(5)°) but the C2—N1—C16—C17 torsion angle is 86.39(14)° so that the rotational orientation of the benzyl group is nearly perpendicular to the plane of the 5-membered ring. A similar orientation is found for the ethyl group as the C1—N2—C23—C24 torsion angle is -98.31(13)°. The mean planes of the C4···C9 and C10···C15 benzene rings are inclined to that of the imidazolidine ring by 68.47(4) and 81.88(4)°, respectively, while the dihedral angle between them is 78.12(3)°. All bond distances and interbond angles appear as expected for the formulation given. In the crystal, C23—H23B···Cg2 (Cg2 is the centroid of the C4···C9 benzene ring) interactions and C9---H9...O1 hydrogen bonds form chains of molecules extending along the c-axis direction. The chains pack with normal van der Waals contacts.

2.2. Optimized Geometries

Figure 2 shows the optimized structures of 4, 5 and 6 obtained at the DFT / B3LYP / 6-311++G(d,p) level. Table 2 presents the pertinent optimized structural parameters of 4, 5 and 6, including bond lengths (Å) and angles (°) with comparisons with their experimental counterparts. Generally, the theoretical values align closely with the experimental ones with the largest differences occurring with the C=O distances for the two carbonyl groups in each molecule. Here the calculated values are up to 0.1 Å longer than the experimental and considering that the remaining distances and the interbond as well as the torsion angles agree quite well suggests that the differences are not primarily due to the phase difference between the calculated molecules and their experimental counterparts but rather to some artifact of the computational method.

2.3. DFT Results

Its frontier molecular orbitals shape the molecular interactions of a compound. Specifically, the LUMO (lowest unoccupied molecular orbital) serves as a potential electron acceptor while the HOMO (highest occupied molecular orbital) typically acts as an electron donor. As a result, the energy level of the HOMO is linked to the ionization potential of the molecule, while the energy of the LUMO corresponds to its electron affinity see Figure 3 and Table 3. The energy gap between the HOMO and LUMO influences the structural stability of a molecule. Molecules with a narrow gap, commonly termed "soft" molecules, often exhibit increased polarization and are expected to have significant chemical reactivity and thus low kinetic stability [15,16,17]. The gap energies calculated for 4-6 are essentially the same and one can predict them to have quite similar reactivities. In molecules, 4, 5 and 6, benzyl, methyl, and ethyl groups differ at N3, respectively. One may conclude that the N-alkylation has weak influences on the electronic and molecular properties.

2.4. Hirshfeld Surface Analyses

The Hirshfeld surfaces of 4- 6 plotted over dnorm are displayed in Figure 4. The shape index of the HS serves as a tool to visualize π···π stacking through the appearance of adjacent red and blue triangles. indicates the absence of π ... π interactions. The absence of adjacent red and/or blue triangles in Figure S1) implies there are no π···π interactions in 4, 5 and 6. The overall two-dimensional fingerprint plots, Figure S2(a—c) and those delineated into H···H, O···O, O···H/H···O, C···C, N···H/H···N, C···O/O···C and C···H/H···C (4), C···C, C···H/H···C, N···H/H···N, O···H/H···O, N···O/O···N, C···O/O···C, O···O and H···H, (5) and H···H, O···H/H···O, C···C, N···H/H···N and C···H/H···C (6) contacts [18] are illustrated in Figure S2 (b—h) (4), Figure S2(b—i) (5), Figure S2 (b—f) (6), together with their relative contributions to the HS. The most important interactions are H···H contributing 56.4% (4), 56.1% (5) and 62.4% (for 6), to the overall crystal packings, which are reflected in Figure 5 (a—c) with the tips at de = di = 1.09 Å (4), de = di = 1.13 Å (for 5) and de + di = 1.22 Å (for 6). In the presence of 4, 5 and 6 for C—H···π interactions, the C···H/H···C contacts contributed 30.4% (4), 28.8% (5) and 22.9% (6) to the overall crystal packings, are reflected in Figure S2 (b) with the tips of the wings at de + di = 2.74 Å (4), de + di = 2.65 Å (5) and de + di = 2.67 Å (6). The symmetric pairs of spikes for the O···H/H···O contacts contributing 9.8% (for 4), 12.9% (5) and 11.90% (for 6) to the overall crystal packings, are reflected in Figure S2 (d) with the tips at de + di = 2.49 Å (4), de + di = 2.35 Å (5) and de + di = 2.41 Å (6). The C···C contacts [ Figure S2 (e) appearing as bullet-shaped distributions of points have contributions of 1.8% (4) and 1.7% (6) to the HSs with the tips at de = di = 1.76 Å (4) and de = di = 1.79 Å (6). The Hirshfeld surface representations with the function dnorm plotted onto the surface are illustrated for the C···H/H···C, H···H and O···H/H···O interactions of 4, 5 and 6 in Figure 5. The HS analyses affirm the significance of H-atom contacts in shaping the crystal packings. The prevalence of C···H/H···C, O···H/H···O and H···H interactions suggests that hydrogen bonding and van der Waals interactions are crucial factors in the crystal packings [19].

2.5. Molecular Docking Studies

As mentioned in the methodology section, the estimated cholesterol oxidase inhibition of 4, 5 and 6 is investigated by determining their binding affinities into the binding site of cholesterol oxidase. All three fit well into the binding site of cholesterol oxidase, and they established stable complexes with its amino acids with binding energies of -10.97, -8.94, and -9.66 kcal/mol. The negative binding energies indicate that cholesterol oxidase inhibition is thermodynamically favorable and spontaneous. Figure 6 displays 2D and 3D binding interactions of 4, 5 and 6 with the amino acids of cholesterol oxidase. Due to the high binding affinity in terms of binding energy and binding interactions of the docked compounds into the binding site of cholesterol oxidase, one may estimate that 4 may have a higher inhibition efficiency towards the cholesterol oxidase than 5 or 6 (Figure 6). Compounds 5 and 6, alkylated with methyl and ethyl groups, exhibited close binding energies and interactions toward cholesterol oxidase. Hence, one concludes that these compounds may have comparable potency to inhibit cholesterol oxidase. The estimated cholesterol oxidase inhibition efficiency of 4, 5 and 6 is lower than the estimated one of the original docked ligands flavin-adenine dinucleotide, which leads to the formation of a stable complex of relative binding energy to 4, 5 and 6 compounds of 6 kcal-mol^-1 approximately.

2.6. Molecular Docking Studies MD Simulation

2.6.1. MD Simulation of Protein-Ligand Complex

An 80 ns MD simulation was performed on the complex of cholesterol oxidase with 4. The objective was to evaluate the stability of the complex and delineate the interactions between the protein and ligand throughout the simulation. Stability was assessed by examining variations in root mean square deviation (RMSD), RMSF, secondary structure elements (SSE) in cholesterol oxidase (including alpha-helices and beta-sheets), radius of gyration (Rg), SASA, and the formation of hydrogen bonds [20]. Figure 7 illustrates the RMSD, RMSF, SSE, Rg, SASA and hydrogen bond plots derived from the molecular dynamic simulation. RMSD measures the average displacement of atoms within the simulated structure relative to a reference structure. Initially, all 800 frames of the protein were aligned to the reference frame's backbone (frame 0), after which RMSD was computed using the C-alpha atoms. RMSD offers insights into the degree of conformational changes experienced by the molecule during the simulation [20]. A heightened RMSD value signifies instability in the protein system. Examination of the RMSD plot (Figure 7a) disclosed an initial sharp increase in both protein and ligand RMSDs at the onset of the simulation, followed by stabilization with consistent fluctuation. After five nanoseconds of simulation, the protein-ligand system reached a state of equilibrium and maintained a consistent RMSD for the remainder of the simulation. The average RMSD for the protein and ligand post-equilibration was determined to be 1.24 ± 0.07 Å and 4.46 ± 0.51 Å, respectively. An RMSD value within the range of 1 to 3 Å is considered acceptable for proteins, suggesting no significant alterations in the protein's structure compared to the reference frame. However, concerning the ligand, significant conformational changes occurred during the simulation, as indicated by a notable increase in RMSD. The rise in ligand RMSD is attributed mainly to the formation and collapse of bonds responsible for anchoring the ligand structure to the protein.

Similar to the RMSD plot, the RMSF plot (Figure 7b) depicts the deviations observed in the residues of the protein throughout the simulation. The RMSF plot offers valuable insights into how individual atoms behave within a molecular system throughout the simulation. A specific atom or residue with a higher RMSF value indicates greater fluctuation. This data is pivotal for discerning which areas of a protein undergo notable structural changes or display heightened flexibility during the simulation [21]. The RMSF plot highlighted notable fluctuations in specific residues, including TYR436, ASN79, THR5, ALA7, GLY9, GLY82, GLY393, GLY258, ASP8, ASP10, and THR394, falling within the acceptable range of 1.5 to 2.4 Å. Among these residues, GLY393 and THR394 exhibited ligand contact. Figure 7c presents the analysis of protein SSE, providing insight into each residue's percentage of SSE contributed. The examination reveals that alpha-helices and beta sheets constitute 20.87% and 17.79%, respectively, totaling 38.66% of the observed SSE. Notably, the SSE had no notable structural alterations throughout the simulation, suggesting a reliable binding of 4 with cholesterol oxidase.

The radius of gyration (Rg) serves as a measure of the compactness or spatial extent of a large molecule such as a protein [21]. In Figure 7d, the Rg plot for cholesterol oxidase, both unbound and complexed with 4, is presented over an 80 ns molecular dynamics (MD) simulation. For the unbound protein, the Rg initially measured 22.14 Å, steadily increasing to a peak of 22.5 Å at 50 ns, followed by a gradual decrease to 22.3 Å at 80 ns. In contrast, for the cholesterol- 4 complex, the Rg started at 22.1 Å, increased to 22.51 Å at 50 ns, and then steadily decreased to 22.2 Å by the end of the simulation. In the context of a protein, a smaller Rg value signifies a more compact and folded structure, whereas a larger Rg value indicates a more extended or unfolded conformation. The fluctuations in observed Rg (ranging from 22.1 to 22.5 Å) for both the native and ligand-bound proteins suggest that the ligand remained firmly attached throughout the simulation, contributing to the stability and compactness of the protein [22].

The SASA is another crucial quantitative metric used to evaluate the interaction between a protein and a ligand in molecular dynamics (MD) simulations [23]. SASA quantifies the accessible surface area of a biomolecule, such as a protein, to solvent molecules. This calculation considers the molecular surface and determines the accessible surface area based on the effective radius of a water molecule. SASA serves as a valuable tool for investigating biomolecular interactions and their environmental context, offering insights into the surface availability for interactions with solvent molecules, ligands, or other molecular entities [24]. In Figure 7e, the SASA plot illustrates the dynamics of native cholesterol and cholesterol bound to 4. For the native protein, SASA exhibited fluctuations within the range of 8827–9934 Å2 throughout the 80 ns simulation. The average observed solvent-accessible surface area was 9310 ± 159 Å2. Upon ligand binding to the protein, the average SASA value experienced a slight decrease to 9197 ± 74 Å2. Furthermore, fluctuations within the range of 8931–9306 Å2 were observed, mirroring the SASA fluctuation pattern seen in the native protein. This marginal reduction in the average SASA value is attributed to the ligand occupying the binding site in cholesterol oxidase.

In Figure 7f, the hydrogen bond formation between cholesterol oxidase and 4 during the MD simulation is presented. This analysis provides insights into the dynamic nature of the protein-ligand complex, revealing key interactions that contribute to the stability or transient nature of the binding [25]. On average, 1.5 ± 0.66 hydrogen bond formations were observed, reaching a peak of 4 hydrogen bonds at 21.7 ns. Towards the conclusion of the simulation, only one hydrogen bond persisted, indicating that most of the hydrogen bonds formed during the simulation were transient. This suggests that another type of non-bonded interaction may also play a crucial role in maintaining the stability of the protein-ligand complex. Figure 8 illustrates the dynamic movements of the protein-ligand complex during an 80 ns MD simulation. Minimal positional and conformational changes were observed in the ligand throughout the simulation. Figure 8a presents relative changes in ligand position and conformation, showcasing 3D snapshots of cholesterol-derivative complexes taken during the simulation, specifically highlighting significant shifts in the ligand's conformation relative to the protein. Compared to the reference snapshot at 1 ns, notable shifts in the ligand's conformation were observed at intervals of 20.2, 28.4, 52.2, 63.3, 75, and 80 ns. These changes primarily arise from rotational movements enabled by the four aromatic rings connected to the ligand's core structure via flexible, rotatable bonds. At the onset of the simulation, a minor rotational adjustment towards the left side was detected in the ligand at 20.2 ns, which was subsequently reversed by 28.4 ns. Following this, another rotational shift occurred towards the left side of the ligand at 52.2 ns. Later, there was a reversal in the ligand's rotation direction at 63.3 ns, with the only noticeable alterations thereafter being slight rotations of two aromatic rings located at the ligand's base during 75 and 80 ns. The changes observed in the ligand’s conformation correspond closely to the RMSD peaks noted throughout the simulation. Figure 8b displays the ligand's position at the binding site in cholesterol oxidase at 1 ns and 80 ns. A 90o angular rotation along the central axis is evident, with one aromatic ring remaining attached to the original binding cleft (residues between 250 – 310) and two attaching to a cleft adjacent to the original site. Given the binding site's capacity for accommodating only three aromatic rings at a time, one of the aromatic rings is consistently observed to be unattached during the simulation.

2.6.2. PCA, DCCM and FEL Analysis

Results of DCCM, PCA and FEL analysis are depicted in Figure 9. Major component analysis is one of the most used methods to study the intricacy of the dynamics of biomolecular systems during MD simulation [26]. To perform PCA, the atomic coordinates representing protein and ligand structure were collected from the MD simulation data and the structures were aligned to remove translational and rotational motions. Subsequently, a covariance matrix that quantifies the atomic motions is constructed using the aligned trajectories. The eigenvalue decomposition of this matrix yields eigenvectors and eigenvalues, where the former describes the principal components representing the major directions of motion and the latter indicates the magnitude of these motions. Projecting the original trajectory onto these eigenvectors yielded principal component trajectories, providing insights into the system's predominant modes of motion [27]. Figure 9 shows the PC1, PC2, PC3, PC9 and PC10 plots derived from atomic coordinates of cholesterol oxidase - 4 complexes with corresponding variance. The PCA trajectory analyses revealed a conformational shift within clusters, transitioning from the initial purple cluster to the intermediate green cluster and subsequently progressing to the final red color cluster. Each data point's color corresponds to the trajectory data at a specific time within the simulation. For instance, the color gradient from purple to blue signifies the period from 0 to 30 ns, while the transition from green to cyan represents the period between 30 to 50 ns. Lastly, the shift from yellow to red indicates the time frame from 50 to 80 ns. In terms of internal motion of trajectories, the protein-ligand complex with the bound ligand exhibited the highest variability in the PC1 cluster, reaching 18.29 %. However, in subsequent analysis, PC2 and PC3 statistics demonstrated a notable reduction in variability to 7.48 % and 4.05 %, respectively. The convergence observed in principal components across PC2 and PC3 suggests a minimized atomic motion variability (periodic movement) within the protein-ligand complex during simulation [28]. Furthermore, a full convergence of principal components into a singular cluster was observed across PC9 and PC10. This finding supports the notion of a stable binding interaction between 4 and the cholesterol oxidase enzyme, aligning well with the snapshots collected during the simulation period.

DCCM is yet another method that serves as a valuable approach for characterizing the dynamics of biomolecular structures. Derived from a covariance matrix that incorporates fluctuations in atomic positions, DCCM precisely quantifies correlations between pairs of atoms, providing a detailed account of the dynamic interplay within the biomolecular system during the simulation [29,30]. Figure 9b displays the dynamic cross-correlation map of the cholesterol oxidase–4 complex throughout the MD simulation. In this representation, the red color signifies a positive correlation (> 0.8), while blue indicates anti-correlation (< -0.4). The intensity of the color reflects the strength of the correlation or anti-correlation observed. Positive correlations imply atoms moving in the same direction, whereas negative correlations suggest motions in opposite directions. The DCCM map reveals a prevalence of highly correlated residues over anti-correlated residues. Notably, residues within the ligand binding site (residue index range 250-310) exhibit substantial correlation, indicating a consistent motion pattern of these binding site residues upon binding with 4.

Lastly, the free energy landscape of the protein-ligand system was constructed using the trajectory data derived from the MD simulation. The FEL analysis in MD simulations serves as a fundamental tool for unraveling the energetics and exploring the conformational space of biomolecular systems [31]. Figure 10 shows the free energy landscape of the cholesterol oxidase –4 complex with Gibb’s free energy measured against the principal components 1 and 2. The figure depicts several energy wells, with particular emphasis on two of them. A well located in the bottom right quadrant demonstrates a global minimum energy of 0 kJ/mol, while another situated in the middle-left quadrant displays a global minimum energy of 0.192 kJ/mol. These low-energy states correspond to the protein-ligand complex at 53.5 and 77.1 ns, respectively. At 53.5 ns, numerous hydrophobic and ionic interactions between the ligand and the protein were detected. Additionally, the ligand established two hydrogen bonds with water molecules. At 77.1 ns, the ligand formed a hydrogen bond (with GLN304), a water bridge (with VAL250), a pi-cation interaction (with LYS296), and various hydrophobic interactions with the protein. The free energy landscape of the cholesterol oxidase–4 complex reveals the achievement of a highly stable, lowest energy configuration, suggesting a well-folded protein structure following the attachment of 4 to the protein [32].

2.6.3. Protein-Ligand Contacts During the Simulation and MM-GBSA Binding Energies

A summary of bonded interactions established between cholesterol oxidase and 4 is illustrated in Figure 11. The histogram depicting the fraction of different interactions made by the residues concerning ligands is shown in Figure 11a. Four distinct interactions were documented during the simulation, encompassing hydrogen bonds, water bridges, ionic bonds, and hydrophobic interactions. Of these relations, the hydrogen bond stands out as the most crucial in terms of affording stability for the attachment of the ligand to the protein. Moreover, hydrogen bonds appear to play a pivotal role in influencing the drug's specificity [33]. Throughout the MD simulation, two hydrogen bonds were noted between cholesterol oxidase and 4. These hydrogen bonds were established by the residues GLN304 and SER300 with one of the oxygen atoms attached to the imidazole group in 4. Among these hydrogen bonds, the bond formed by GLN304 persisted for more than 30% of the time during the simulation (Figure 11b).

During the simulation, hydrophobic interactions like π-π and π-cation interactions were detected, involving residues HIS306 (π-π), LYS296 (π-cation), and ARG249 (π-cation). Additionally, non-specific hydrophobic interactions were formed by LEU297, MET301, and PHE390. Water bridges involving residues MET42, ARG249, VAL250, GLU266, SER300, GLN304, and GLY395 were also noted. No ionic interactions were detected throughout the simulation. Figure 11c displays superimposed images of the protein and ligand structure taken at 10 ns intervals during the simulation. While the protein structure maintains its relative position, distortions in the ligand structure are evident, corroborating the RMSD results. Despite these distortions, the ligand consistently occupies the same cleft or binding region in the protein. Moreover, a comparison of bonded interactions between the initial and final frames of the simulation indicates the preservation of many original interactions (enlarged image in Figure 11c). Post-simulation analysis of interactions between the protein and ligand (Figure 11d) indicates the presence of one hydrogen bond (formed by GLN300) and numerous hydrophobic interactions that endure until the end of the simulation. This suggests that the complex formed by cholesterol oxidase and 4 is highly stable and characterized by numerous non-bonded interactions. Furthermore, the computed MM-GBSA binding energies for the protein-ligand complex highlight the significant binding affinity of 4 towards cholesterol oxidase, with a total MM-GBSA binding energy (ΔGbind) of -51.5 ± 10.6 kcal/mol.

Figure 12 depicts the varied binding free energies calculated for the protein-ligand complex using the MM-GBSA method. Notably, the energies associated with van der Waals and electrostatic interactions contribute the most to the overall binding energy, underscoring their pivotal role in sustaining contact between the protein and ligand. In summary, the results from MD simulation and MM-GBSA analysis collectively suggest that 4 holds promise as a potential inhibitor of cholesterol oxidase.

3. Materials and Methods

3.1. Spectroscopic Characterization

The structures of the products were deduced from their nuclear magnetic resonance (NMR) spectra which were obtained on a Bruker Advance II 300 Fourier transform spectrometer, operating at 300 MHz for 1H nuclei and 75 MHz for 13C nuclei, with tetramethylsilane (TMS) used as the internal reference. In addition, liquid chromatography-mass spectrometry (LCMS) was conducted using an Acquity UPLC-BEH C18 in EIS mode, with information obtained in positive mode. IR spectra were obtained using a Nicolet IS50 instrument with potassium bromide (KBr) beam splitters, measuring frequencies between 450 and 4000 cm-1. UV/VIS spectra were recorded on a s Biobased Bk-D580 Double Beam Scanning UV/Vis spectrophotometer. Melting points were determined using a Digital Melting Point Apparatus.

2.2. Preparation of N-alkylated Phenytoin Derivatives (3-13)

5,5-Diphenylimidazolidine-2,4-dione 2 was prepared using the Biltz method [34] which involves a base-catalyzed condensation between benzyl and urea. The phenytoin, 2, was alkylated with benzyl chloride to obtain compounds 3 and 4. Compound 3 was then alkylated under solid-liquid phase transfer catalysis (PTC) in DMF containing tetrabutylammonium bromide and K₂CO₃ with the alkylating agents iodomethane, bromoethane, 2-bromobutane, 1-bromooctane, 1-bromononane, 1-bromodecane, 1-bromododecane, 1-bromohexadecane, and 1-bromooctadecane. This led to the formation of new phenytoin derivatives (5-13) in high yields, as shown in Scheme 1.

Scheme 1. Syntheses of phenytoin derivatives.

Scheme 1. Syntheses of phenytoin derivatives.

3.2.1. Synthesis Procedure of 2.

Benzyl (9.51 mmol, 2 g) and urea (15.98 mmol, 0.96 g) were dissolved in 50 mL of ethanol by stirring for 10 minutes in a suitably equipped 200 mL flask after which 5 mL of a 10 mol/L aqueous potassium hydroxide solution was added and the mixture was heated at 78 °C for 2 hours. The reaction mixture was then cooled to room temperature and 150 mL of ice water was added. The pale beige solid byproduct was removed by filtration and the transparent solution was acidified with 6 M hydrochloric acid full precipitation of phenytoin occurred. This was collected by vacuum filtration and purified by recrystallization from ethanol.

5,5-diphenylimidazolidine-2,4-dione (2)

Yield: 80%; Colorless crystals (recrystallized from ethanol); MP = 569 K; ¹H NMR (DMSO-d₆ - 300.13 MHz) δ (ppm): 9.33, 11.13 (s, 1H, NH); 7.34 -7.45 (m, 10H, HAr); ¹³C NMR (DMSO d6 - 75.47MHz) δ (ppm): 156.49, 140.27, 70.74 (Cq); 174.33, 162.81 (C=O); 127.01-129.03 (CHAr); LCMS (ESI): 253.09 [M+H+]; (Figures S3-S5).

3.2.2. Alkylation of Phenytoin with Benzyl Chloride

A mixture of 2 (1.98 mmol, 0.5 g) and an equimolar quantity of benzyl chloride in 20 mL of DMF, along with potassium carbonate (2.97 mmol, 0.41 g) and a catalytic quantity of tetra-n-butylammonium bromide (0.2 mmol, 0.06 g), was stirred at room temperature for 48 hours. The progress was monitored by thin-layer chromatography. Upon completion, the solvent was evaporated under reduced pressure, and the salts were eliminated via liquid extraction with water and dichloromethane. The two products, 3 and 4 were separated by chromatography column using hexane/ethyl acetate (5: 1) as eluent.

3-benzyl-5,5-diphenylimidazolidine-2,4-dione 3:

Yield: 80%; Colorless crystals (recrystallized from methanol); MP = 480 K; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm): 4.74 (s, 2H, CH2); 6.74 (s, 1H, NH), 7.26-7.41 (m, 15H, HAr); ¹³C NMR (CDCl₃ - 75.47MHz) δ (ppm): 173.63, 156.39 (C=O); 139.10, 135.85, 70.21 (Cq); 126.84-128.83 (CHAr); 42.65 (CH2); FTIR (KBr, ν cm -1): 3274 (NH); 1760, 1704 (C=O); 1494, 1447 (C=C) Ar; 1139 (C-N); 750, 696 (CH) Ar; LCMS (ESI): 323.14 [M+H+]; UV–Visible λ_max(nm): 288 in dichloromethane (Figures S6-S10).

1,3-dibenzyl-5,5-diphenylimidazolidine-2,4-dione 4:

Yield: 8%; White crystals (recrystallized from dichloromethane-hexane); MP = 395 K; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm): 4.81, 4.58 (S, 2H, CH2); 6.79-7.34 (m, 20H, HAr); ¹³C NMR (CDCl₃ -75.47MHz) δ (ppm): 173.55, 156.33 (C=O); 136.69, 136.57, 136.12, 75.52 (Cq); 126.99-128.70 (CHAr); 45.42, 42.99 (CH2); FTIR (KBr, ν cm -1): 1763, 1704 (C=O); 1494, 1447 (C=C) Ar; 1135 (C-N); 751, 693 (C-H) Ar; LCMS (ESI): 433.19 [M+H+]; UV–Visible λ_max(nm): 288 in dichloromethane (Figures S11-S15).

3.2.3. Synthesis of Compounds (5-13):

A mixture of 3 (0.87 mmol, 0.3 g), tetrabutylammonium bromide (0.15 mmol, 0.05 g), and potassium carbonate (1.73 mmol, 0.24 g) in 15 mL of DMF was stirred at room temperature for 30 minutes. Following this, 1.5 equivalents of the alkyl halide (iodomethane, bromoethane, 2-bromobutane, 1-bromooctane, 1-bromononane, 1-bromodecane, 1-bromododecane, 1-bromohexadecane, and 1-bromooctadecane) were added. The mixture was agitated at room temperature for 48 hours and the progress of the conversion was monitored by thin-layer chromatography. Upon completion, the solvent was evaporated under reduced pressure and the salts were extracted with water and dichloromethane.

3-benzyl-1-methyl-5,5-diphenylimidazolidine-2,4-dione (5):

Yield: 95%; White crystals (recrystallized from dichloromethane-hexane); MP = 397 K; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm): 2.81 (S, 3H, CH3); 4.78 (s, 2H, CH2); 7.19-7.41 (m, 15H, HAr); ¹³C NMR (CDCl₃ -75.47MHz) δ (ppm): 173.44, 155.91 (C=O); 136.41, 136.15, 74.42 (Cq); 127.92-128.83 (CHAr); 43 (CH2-N); 26.68 (CH3); FTIR (KBr, ν cm-1): 1777, 1703 (C=O); 1489, 1445 (C=C)Ar; 1142 (C-N); 743, 696 (C-H)Ar; LCMS (ESI): 357.16 [M+H+]; UV–Visible λ _max (nm): 236 in dichloromethane (Figures S16-S20).

3-benzyl-1-ethyl-5,5-diphenylimidazolidine-2,4-dione (6):

Yield: 93%; White crystals (recrystallized from dichloromethane-hexane). MP = 140-143 °C; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm): 4.78 (s, 2H, CH2); 7.24-7.42 (m, 15H, HAr); 3.45 (q, 2H, CH2, 3JH-H=7.2MHz); 0.64 (S, 3H, CH3, 3JH-H=7.2MHz). ¹³C NMR (CDCl3 - 75.47MHz): δ (ppm) 173.63, 155.50 (C=O); 137.3, 136.26, 74.82 (Cq); 127.85-128.86 (CHAr); 42.83, 37.02 (CH2-N); 13.37 (CH3); FTIR (KBr, ν cm-1): 2917, 2850 (C-H); 1777, 1703 (C=O); 1489, 1445 (C=C) Ar; 1142 (C-N); 743, 696 (C-H)Ar; LCMS (ESI): 371.17 [M+H+]; UV–Visible λmax(nm): 232 in dichloromethane. (Figures S21-S25).

3-benzyl-1-butyl-5,5-diphenylimidazolidine-2,4-dione (7):

Yield: 95%; Viscous oil; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm): 4.76 (s, 2H, CH2); 7.22-7.40 (m, 15H, HAr); 3.33 (m, 2H, CH2); 1.02 (qt, 2H, CH2); 0.86-0.94 (m, 2H, CH2); 0.67 (t, 3H, CH3); ¹³C NMR (CDCl₃ - 75.47MHz): δ (ppm) 173.58, 155.65 (C=O); 137.29, 136.25, 74.91 (Cq); 127.81-128.84 (CHAr); 42.81-42 (CH2-N); 29.98, 19.97 (CH2); 13.41 (CH3); FTIR (cm-1): (KBr, ν cm-1): 2917, 2850 (C-H); 1777, 1703 (C=O); 1489, 1445 (C=C)Ar; 1142 (C-N); 743, 696 (C-H)Ar; LCMS (ESI): 399.20 [M+H+]; UV–Visible λ_max(nm): 288 in dichloromethane (Figures S26-S30).

3-benzyl-1-octyl-5,5-diphenylimidazolidine-2,4-dione (8):

Yield: 87%; viscous oil; ¹H NMR (CDCl₃ - 300.13MHz): δ (ppm) 7.22-7.41 (m, 15H, HAr); 4.76 (s, 2H, CH2); 3.32 (t, 2H, CH2); 0.90-1.25 (m, 12H, CH2); 0.86 (t, 3H, CH3); 13C NMR (CDCl₃ - 75.47MHz) δ (ppm): 173.60, 155.62 (C=O); 137.30, 136.26, 74.89 (Cq); 127.81-128.84 (CHAr); 42.83-42.28 (CH2-N); 22.62-31.79 (CH2); 14.08 (CH3); FTIR (KBr, ν cm-1): 1774, 1711 (C=O); 1493, 1440 (C=C)Ar; 1139 (C-N); 740, 694 (C-H)Ar; LCMS (ESI): 455.26 [M+H+]; UV–Visible λ_max(nm): 288 in dichloromethane (Figures S31-S35).

3-benzyl-1-nonyl-5,5-diphenylimidazolidine-2,4-dione (9):

Yield: 83%; Viscous oil; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm): 4.76 (s, 2H, CH2); 7.22-7.39 (m, 15H, HAr); 3.33 (t, 2H, CH2); 1.16-1.29 (m, 14H, CH2); 1.01 (t, 3H, CH3); ¹³C NMR (CDCl₃ - 75.47MHz): δ (ppm) 173.59, 155.62 (C=O); 137.32, 136.26, 74.90 (Cq); 127.80-128.83 (CHAr); 42.83-42.28 (CH2-N); 22.62-31.79 (CH2); 14.08 (CH3); FTIR (KBr, ν cm-1): 2921, 2851 (C-H); 1768, 1711 (C=O); 1494, 1447 (C=C)Ar; 1139 (C-N); 753, 696 (C-H)Ar; LCMS (ESI): 399,2085 [M+H+]; UV–Visible λ_max(nm): 288 in dichloromethane (Figures S36-S40).

3-benzyl-1-decyl-5,5-diphenylimidazolidine-2,4-dione (10):

Yield: 91%; Viscous oil; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm): 4.77 (s, 2H, CH2); 7.40-7.24 (m, 15H, HAr); 3.34 (t, 2H, CH2); 0.92-1.32 (m, 16H, CH2); 0.9 (t, 3H, CH3); 13C NMR (CDCl3 - 75.47MHz) δ (ppm): 173.58, 155.63 (C=O); 137.34, 136.28, 74.91 (Cq); 127.76-128.85 (CHAr); 42.84-42.30 (CH2-N); 22.68-31.87 (CH2); 14.01 (CH3); FTIR (KBr, ν cm-1): 2921, 2851 (C-H); 1768, 1711 (C=O); 1494, 1447 (C=C)Ar; 1139 (C-N); 753, 696 (C-H)Ar; LCMS (ESI): 483.29 [M+H+]; UV–Visible λ _max(nm): 288 in dichloromethane (Figures S41-S45).

3-benzyl-1-dodecyl-5,5-diphenylimidazolidine-2,4-dione (11):

Yield: 89%; Viscous oil; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm) : 4.76 (s, 2H, CH2); 7.22-7.40 (m, 15H, HAr); 3.33 (t, 2H, CH2); 1.01-1.33 (m, 20H, CH2); 0.93 (t, 3H, CH3); ¹³C NMR (CDCl₃ - 75.47MHz) δ (ppm): 173.60, 155.62 (C=O); 137.30, 136.26, 74.89 (Cq); 127.81-128.84 (CHAr); 42.83-42.28 (CH2-N); 22.70-31.93 (CH2); 14.14 (CH3); LCMS (ESI): 511.33 [M+H+]; UV–Visible λ _max(nm): 288 in dichloromethane (Figures S46-S49).

3-benzyl-1-hexadecyl-5,5-diphenylimidazolidine-2,4-dione (12):

Yield: 80%; Viscous oil; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm): 4.76 (s, 2H, CH2); 3.32 (t, 2H, CH2); 0.92-1.29 (m, 28H, CH2); 0.90 (t, 3H, CH3); 7.21-7.40 (m, 15H, HAr); ¹³C NMR (CDCl₃ - 75.47MHz) δ (ppm): 173.60, 155.61 (C=O); 137.29, 136.24, 74.89 (Cq); 127.81-128.84 (CHAr); 42.83-42.27 (CH2-N); 22.70-31.94 (CH2); 14.13 (CH3); FTIR (KBr, ν cm-1): 2921, 2851 (C-H); 1768, 1711 (C=O); 1494, 1447 (C=C)Ar; 1139 (C-N); 753, 696 (C-H)Ar; LCMS (ESI): 567.39 [M+H+]; UV–Visible λ max(nm): 288 in dichloromethane (Figures S50-S54).

3-benzyl-1-octadecyl-5,5-diphenylimidazolidine-2,4-dione (13):

Yield: 86%; Viscous oil; ¹H NMR (CDCl₃ - 300.13MHz) δ (ppm): 4.77 (s, 2H, CH2); 3.34 (t, 2H, CH2); 1-1.30 (m, 32H, CH2); 0.91 (t, 3H, CH3); 7.21-7.42 (m, 15H, HAr); ¹³C NMR (CDCl₃ - 75.47MHz) δ (ppm): 173.59, 155.62 (C=O); 137.32, 136.27, 74.90 (Cq); 128-129 (CHAr); 42.83-42.28 (CH2-N); 22.72-31.95 (CH2); 14.15 (CH3); FTIR (KBr, ν cm-1): 2921, 2851 (C-H); 1768, 1711 (C=O); 1494, 1447 (C=C)Ar; 1139 (C-N); 753, 696 (C-H)Ar; UV–Visible λ _max(nm): 288 in dichloromethane (Figures S55-S58).

3.3. X-ray and Theoretical Studies

The X-ray intensity data for 4-6 were collected using a Bruker D8 VENTURE PHOTON 3 CPAD diffractometer equipped with an INCOATEC IμS-Cu microfocus source. The frames were integrated with the Bruker SAINT software package [35] using a narrow-frame algorithm, and empirical absorption corrections and merging of equivalent reflections were carried out with SADABS [36]. The structures were solved using SHELXT, and all non-hydrogen atoms were anisotropically refined on F2 by full-matrix, least-squares methods with SHELXL [37,38]. The crystallographic and experimental data for 4, 5 and 6 are shown in Table 1.

3.4. Hirshfeld Surface Analysis

The Hirshfeld surfaces of solid crystals are generated using Crystal Explorer 17.5 [39,40]. The HS analysis helps in identifying the intercontacts within crystals [41]. In the Hirshfeld surface plots, the white surface represents contacts with distances equal to the sum of van der Waals radii, and the blue and red colors indicate distances longer (distant contact) or shorter (in close contact) than the sum of the van der Waals radii. The dnorm plots were calculated over the ranges -0.0600 to 1.3750 a.u. (4), -0.1545 to 1.3476 a.u. (5), and -0.1308 to 1.6040 a.u. (6) respectively. The bright-red spots on the surface over atoms in the molecule indicate they are either donors or acceptors, with corresponding red and blue regions representing negative (hydrogen-bond acceptors) and positive (hydrogen-bond donors) potentials on the HS mapped over the electrostatic potential [42].

3.5. DFT Calculations

Optimization and frequency calculations were conducted using the DFT method at the B3LYP/6-3111++G(d,p) level of theory as implemented in Gaussian 16W [43]. The tight convergence in the SCF procedure (Default in Gaussian 16), which uses a combination of EDIIS and CDIIS algorithms is used as the criterion for molecular structure convergence. X-ray atom coordinates of 4, 5 and 6 were used as starting input coordinates in DFT calculations. Molecular and electronic properties of 4, 5 and 6 were calculated at the same level of theory. Frequency calculations confirm that all optimized structures are true minima, i.e., all imaginary frequencies are positive.

3.6. Molecular Docking Studies

A molecular docking study was performed to estimate the binding affinity of 4, 5 and 6 into the binding site of cholesterol oxidase using the AutoDock 4 package. The X-ray coordinates of cholesterol oxidase complexed with steroid substrate flavin-adenine dinucleotide were downloaded from the RCSB website using the 3COX pdb file [44] and its molecular docking is validated by re-docking the original ligand flavin-adenine dinucleotide into the binding site of cholesterol oxidase. The flavin-adenine dinucleotide fits well into the binding site of cholesterol oxidase and forms a stable complex with its amino acids of binding energy of -16.19 kcal/mol.Å. The superposition between the geometrical coordinates of the original ligand before and after it yields an RMSD value of 1.69Å. Molecular docking stepwise has been reported previously. Protein-ligand interactions were visualized using the Discovery Studio Visualizer Client v21.1.0.20298.

3.7. Molecular Dynamics Simulation

The highest binding affinity complex was selected for the MD simulation. For this purpose, the Desmond v2022-1 module within the Schrödinger package was used [21,45]. To prepare the protein-ligand complex, the Maestro v2022-1 software was used, and the protein preparation wizard was employed. This process involved removing water molecules that were more than 5Å away from the ligand, adding any missing hydrogen atoms, determining the bond order of the ligand, and identifying the optimal protonation states for both the ligand and the histidine residues within the protein. Additionally, the structure refinement was accomplished through a restrained minimization using the Prime module and the OPLS 2005 forcefield [46]. After preparing the protein-ligand complex, a cubic box with 10 Å dimensions was created using the SPC (simple point charge) water model to solvate the system. To neutralize the system, counter ions like sodium were added [47]. Additionally, the assembled system underwent an energy minimization process with a convergence threshold for the energy gradient set at 1 kcal/mol. Subsequently, pre-equilibration was conducted by Desmond's default six-step relaxation protocol. The NPT ensemble was employed for all runs, with the pressure set at 1 bar and the temperature set at 300 K. Temperature control used the Nosé–Hoover chain coupling scheme, while pressure control was managed using the Martyna–Tuckerman–Klein chain coupling scheme with a coupling constant of 2.0 ps [48]. The simulation was set to run for 80 ns, saving trajectories every 4.5 ps. non-bonded interactions were computed using the r-RESPA integrator, where long-range interactions were calculated every three steps, and short-range interactions were updated at each step [49]. Following the simulation, the interaction diagram module in Maestro was employed to analyze the trajectories. Plots depicting root mean square deviation (RMSD), solvent accessible surface area (SASA), radius of gyration (Rg), root mean square fluctuation (RMSF), and protein-ligand interactions were generated. Moreover, the analysis included evaluations of principal components (PCA), dynamic cross-correlation matrix (DCCA), and free energy landscape (FEL). These assessments were conducted using the trajectory file obtained from the MD simulation. The PCA and DCCA procedures employed the Bio3D (v2.4-4) package within the R (v4.2.2) statistical programming software [50]. Similarly, for FEL analysis, the Geo Measures (v0.9) plugin of PyMOL (v2.5) was employed [51].

3.8. MM-GBSA Analysis

"Molecular Mechanics-Generalized Born Surface Area Analysis" is a method used for studying the interactions between different molecules. It combines the principles of molecular mechanics and generalized Born models [19]. By offering insights into the energetics of molecular binding, this approach aids in the understanding of ligand-receptor interactions [52]. The binding free energy ($\mathsf{\Delta }{\mathrm{G}}_{\mathrm{b}\mathrm{i}\mathrm{n}\mathrm{d}}$) was determined by applying the specified equation to the last 300 frames extracted from the MD simulation:

$$\mathsf{\Delta }{\mathrm{G}}_{\mathrm{b}\mathrm{i}\mathrm{n}\mathrm{d}}=\mathsf{\Delta }{\mathrm{G}}_{\mathrm{s}\mathrm{o}\mathrm{l}\mathrm{v}}+\mathsf{\Delta }{\mathrm{G}}_{\mathrm{M}\mathrm{M}}+\mathsf{\Delta }{\mathrm{G}}_{\mathrm{S}\mathrm{A}}$$

(1)

where ΔG_solv is the solvation energy of the ligand-receptor complex, ΔGMM indicates free energy differences of ligand + protein complex, and the total energies of protein and ligand in isolated form. Additionally, ΔGSA denotes the area energy differences of the surface between protein and ligand. All MM-GBSA calculations were executed by the Prime-MMGBSA module in the Schrödinger package [43].

4. Conclusions

This article discusses the synthesis of new derivatives of N-alkylated phenytoin (4-13). Their structures were determined by ¹H and ¹³C NMR spectral data, and X-ray diffraction studies confirmed the structures of 4, 5 and 6. To accurately replicate the experimental spectral data, the energy gap between the HOMO and LUMO orbitals is a key indicator of the structural stability of molecules. The results obtained show a clear trend, with the following stability order: 6 > 5 > 4. Hirshfeld surfaces were used to confirm the presence of intermolecular interactions. In addition, molecular docking investigations of 4, 5 and 6 demonstrated the expected binding pose. The study's findings can guide drug design or experimental efforts. Molecular dynamics simulation studies showed that protein-4 complex is highly stable and included an RMSD plot of protein Cα atoms and ligand fit the protein, an RMSF plot where green lines indicate ligand contact and a distribution of protein secondary structure elements by residue index.