Abstract

A new nicotinamide derivative, (E)-N-(4-(1-(2-(4-benzamidobenzoyl)hydrazone)ethyl)phenyl)nicotinamide, was designed as a VEGFR-2 inhibitor. Utilizing the density functional theory (DFT) calculations, the three-dimensional structure of the designed compound was determined, shedding light on its stability and reactivity. Molecular docking revealed its capability to inhibit VEGFR-2, which was further supported by molecular dynamics (MD) simulations confirming its binding to the target protein. In addition, molecular mechanics-generalized born surface area (MM-GBSA), protein-ligand interactions profiler (PLIP), and essential dynamics studies provided further validation of the compound’s precise binding with optimal energy. Then, the “compound 10” was synthesized and subjected to in vitro assays. Compound 10 inhibited VEGFR-2 with an IC50 value of 105.4 ± 0.896 nM, comparing sorafenib’s IC50 value of 61.65 ± 0.934 nM. Besides, it exhibited cytotoxicity against HepG2 and MCF-7 cancer cell lines, with IC50 values of 35.78 ± 0.863 μM and 57.62 μM ± 0.871, comparing sorafenib’s IC50 values of 5.95 ± 0.917 μM and 8.45 ± 0.912 μM. Furthermore, compound 10 demonstrated a lower level of toxicity towards Vero cell lines, with an IC50 value of 127.3 μM. Likewise, compound 10 induced apoptosis in HepG2 cell lines through a flow cytometric analysis in addition to an increase in the levels of caspase-3 and caspase-9. Moreover, compound 10 hindered the migration and healing abilities of HepG2 cells. In conclusion, our study positions compound 10 as a promising candidate for further chemical modifications and biological evaluations.

1. Introduction

Cancer, a serious health condition and the second leading cause of death globally after cardiovascular diseases, continues to pose significant challenges [1]. The incidence and mortality rates of cancer have been steadily increasing in recent decades [2]. Despite efforts to develop safer and more effective antitumor drugs, traditional chemotherapy treatments are still associated with systemic toxicity, resulting in severe side effects [3]. Consequently, researchers are actively exploring the development of novel anticancer drugs that offer improved selectivity, efficacy, and safety. Cancer is characterized by uncontrolled cell growth, angiogenesis (formation of new blood vessels), and metastasis (spread to other parts of the body) [4]. VEGFR-2, a transmembrane tyrosine kinase receptor, plays a crucial role in cancer progression by regulating various cellular processes, including angiogenesis and metastasis [5, 6]. Inhibiting VEGFR-2 signaling has emerged as a promising approach for cancer treatment. By impeding the formation of new blood vessels and disrupting angiogenesis, VEGFR-2 inhibitors can limit the supply of nutrients and oxygen to tumors, resulting in reduced tumor growth and metastasis [7]. Numerous types of cancer, such as breast, prostate, colon, cervix, brain, and ovarian cancer, exhibit high levels of VEGFR-2 expression, making it an attractive target for therapeutic interventions [8]. Consequently, extensive research is underway to develop VEGFR-2 inhibitors with enhanced efficacy, specificity, and safety, which could significantly improve cancer treatment outcomes for patients.

Our laboratory’s research has presented a variety of anticancer compounds targeting VEGFR-2. These compounds include quinoline [9], naphthalene [10], nicotinamide [11], theobromines [12, 13], thiazolidine [14], pyridine [15], indole [16], and isatin [17] derivatives.

1.1. Rational

Considering sorafenib I (Figure 1) as a representative example of FDA-approved VEGFR-2 inhibitors, it becomes evident that it possesses four essential pharmacophoric features necessary for optimal binding to the active site of VEGFR-2. These features include the following: (i) a heteroaromatic ring system that aligns with the hinge region [18], (ii) a central aromatic linker occupying the gatekeeper region [19], (iii) a pharmacophore (hydrogen-bonding center) facilitating hydrophilic interactions with Glu883 and Asp1044 at the DFG-motif region [20], and (iv) a hydrophobic tail occupying the allosteric pocket [2123].

Recently, our research team designed and synthesized a nicotinamide derivative, compound II, as a potential anticancer agent targeting VEGFR-2 (Figure 2). This compound exhibited promising inhibitory activity against VEGFR-2 with an IC50 value of 51 nM. To further investigate its potential, compound II underwent comprehensive in silico and in vitro evaluations [24]. The encouraging outcomes prompted us to modify compound II to enhance its effectiveness as a VEGFR-2 inhibitor. In this study, we employed a substitution pattern variation strategy on the central phenyl ring, transforming it from a metadisubstituted structure to a paradisubstituted structure. This modification aimed to induce a change in the orientation of the new compound, potentially facilitating an ideal binding pattern at the active site.

