2022 Volume 13 Issue 2

Structure-based Multi-targeted Molecular Docking and Molecular Dynamic Simulation Analysis to Identify Potential Inhibitors against Ovarian Cancer


Bandar Hamad Aloufi*


Ovarian Cancer (OC) is among the most prevalent cancers in females. OC is one of the deadliest and worst prognosis diseases. Currently, there are no approved OC screening tests or early detection methods. Hence, new screening, prevention, and early detection strategies are still highly demanded. Plant compounds have recently become increasingly important in developing new, effective, and affordable anti-cancer drugs. To investigate the protein-ligand interaction, molecular docking was used with the Molecular Operating Environment (MOE) tool to find the best inhibitor for the target proteins. Compound spatial affinity for the active sites of the NY-ESO-1, RUNX3, and UBE2Q1 proteins was calculated using molecular docking. ADMET analysis was used to determine the drug-likeness of the selected compounds, while MD simulation and MMGBSA/MMPBSA experiments were used to further understand the binding behaviors. Pre-clinical tests can help confirm the validity of our in silico studies and determine whether the compound can be used as an anti-cancer drug to treat OC.

Keywords: Ovarian cancer, Anti-cancer drugs, NY-ESO-1, RUNX3, UBE2Q1


OC is the deadliest gynaecological cancer and the fifth most common cause of cancer death among females. During her lifetime, 1 out of every 54 women will develop OC (Ghilardi et al., 2022). OC is classified into three types based on the origin of the cells: epithelial, stromal, and germ cell. Despite advances in therapeutic and diagnostic procedures, OC has the lowest survival rate of any gynaecological malignancy in developed nations (Siegel et al., 2016; Torre et al., 2018; Manchanda, 2022). Furthermore, it is critical to investigate the role of the tumor-causing microenvironment during the proliferation, early stages, and metastasis. Hence, it is critical to comprehend the root cause from various perspectives, such as molecular pathogenesis, histological subtypes, hereditary factors, epidemiology, treatment methods, and diagnostic perspectives (Auner et al., 2010; Ledermann et al., 2013; Kbirou et al., 2022). 

OC has a low survival rate and poor prognosis due to a lack of early screening methods and ineffective treatments for advanced stages of the disease. The high mortality rate from OC is because 30% of advanced-stage tumors do not respond to standard chemotherapy, and most responders relapse over time (Auner et al., 2010; Haque et al., 2022). Patients with OC typically experience recurrence and progression of the disease due to resistance to current chemotherapy treatment. Most chemotherapy drugs prescribed are synthetically derived and toxic to cancer cells and normal cells (Gupta et al., 2001; Pahwa et al., 2022). On the other hand, recent findings show that naturally extracted phytochemicals from plants have significant selective cytotoxicity for cancer cells while having minimal toxicity for normal cells. This may be a viable cancer treatment option (Devi et al., 2015; Roy et al., 2022).

Therefore, in our study, 2500 natural compounds of plant origin were collected and selected for their inhibitory properties. The compounds screened for docking were first put via in silico docking and high-throughput virtual screening before being chosen as the best inhibitor for the target proteins (NY-ESO-1, RUNX3, and UBE2Q1). The study of their interactions with the target proteins may assist in developing novel drugs for OC treatment.

Materials and Methods

Retrieval and Refinement of OC Proteins

Crystallized structures of NY-ESO-1, RUNX3, and UBE2Q1 were retrieved from Protein Data Bank (PDB) (Bittrich et al., 2022). Molecular Operating Environment (MOE) was used to prepare the retrieved structures for docking (Vilar et al., 2008; Wang et al., 2021).

Ligand Database Preparation

An extensive literature search was carried out to find phytochemicals that have been reported to be effective against OC. Phytochemical chemical structures were obtained from the MPD3 database, the Pubchem database, the MAPS database, and the Zinc database (Irwin & Shoichet, 2005; Kim et al., 2016; Mumtaz et al., 2017) in various ligand file formats, including mol, sdf, and mol2. These ligand structures were optimized with the Protonate3D module by adding partial charges in MOE. The MMFF94X force field was utilized to minimize each ligand’s energy. The ligands were then individually added to the MOE ligand database for docking purposes.

Binding Site Analysis