2. Results and Discussion

2.1. Molecular Docking

Using the MOE software, molecular docking simulations were conducted to investigate the potential binding pattern of the targeted nicotinamide derivative within the active site of VEGFR-2. The docking simulation utilized a target enzyme with the PDB ID of 2OH4. This protein was specifically selected for its clear binding mode with sorafenib, as well as its accurate manifestation of the pharmacophoric features. The binding mode of sorafenib to the VEGFR-2 active site cavity was consistent with previously observed interactions [25, 26]. Sorafenib occupied key structural components of the receptor, including the hinge region (forming a hydrogen bond with Cys917) and the DFG-motif region (establishing hydrogen bonds with Glu883 and Asp1044). In addition, it formed two networks of hydrophobic interactions in the linker region and the terminal hydrophobic regions (Figure 3).

The binding potential of the targeted nicotinamide derivative within the ATP binding pocket of VEGFR-2 was investigated. It was observed that the pyridine moiety of the compound effectively occupied the hinge region, forming a hydrogen bond with the critical amino acid Cys917. Furthermore, the central phenyl moiety (PhI) interacted with the linker region, establishing hydrophobic interactions with Ala864, Val914, and Val846. In addition, the amide I group of the targeted nicotinamide derivative formed two hydrogen bonds in the DFG-motif region, specifically with Glu883 and Asp1044. The compound fitted well into the hydrophobic region and established interactions with its hydrophobic anchors (Ile886 and Leu887) through the second phenyl moiety (PhII) (Figure 4). These findings provide a potential inhibitory activity of the designed candidate against VEGFR-2.

2.2. MD Simulations

The analysis of the production run (100 ns) reveals that the system exhibits a consistent and stable binding conformation. Both the apo VEGFR-2 protein (represented by the blue curve) and the holo VEGFR-2 protein (represented by the red curve) display average root mean square deviation (RMSD) values of around 2.7 Å when considering all alpha carbon atoms (Figure 5(a)). In contrast, the RMSD of the ligand demonstrates minimal variation, with an average value of 1.1 Å (Figure 5(b)). Comparable patterns are observed in the radius of gyration (RoG) (Figure 5(c)) and solvent-accessible surface area (SASA) (Figure 5(d)), where the averages for both systems remain consistently at 20.5 Å and 17500 Å2, respectively. The hydrogen bond analysis (Figure 5(e)) shows a consistent variation with an average of 70 bonds for both systems. In the RMSF (root mean square fluctuation) plot (Figure 5(f)), except for specific regions such as the free N-terminal (11.7 Å in the holoprotein), Tyr994: Asp996 loop in the apo system (2.5 Å), Gly1046: Leu1065 loop in the apo system (6.6 Å), Lys1053: Asp1062, and the C-terminal (10.3 Å in the apoprotein and 8.5 Å in the holoprotein), the fluctuations of amino acids are minimal, with values less than 2 Å.

Throughout the simulation, the ligand maintains an average distance of 7.2 Å from the protein’s center of mass. Overall, the binding of the ligand to the protein remains stable, as demonstrated in Figure 5(g).

2.3. MM-GBSA Analysis

MM-GBSA is a computational approach employed in estimating the binding free energy of biomolecular complexes, such as protein-ligand interactions. It integrates molecular mechanics (MM) to calculate internal energy and generalized born surface area (GBSA) for accounting solvation effects. This technique holds great significance in drug discovery as it aids in predicting binding affinity and choosing promising drug candidates for subsequent evaluation [27]. The MM-GBSA binding free energy study between the targeted nicotinamide derivative and VEGFR-2 provides insights into the different components contributing to the binding interaction (Figure 6). The typical van der Waals and electrostatic interactions for the nicotinamide derivative are estimated to be approximately −58.85 kcal/mol and −27.94 kcal/mol, respectively, resulting in an average total binding energy of −43.73 kcal/mol. The decomposition analysis depicted in Figure 7 identifies the specific amino acids that are in close proximity to the ligand and have the most significant impact on the interaction. A value larger (or lesser) than −1 kcal/mol is attributed to the following amino acids: Ile886 (−1.02 kcal/mol), Val897 (−1.43 kcal/mol), Phe916 (−1.69 kcal/mol), Cys917 (−2.02 kcal/mol), Cys1022 (−2.15 kcal/mol), Ile1023 (−2.44 kcal/mol), Leu1033 (−1.36 kcal/mol), Cys1043 (−4.61 kcal/mol), and Asp1044 (−1.12 kcal/mol).

2.4. PLIP Analysis

PLIP, protein-ligand interaction profiler, is a computational tool utilized in molecular biology to scrutinize the interactions between proteins and ligands. Its purpose is to identify and elucidate various types of interactions occurring within protein-ligand complexes, including hydrogen bonds, hydrophobic interactions, and π-π stacking. By offering a comprehensive profile of these interactions, PLIP plays a vital role in comprehending the binding mode and strength between a protein and its ligand. This insight proves invaluable in drug discovery and molecular biology, enabling researchers to tailor compounds for specific biological targets. Moreover, PLIP can be employed to investigate protein-protein interactions, furnishing valuable information on cellular signaling pathways and protein functions [28]. Subsequently, the trajectory underwent clustering to select representative frames for each cluster. As outlined in the methods section, the elbow method was employed to automatically determine the optimal number of clusters, resulting in four clusters. The PLIP website was utilized to evaluate the quantity and nature of interactions between the targeted nicotinamide derivative and VEGFR-2 for each cluster representative (Table 1). Among all the cluster representatives, hydrophobic interactions were the most prevalent, accounting for 34 interactions, compared to 16 hydrogen bonds. This observation aligns with the discrepancy observed in the electrostatic and van der Waals energy levels obtained from MM-GBSA. Cys917, Ile1023, and Asp1044 were the most common amino acids involved in hydrogen bonding across all cluster representatives, while Leu838, Val897, Ile1023, Leu1033, and Asp1044 were the most frequently occurring amino acids involved in hydrophobic contacts across all four cluster representatives. In addition to providing information on interaction types and frequencies, PLIP also generated a .pse file to visualize the three-dimensional conformation of the ligand and its interaction with the protein (Figure 8).

2.5. Essential Dynamics Analysis

Essential dynamics is a computational method utilized in molecular dynamics simulations to pinpoint the primary and crucial movements within a biomolecular system, such as a protein or a protein-ligand complex. This technique simplifies the analysis of large-scale molecular systems by focusing on the most significant motions. These motions encompass collective movements or structural changes that play a vital role in the biological function of the system. By isolating these key motions, valuable insights into the dynamic behavior of biomolecules can be obtained. This knowledge is particularly pertinent in comprehending the functional mechanisms of proteins and their interactions with ligands, which holds relevance in areas such as drug development and protein engineering [29, 30]. The principal component analysis was employed to determine the underlying cause of the trajectory’s highly synchronized motion. The selection of the essential subspace was based on criteria such as the scree plot, eigenvector distribution, and variance analysis, as outlined in the methods section. Examination of the scree plot revealed a significant flattening of the slope after the second principal component. Among the eigenvectors, the top three accounted for nearly 78.5% of the total variance, with the first eigenvector alone capturing 69.2% of the variance (Figure S.1). Notably, the first two eigenvectors deviated from a Gaussian distribution (Figure S.2). Consequently, the reduced subspace was represented by the first three eigenvectors.

To assess the randomness of motion represented by the first 10 eigenvectors, the cosine content was calculated for both apo and holo VEGFR-2 simulations. The cosine content values for the first 10 eigenvectors of both the apoprotein and holoprotein were found to be less than 0.3 (Figure S.3). Upon comparing the two trajectories, it became evident that they were sampled differently, as indicated by the limited overlap observed in the first three eigenvectors (1.8% according to the root mean square inner product, RMSIP). Furthermore, the RMSIP analysis revealed a mere 13.2% similarity between the apo and holo C matrices.

2.6. Bidimensional Projection Analysis

The outcomes obtained by projecting each trajectory onto the first three eigenvectors of the new C matrix are depicted in Figure 9. Figure 9(a) illustrates the projection on the first two eigenvectors, demonstrating distinct average structures for the two trajectories and a moderate level of overlap in their sampling. In Figure 9(b), the average structures projected onto the first and third eigenvectors exhibit a significant degree of similarity, and their trajectories also show substantial overlap. Figure 9(c) reveals a similar level of overlap between the two trajectories as seen in Figure 9(a) when projecting onto the second and third eigenvectors.