Active sites of the targeted receptors were predicted through the CPORT tool present in the Haddock interphase (de Vries & Bonvin, 2011). These active site pockets are where the ligand will bind to the target receptors to inhibit its activity. The sites consist of hydrophilic, hydrophobic, donor, acceptor and metal-binding domains. The best predicted binding pocket was chosen for performing molecular docking.

Molecular Docking

Within the specified docking sites of the NY-ESO-1, RUNX3, and UBE2Q1 proteins, the MOE Dock tool was utilized to dock a ligand database of 2500 phytochemicals. The default ligand placement method was used to find 1000 best poses of docked molecules by the triangular matcher algorithm (Podvinec et al., 2010; Alhumaydhi et al., 2021). The London dG scoring function was used to rescore simulated poses. The Force field refinement algorithm that applies the Generalized Born solvation model to calculate final binding energy while keeping receptor residues rigid was used to further minimize the top ten ranked poses per molecule created by London dG (Li et al., 2022). Binding affinity, S-score, and Root-Mean-Square Deviation (RMSD) values were used to rank all compounds. The goal was to select only those compounds for further investigation from the top-ranked poses that bind to active OC protein residues with a good dock score.

Analysis of Ligand- Receptor Interaction

To clearly understand the best-docked complex receptor-ligand interactions, MOE's LigX tool was utilized for analyzing ligand-receptor interactions in 2D plots (Khalifa et al., 2020; Alhumaydhi et al., 2021).  It creates electrostatic, hydrophobic, Van der Waals forces in the active site of OC proteins and a two-dimensional graph of hydrogen bonding, contributing to drug-like molecule affinity. Using UCSF Chimera, the 3D images of OC protein inhibitor complexes were created (Pettersen et al., 2004; Pettersen et al., 2021).

Drug scan/ADME Toxicity

A computational approximation of the drug probability of docked phytochemicals was found through the drug scanning tool on the Molinspiration server within limits set by Lipinski's Rule. Using the ADMETsar server, the deposition, absorption, excretion, metabolism, and toxicity profiles of these hits were virtually evaluated (Cheng et al., 2012). Carcinogenic potential and AMES toxicity of inhibitors were also evaluated.

Free Energy Calculation Binding

The free energies binding of receptor-ligand complexes was calculated by the Schrodinger suite's Prime module and the OPLS-2005 force field (Mishra & Singh, 2022). The free energy binding was calculated by subtracting the complex free energy from the ligand and protein free energies sum (Ikot et al., 2020).

ΔG(binding) = ΔG(complex) – [ΔG(ligand)+ ΔG(protein)]


where ΔG(binding) denotes the binding free energy and ΔG(complex), ΔG(ligand), and ΔG(protein) denote the free energies of the complex, ligand, and protein, respectively. Using the MM-GBSA method and the MM-PBSA.py program, with 250 snapshots taken at equal intervals over the last 20 ns of simulation, binding free energies were estimated.

MD Simulation

The Desmond software was used to run MD simulations (Release, 2017). To determine protein-ligand interactions, which were solvated using the simple point charged (TIP4P) water model, the optimized potentials for liquid simulations (OPLS)-2005 force field was used in this system. A 10Å buffer region was created among the box sides and protein atoms using the orthorhombic water box. Na+ ions neutralized the systems after overlapping water molecules were removed. The energy was calculated using the OPLS-2005 force field. The temperature was kept constant at 300 K during the integration step, yielding a 2.0 fs value. Finally, for monitoring the consistency of the OC proteins in their native motion, the Root Means Square Deviation (RMSD), Root Mean Square Fluctuations (RMSF), Solvent Accessible Surface Area (SASA), and Radius of Gyration (RoG) were calculated. For up to 100 ns, the synchronize file was saved every 5000 ps, and the result was examined using Nagasundaram et al.'s method (Kabir et al., 2022).

Results and Discussion

Proteins Preparation

The study began with retrieving the crystal structures of NY-ESO-1, RUNX3, and UBE2Q1 proteins from the PDB database with PDB IDs of 1S9W, 5W69, and 2QGX, respectively. X-ray diffraction was used to determine the structure up to a resolution of 2.20 Å, 2.80 Å, and 2.56 Å, respectively.  MOE was used to prepare the structures for docking. 3D protonation, removal of water molecules, and energy minimization were performed using MOE.

Binding Site Analysis