In Figure 10, we utilized Porcupine diagrams to provide a visual representation of the movement associated with the first three eigenvectors. Among these, the most prominent motion is observed in the Gly1046: Leu1065 loop. In more detail, the first eigenvector of each trajectory delineates a specific motion pattern. The apoprotein, represented by the green cartoon, displays a motion characterized by a loop opening. Conversely, the holoprotein, depicted in red, exhibits a motion pattern indicative of loop closure. Similarly, the second eigenvector highlights a comparable motion trend. Moving on to the third eigenvector, both the apoprotein and holoprotein showcase a motion primarily characterized by loop opening. These observations provide crucial insights into the dynamic behavior of the system, particularly with regard to the Gly1046: Leu1065 loop region.

2.7. DFT Calculation
2.7.1. Geometry Optimization and Mulliken Charge

The geometry optimization of the targeted nicotinamide derivative was performed using Gaussian 09 software at the B3LYP/6-31G++(d, p) level of theory. The resulting optimized structure, depicted in Figure 11(a), provides information on the imine bond length (C16-N18) and the bond angles at the ends of the created bond. The chemical system comprises 59 atoms and 250 electrons, with a ground energy or total energy (TE) of −43013.1 eV, indicating a stable structure. The predicted dipole moment magnitude of 5.451 Debye suggests a highly reactive and polarizable nature of the compound. Figure 11(b) illustrates the polarizability and distribution of Mulliken charges across the chemical structure, offering insights into the potential for charge transfer within the targeted nicotinamide derivative. The blue and orange dashed circles in Figure 11(b) highlight the most electronegative and electropositive atoms, respectively.

2.7.2. Frontier Molecular Orbital (FMO) Analysis

The HOMO (highest occupied molecular orbital) and LUMO (lowest unoccupied molecular orbital), known as the frontier molecular orbitals, provide quantitative information about the affinity of the chemical system to interact with the target. In this study, the HOMO and LUMO were computed for the targeted nicotinamide derivative. Figure 12 illustrates that the LUMO density lobes are distributed throughout the structure, except for the N-phenylbenzamide terminal, while the localized area is reduced in the case of the HOMO density. The HOMO energy (EHOMO), LUMO energy (ELUMO), and the HOMO-LUMO energy gap (Egap) are displayed. The Egap value is a significant stability parameter that describes the reactivity and stability of a molecule. A structure with a short Egap indicates a more reactive and polar structure, where electrons can be readily transferred to an acceptor. The calculated theoretical Egap suggests a significant inhibitory reactivity due to the molecule’s charge transfer capacity [31]. Table 2 provides the ionization potential (IP) and electron affinity (EA) values, with a short Egap indicating a high level of interaction with the target, resulting from electron donation by the developed anticancer drug. Table 2 shows the ionization potential (IP) and electron affinity (EA) values, with a short Egap indicating high interaction with the target due to electron donation from the targeted nicotinamide derivative.

2.7.3. Chemical Reactivity Descriptors and Total Density of State (TDOS)

Koopman’s theorem was applied to compute the electronic properties of the targeted nicotinamide derivative, utilizing various global reactivity parameters such as ionization potential (IP), maximal charge acceptance (Nmax), electron affinity (EA), chemical potential (μ), energy change due to charge transfer (∆E), chemical hardness (η), chemical electronegativity (χ), electrophilicity (ω), and global softness (δ). The results in Table 2 indicate high values for parameters such as σ (electrophilicity), ω (chemical potential), and η (chemical hardness), suggesting a significant bioreactivity. These parameters have an impact on the compound’s ability to acquire additional charge (∆Nmax) from neighboring chemical systems or targets [32]. The high inhibitory reactivity of the synthesized medication can be attributed to the favorable acquisition of additional charge and the low value of Egap. When a compound is in an equilibrium state, the number of occupied states per unit volume can be estimated by multiplying the probability function by the density. The findings obtained from this analysis can be utilized to examine various properties of the molecule. Figure 13 displays the spectrum of the total density of all states, with the highest density observed in the empty orbitals according to the TDOS (total density of state) spectrum.

2.7.4. Electrostatic Potential (ESP) Maps

Electrostatic attraction is a crucial factor influencing the binding of molecules to proteins. Quantitative analysis of these attractive forces allows for the calculation of the energy of the protein-drug complex. Steric and electrostatic forces are the primary driving factors in biomolecules. The electrostatic potential (ESP) surface of a molecule exhibits various colored patches, which correspond to different electron density levels. Atoms with partial negative charges (typically highly electronegative atoms) appear as red patches and can act as hydrogen-bond acceptors. Conversely, atoms with partial positive charges and low electronegativity are depicted as blue patches and can serve as hydrogen-bond donors. Neutral atoms, which have no charge, range in color from green to yellow and are involved in π-stacking interactions and other types of interactions. These distinct colored regions facilitate the prediction of how molecules will interact with the target and aid in understanding reaction mechanisms.

In the structure of the targeted nicotinamide derivative, blue regions indicate the presence of hydrogen-bond donors at hydrogen atoms, while red patches suggest the potential for hydrogen-bond formation with polar amino acids at oxygen atoms. The clouds surrounding the phenyl rings, which mostly appear greenish in color, facilitate the formation of π-stacking bonds with aromatic amino acids (Figure 14).

2.8. In Silico ADMET Analysis

The evaluation of both pharmacokinetic characteristics and biological activity is crucial during the approval process of a new medicine. To prevent delays or potential drug withdrawal, it is necessary to assess the pharmacokinetic features of a new compound early in the drug discovery phase [33]. In this study, Discovery Studio 4.0 software was employed to compute ADMET (absorption, distribution, metabolism, excretion, and toxicity) parameters for the targeted nicotinamide derivative in comparison to sorafenib. Interestingly, the ADMET results of the targeted nicotinamide derivative exhibited a high degree of similarity to those of sorafenib (Figure 15 and Table 3). Both compounds demonstrated a low ability to cross the blood-brain barrier (BBB), exhibited favorable levels of intestinal absorption (LIA), did not inhibit CYP2D6, and displayed plasma protein binding (LPPB) abilities exceeding 90%. Notably, the targeted nicotinamide derivative showed good aqueous solubility (LAS), which was found to be superior to the predicted low solubility of sorafenib.

2.9. Toxicity Studies

The early assessment of toxicity is crucial in order to minimize drug approval failures [34]. However, traditional in vitro and in vivo research methods for toxicity evaluation are limited by ethical restrictions, budget constraints, and time consumption. Consequently, in silico techniques have gained significant advantages and are being utilized for toxicity prediction [35]. In this particular study, the toxicity of the targeted nicotinamide derivative was estimated using nine toxicity models available in the Discovery Studio program, and the results were compared to those of sorafenib. The chemical structure of the targeted nicotinamide derivative underwent a comprehensive analysis by comparing it with a large dataset of previously studied drugs in order to predict its toxicity and safety profile. This analysis involved the calculation of various essential parameters, such as probability, enrichment, Bayesian score, Mahalanobis distance, and Mahalanobis distance value, for each toxicity model. For a detailed and comprehensive understanding of the results, refer to the detailed toxicity report in the Supplementary data (available here).

Fortunately, the results revealed that the targeted nicotinamide derivative exhibited safe values across all the models employed (Table 4).

2.10. Chemistry

According to our design and the conducted studies, the targeted nicotinamide derivative exhibits a strong ability to bind to and inhibit VEGFR-2. Computational investigations have provided further support for the proposed binding affinity. Following these findings, the synthesis of the targeted nicotinamide derivative was carried out as outlined in Scheme 1.

The synthesis process involved several steps. First, nicotinoyl chloride 3 was obtained by acylating nicotinic acid 2 with thionyl chloride. Subsequently, the appropriate nicotinamide derivative 4 was synthesized by reacting nicotinoyl chloride 3 with 4-aminoacetophenone. On the other hand, 4-aminobenzoic acid 5 was esterified by refluxing it in a solution of methanol and sulfuric acid, resulting in the formation of the corresponding ester 6. In accordance with the described procedures, ester 6 underwent a benzoylation reaction with the addition of benzoyl chloride 7 drop-by-drop in a mixture of dichloromethane (DCM) and triethylamine (TEA) at 0°C. This reaction produced the desired benzoyl derivative 8. Furthermore, the acid hydrazide derivative 9 was obtained by reacting compound 8 with hydrazine hydrate in pure ethanol. Finally, compound 10, the final nicotinamide derivative, was synthesized by condensing compound 9 with compound 4 in an ethanol/glacial acetic acid mixture.