Binding pockets against the target receptors were chosen from Cport. For receptor NY-ESO-1, the reported active residues were Tyr A27, Asp A30, Thr A31, Gln A32, Arg A48, Ser B52, Asp B53, Ser B55, Tyr B63, Leu B65, Tyr B67, Phe A241. Similarly, for RUNX3, it was Tyr D26, Tyr C27, Gln C32, Ser D52, Asp D53, Ser D55, Tyr D63, Leu D65 and Phe C241. While Val A32, Asp A38, Asn A71, Ser A73, Arg A 86, Ser C147, Leu A148, and Thr C149 were the residues reported for UBE2Q1.

Molecular Docking

OC proteins were docked against the phytochemical ligand database. A strict filter that included four factors was used to rank docked compounds: hydrogen bonding strength, maximum binding pocket occupancy with lowest Gibbs free energy, and other potential non-covalent interactions, all calculated and represented using an S-score function. Top docking poses were chosen from among 10,326 docked molecules. The ranking criteria included threshold-based criteria that required a ligand to exhibit the desired S-score values (stronger the affinities and interaction, lower the score) and bind to the target OC proteins using all of the hotspot conserved residues in the binding pocket (Table 1).



Table 1. Docking statistics of target proteins against plant compounds

Compounds I’D

Compounds Name

Binding Affinity


Interacting residues




-17.65 kcal/mol


Asp A30, Thr A31, Gln A32, Leu B65, Ser B52, Tyr A27, Tyr B63



-15.34 kcal/mol


Tyr A27, Ser B55, Asp B53, Ser  B52, Tyr B63, Phe A241, Gln A32



-14.76 kcal/mol


Tyr A27, Gln A32, Arg A48, Tyr B67, Ser B52, His B51, Asp A30




-16.65 kcal/mol


Tyr C27, Ser D55, Tyr D63, Phe C241, Asp D53, Ser D52, Gln C32, Leu D65, Tyr D26


Buddledin C

-15.76 kcal/mol


Tyr D26, Thr C233, Ala C211, Tyr D63, Pro C235



-13.87 kcal/mol


Arg C48, Tyr D67, Phe C241, Tyr C27







Asp A38, Val A32, Leu A148, Ser A73, Asn A71, Arg A 86






Lys A75, Asp A38, Ser A73, Ser C147, Thr C149


Vibsanin H




Ser A73, Leu C148, Asp B34, Asn A33, Arg C86



Humulene, Thiamine, and Deoxyactein were found to bind with high affinity within the active sites of the NY-ESO-1 (Figure 1). Humulene was bound to NY-ESO-1 protein, having a score of -17.65 KJ / mol creating hydrogen bonds with Asp A30, Thr A31, Gln A32, Leu B65, Ser B52, Tyr A27, and Tyr B63 side chains (Figure 1a) and Thiamine D was bound having a binding score of -15.34 kcal/mol creating hydrogen bonds with Tyr A27, Ser B55, Asp B53, Ser B52, Tyr B63, Phe A241, and Gln A32 side chains (Figure 1b). In contrast, Deoxyactei was bound to NY-ESO-1 protein with a -14.76 KJ / mol, forming hydrogen bonds with the side chains of Tyr A27, Gln A32, Arg A48, Tyr B67, Ser B52, His B51, and Asp A30 (Figure 1c).







Figure 1. Interaction mechanisms and binding modes of new NY-ESO-1 protein inhibitors Interaction analysis of Humulene (a), Thiamine (b), and Deoxyactein (c)



Likewise, (3-(Benzyloxy) isoxazol-5-yl) methanol, Buddledin C, and Terthiophene were observed to have a high affinity within the active sites of the RUNX3 (Figure 2). (3-(Benzyloxy) isoxazol-5-yl) methanol was bound to RUNX3 protein, having a score of -16.65 KJ / mol creating hydrogen bonds with the side chains of Tyr C27, Ser D55, Tyr D63, Phe C241, Asp D53, Ser D52, Gln C32, Leu D65, and Tyr D26 (Figure 2a) and Buddledin C was bound with a binding score of -15.76 kcal/mol creating hydrogen bonds with side chains of Tyr D26, Thr C233, Ala C211, Tyr D63, and Pro C235 (Figure 2b). While Terthiophene was bound to RUNX3 protein with a -13.87 KJ / mol score, creating hydrogen bonds with the side chains of Arg C48, Tyr D67, Phe C241 and Tyr C27 (Figure 2c).