The determination of the molecular formula (C28H23N5O3) was supported by EI-MS and elemental analysis (Supplementary data (available here)). The IR spectrum of compound 10 exhibited absorption bands at 1670 cm−1 for C=O and 3203 cm−1 for NH, confirming the presence of these functional groups. The 1H NMR spectrum showed significant singlet signals for the three amidic protons at δ 10.70, 10.63, and 10.54 ppm (Supplementary data (available here)). In addition, a characteristic singlet signal corresponding to CH3 protons appeared at δ 2.40 ppm. The 13C NMR spectrum further confirmed the chemical structure of the compound (Supplementary data (available here)).

2.11. Biology
2.11.1. In Vitro VEGFR-2 Inhibition

In the pursuit of developing an inhibitor targeting VEGFR-2, a key regulator of angiogenesis, compound 10 was specifically designed. Through in silico investigations, compound 10 exhibited promising inhibitory properties. Building upon these findings, we further explored the inhibitory potential of compound 10 against the VEGFR-2 protein through in vitro studies. Remarkably, as shown in Figure 16, compound 10 demonstrated promising inhibition of the VEGFR-2 protein, as indicated by an IC50 value of 105.4 ± 0.896 nM, in comparison to sorafenib’s IC50 value of 61.65 ± 0.934 nM. These results intriguingly aligned with the promising computational findings, providing strong evidence of compound 10’s robust suppressive capabilities.

2.11.2. Cytotoxicity and Safety

Compound 10 has emerged as a promising agent targeting VEGFR-2, with its potential demonstrated through comprehensive in silico and in vitro investigations. This study aimed to assess the cytotoxic effects of compound 10 against HepG2 and MCF-7 cancer cell lines, using sorafenib as a reference drug. The results (Figure 17) revealed the cytotoxic activity of compound 10 against both cell lines, with IC50 values of 35.78 ± 0.863 μM and 57.62 μM ± 0.871, respectively, comparing sorafenib’s IC50 values (5.95 ± 0.917 μM and 8.45 ± 0.912 μM) against the same 2 cell lines, respectively. Moreover, compound 10 displayed a higher IC50 value of 127.3 μM against Vero cell lines, resulting in selectivity index values of 3.5 and 2.2 against HepG2 and MCF-7 cell lines, respectively.

2.11.3. Apoptosis Assay

Flow cytometry analysis utilizing Annexin V and PI double stains was employed to examine the apoptotic effects of compound 10, at a concentration of 35.78 µM, on HepG2 cells, which were found to be the most sensitive cell line. The results revealed a significant increase in the percentage of HepG2 cells undergoing early apoptosis (from 0.71% to 3.63%) and late apoptosis (from 0.13% to 11.84%) following treatment with compound 10 for 48 hours, as compared to the control group. Furthermore, the percentage of necrotic cells also increased to 4.24% in compound 10-treated cells, compared to 2.22% in control cells (Table 5 and Figure 18). These findings indicate that compound 10 is capable of inducing apoptosis and disrupting the cell cycle in HepG2 cells warranting further investigation to understand and enhance its activity.

2.11.4. Apoptotic Proteins Assay

The aim of this study was to assess the impact of compound 10 on the expression of apoptosis-related proteins, specifically caspase-3 and caspase-9, in HepG2 cells. The experimental approach involved treating the cells with compound 10 and comparing the expression levels of these proteins with a control group of untreated cells using quantitative real-time polymerase chain reaction (qRT-PCR). The results revealed significant modulations of caspase-3 and caspase-9 expressions in response to compound 10 treatment, indicating its potential as an anticancer agent with apoptotic properties. Compound 10 was observed to upregulate the levels of caspase-3 and caspase-9, both of which play crucial roles in the apoptotic pathway. Caspase-3 acts as an initiator caspase, initiating the apoptotic cascade [36], while caspase-9 functions as an executioner caspase downstreaming in the apoptotic signaling pathway [37]. The substantial increase of 2.61-fold in caspase-3 levels and 3.66-fold in caspase-9 levels in compound 10-treated HepG2 cells compared to the control group (Table 6) further supports its potential as a promoter of apoptosis.

2.11.5. The Effect of Compound 10 on HepG2’s Migration and Healing

The assessment of cancer cell migration and healing potential in vitro can be achieved by using the wound healing (scratch or migration) assay, a simple and cost-effective technique [38]. This method involves creating a scratch in a monolayer of cancer cells, measuring the initial width of the scratch, and monitoring the closure of the scratch over time in both treated and untreated cells. In this study, images of the scratch area were captured for treated and untreated HepG2 cell lines at 0 and 48 hours. The results of the scratch assay (Table 7 and Figure 19) demonstrated that the untreated HepG2 cells exhibited a significant closure of the scratch by 65.9% after 48 hours. However, treatment with compound 10 only resulted in a 2% reduction in the width of the scratch, indicating a substantial inhibitory effect on the migration and healing ability of the cancer cells. These findings suggest that compound 10 holds promise as a therapeutic agent for inhibiting cancer cell migration and warrants further investigation into its potential as well as future modifications.

3. Experimental Work

3.1. Docking Studies

Docking was conducted for compound 10 against VEGFR-2 by MOE2014 software [39]. Additional detailed information can be found in the Supplementary section (available here), providing further insights and clarification on the experimental work presented in this study.

3.2. MD Simulations Studies

To assess the stability of the VEGFR-2_compound 10 complex and explore the interactions and distinctions between the holo and apo structures, a 100-ns classical molecular dynamics (MD) simulation was conducted using GROMACS 2021 [40]. The input data were obtained by the CHARMM-GUI web server solution builder module [41, 42]. Additional detailed information can be found in the Supplementary section (available here), providing further insights and clarification on the experimental work presented in this study.

3.3. Binding Free Energy Calculation Using MM-GBSA

The molecular mechanics-generalized born surface area (MM-GBSA) approach was employed to assess the binding strength of the VEGFR-2_compound 10 complex, utilizing the gmx_MMPBSA program [43, 44]. Additional detailed information can be found in the supplementary section, providing further insights and clarification on the experimental work presented in this study.

3.4. ED Analysis

Mass-weighted covariance matrix (C) of alpha carbons in amino acids (Glu826: Leu1161) correlated mobility along MD trajectories by employing PCA analysis [45].

3.5. Bidimensional Projection Analysis

In order to directly compare the frames in the reduced subspace, a specific procedure was implemented. First, the trajectories of the apoprotein and complex were combined. Then, they were aligned with the apoprotein configuration obtained during the equilibration phase. Subsequently, a new C matrix was generated for the merged trajectories. Finally, each trajectory was projected onto the new C matrix. The similarity between the two trajectories was assessed by plotting the projection on the first three eigenvectors, utilizing different combinations of eigenvectors [46].

3.6. DFT

DFT studies were performed for compound 10 by Gaussian 09 and GaussSum3.0 programs [31]. Additional detailed information can be found in the Supplementary section, providing further insights and clarification on the experimental work presented in this study.

3.7. ADMET Studies

ADMET studies were performed for compound 10 by Discovery Studio 4.0 [47]. Additional detailed information can be found in the Supplementary section (available here), providing further insights and clarification on the experimental work presented in this study.

3.8. Toxicity Studies

Toxicity studies were conducted for compound 10 by Discovery Studio 4.0 [47]. Additional detailed information can be found in the Supplementary section (available here), providing further insights and clarification on the experimental work presented in this study.

3.9. Synthesis of Compound 10

N-(4-(Hydrazinecarbonyl)phenyl)benzamide 9 and N-(4-acetylphenyl)nicotinamide 4 were mixed and thoroughly dissolved in a round bottom flask containing absolute ethanol (25 mL). Later, the entire mixture was refluxed for 6 h while being catalyzed by drops of glacial acetic acid. The mixture was concentrated and cooled following the reaction. Crystallization from methanol was used to filter and purify the collected product.

3.9.1. (E)-N-(4-(1-(2-(4-Benzamidobenzoyl)hydrazone)ethyl)phenyl)nicotinamide

Yellowish white crystal (yield, 70%); m. . = 210–212°C; IR (KBr) ν cm−1: 3203 (NH), 3064, 3026 (CH aromatic), 1670 (C=O); 1H NMR (400 MHz, DMSO-d6) δ 10.70 (s, 1H), 10.63 (s, 1H), 10.54 (s, 1H), 9.15 (s, 1H), 8.79 (d, J = 4.8 Hz, 1H), 8.33 (d, J = 7.9 Hz, 1H), 8.01 (d, J = 7.5 Hz, 2H), 7.95 (s, 4H), 7.89 (s, 4H), 7.60 (dt, J = 15.1, 7.6 Hz, 4H), 2.40 (s, 3H); 13C NMR (101 MHz, DMSO-d6) δ 166.34, 164.67, 163.69, 155.22,152.69, 149.21, 142.65, 140.48, 136.02, 135.14, 134.01, 132.33, 130.97, 129.24, 129.17, 128.95, 128.26, 127.44, 124.01, 120.29, 119.90, 14.83. Mass (m/z): 477 (M+, 59%); Anal. Calcd. For C28H23N5O3 (477.18): C, 70.43; H, 4.86; N, 14.67; Found: C, 70.29; H, 5.02; N, 14.89%.

3.10. In Vitro VEGFR Inhibition

In Vitro VEGFR inhibition analysis was conducted for compound 10 with a human VEGFR ELISA kit [48]. Additional detailed information can be found in the Supplementary section, providing further insights and clarification on the experimental work presented in this study.

3.11. In Vitro Antiproliferative Activity

The analysis of antiproliferative activity was conducted for compound 10 against two cancer cell lines (HepG2 and MCF) by the MTT procedure [49, 50]. Additional detailed information can be found in the Supplementary section, providing further insights and clarification on the experimental work presented in this study.

3.12. Safety Assay

Safety assay analysis was conducted for compound 10 by MTT procedure utilizing Vero cell lines. Additional detailed information can be found in the Supplementary section, providing further insights and clarification on the experimental work presented in this study.

3.13. Apoptosis Analysis

Apoptosis analysis was conducted for compound 10 by flow cytometry analysis technique [51]. Additional detailed information can be found in the Supplementary section, providing further insights and clarification on the experimental work presented in this study.

3.14. Apoptotic Proteins Assay

The analysis of apoptotic proteins was conducted using qRT-PCR [52] employing the Qiagen RNA extraction kit and Bio-Rad SYBR Green PCR MMX kit. Additional detailed information can be found in the Supplementary section (available here), providing further insights and clarification on the experimental work presented in this study.

3.15. Wound Healing Assay

The potentialities of HepG2 cancer cells to migrate as well as to heal were determined after and before treatment of compound 10 as described before [53]. Additional detailed information can be found in the Supplementary section (available here), providing further insights and clarification on the experimental work presented in this study.

4. Conclusion

In summary, this study showcases the successful application of computer-assisted drug design in identifying a novel nicotinamide derivative, compound 10, as a potential apoptotic and cytotoxic agent. Through various computational techniques such as DFT calculations, molecular docking, and molecular dynamics simulations, compound 10’s stability, reactivity, and binding potential with VEGFR-2 were assessed. The results from MM-GBSA, PLIP, and essential dynamics studies further supported the precise binding and favorable energy of compound 10. In vitro assays confirmed its effectiveness in inhibiting VEGFR-2 and exerting cytotoxic effects on HepG2 and MCF-7 cancer cell lines, while exhibiting a lower toxicity towards normal cells. In addition, compound 10 demonstrated the induction of apoptosis in HepG2 cells and reduced their migration and healing abilities. Collectively, these findings suggest that compound 10 holds promise as a lead candidate for further modifications and research seeking the development of an apoptotic anticancer lead compound.

Data Availability

The data used to support the findings of the study are available from the corresponding author upon request.

Disclosure

An earlier version of this manuscript has been presented as a preprint (Eissa_ et.al 2023) [54].

Conflicts of Interest

The authors verify that they have no conflicts of interest associated with this publication involving any party or among the co-authors.

Acknowledgments

This research was funded by the Princess Nourah bint Abdulrahman University Researchers Supporting Project number PNURSP2024R142, Princess Nourah bint Abdulrahman University, Riyadh, Saudi Arabia. The authors would like to express their sincere gratitude to AlMaarefa University, Riyadh, Saudi Arabia, for supporting this research.

Supplementary Materials

Figure S.1: change in the eigenvalues with increasing the eigenvectors (blue). The cumulative variance retained in the eigenvectors (red). Figure S.2: 1st ten eigenvectors’ distribution. Figure S.3, the cosine content of the 1st ten eigenvectors for the two trajectories. Also, additional detailed information providing further insights and clarification on the experimental parts of molecular docking, MD simulations, MM-GBSA, PLIP, essential dynamics, ADMET, DFT, synthesis, and in vitro studies are presented in this study. Also, the Supplementary materials contain 1H, 13C, EI-MS, and IR spectra of compound 10. Finally, a detailed report for the conducted in silico toxicity studies has been provided in the Supplementary material. (Supplementary Materials)