Targeting Ebola Virus VP30 with Quercetin Derivatives: An Integrated Molecular Docking, MD Simulation, ADMET, and DFT Study

Abstract

Ebola virus disease still poses a major threat to global public health because of its high case-fatality ratio and lack of antiviral medicines. Therefore, in this study, seventeen quercetin derivatives (D1-D17) were computationally screened against the Ebola virus VP30 protein, which acts as a transcription activator in the Ebola virus life cycle. D1 and D2 derivatives displayed high binding energies (-9.1 and -8.7 kcal mol^-1, respectively), showing the formation of hydrogen bonding, π-sigma, and π-alkyl interactions with the vital amino acid residues in VP30 protein. Stability of the complexes was further assessed by conducting molecular dynamics simulation for 100 ns based on root mean square deviation, root mean square fluctuation, radius of gyration, principal component analysis, and molecular mechanics/poisson-Boltzmann surface area calculations. The results revealed that VP30-D1 and VP30-D2 complexes exhibited high binding free energies with the values of -22.1562 and -19.7180 kcal mol^-1, respectively. The pharmacokinetic and safety analysis predicted that D1 and D2 have good drug-like and safety profiles, while their bioavailability is moderate with the possibility of nephrotoxicity. Density functional theory calculations showed that both analogues have favorable electronic structures, such as a low HOMO-LUMO gap (highest occupied molecular orbital and lowest unoccupied molecular orbital, respectively), electrophilicity index, and high electron delocalization.

Keywords:
Ebola virus; VP30; molecular docking; hemorrhagic disease; quercetin

Introduction

In 1976, the world attention was drawn to an outbreak of a mysterious virus known as Ebola in Sudan (now the Democratic Republic of the Congo).¹ In 1989, Ebola Reston was first reported in monkeys imported from the Philippines to the USA.² The Ebola virus (EBOV) belongs to the Filoviridae class, a name derived from the Latin term for thread-like.³ The number of cases during the 2014 outbreak surpassed the total of all previously reported EBOV outbreaks combined.⁴ As a result of the infection, the most commonly reported symptoms included fever, diarrhea, anorexia, fatigue, nausea, vomiting, joint pain, muscle pain, and headache.⁵ Initially, patients experience symptoms such as headaches, muscle pains, weakness, and problems related to the digestive system. Afterward, the patient might also develop liver problems.⁶ EBOV is 75 80 nm in diameter with a smooth surface and has filamentous structures of about 970 nm in length, which are highly infective. The viral covering, which is from the host cell membrane, surrounds a nucleocapsid of 45 60 nm diameter, and spike features originate from the viral envelope.⁷ Bats have been identified as one of the most common natural hosts of EBOV.⁸ Observation of the virus through electron microscopy has revealed that the viral particles have a typical filamentous shape, with some of them being configured into the shape of the number 6. The diameter of these filaments is about 80 nm, and their length may vary from 800 to 1,000 nm. It is significant to mention that only a small fraction, 1%, of their total mass is made up of ribonucleic acid (RNA).⁹ Understanding of virus-host interaction and the mechanisms of evasion of host antiviral response by the virus are of utmost importance in elucidating the pathogenesis of EBOV disease. Viral protein 35 (VP35) is considered a key component of the EBOV RNA replication complex.¹⁰ During viral reproduction, the polyprotein is sliced into three crucial essential proteins: capsid, membrane, and envelope proteins, which are integral to the assembly of viral units. Additionally, seven non-structural proteins (NS1, NS2A, NS2B, NS3, NS4A, NS4B, and NS5) contribute to viral reproduction, virion assemblage, and evasion of host defence mechanisms.¹¹ Little is known about viral accumulation and the control of EBOV replication due to the limited number of EBOV open reading frames.¹² Structure-guided EBOV viral protein 30 (VP30) variants have shown changes in resemblance to the natriuretic peptides (NPs). The relationship between VP30-NP interaction and the activity of these variants in a cell based tiny replicon assessment suggests that the VP30-NP interaction plays both vital and hindering roles in EBOV RNA synthesis.¹³

EBOV VP30 is a critical factor in viral transcription and is strongly associated with the nucleic acid components in viral particles.¹⁴ The translation activator VP30 is thought to be crucial in EBOV replication, probably by stabilizing the initiation of mRNA synthesis. Structural and mutagenesis studies reveal a potential target for small-molecule inhibitors that could block VP30 action and thereby impede viral spread.¹⁵

Seven proteins encoded by EBOV include nucleoprotein (NP), glycoprotein (GP), RNA-dependent RNA polymerase (L protein), and four viral proteins: VP24, VP30, VP35, and VP40.¹⁶ To fully understand the pathogenesis of EBOV infections, it is essential to grasp how the virus interacts with its host and how it suppresses the host antiviral response.¹⁰ Because of the stringent biosafety requirements for EBOV, the direction of drug discovery and other developmental processes with such a virus is highly challenging. Either way, EBOV research has become more and more attractive due to recent breakthroughs.¹⁷ Due to the high pathogenicity of EBOV, live virus vaccines are ideally not considered a very appropriate choice for vaccination.⁷ While effective drugs have been identified for other viral diseases, there are no recently approved therapeutics for curing EBOV infections.¹⁸ Most likely, during host cell infection, the interaction between EBOV proteins and several host proteins results in viral replication.¹⁹

Quercetin, first identified by Nishino et al.²⁰ as a potent anticancer compound, is one of the most broadly distributed flavonoids throughout the plant world. It has several biological activities, including anti-inflammatory, antioxidant, antibacterial, radical scavenging, antiviral, immune modulatory, and gastroprotective properties. Besides, it also has the potential to be used for the treatment of obesity, cardiovascular problems, and diabetes.²¹ On the one hand, animal experiments reveal that quercetin performs a wide spectrum of biological functions both in vitro and in vivo. On the other hand, it contributes to lowering blood pressure and triggers vasodilation both systemically and in the coronary arteries.²² As a dietary flavonoid, quercetin is of considerable interest as a novel medicinal biomolecule with versatile therapeutic properties.²³ Quercetin, a type of flavanol, is widely studied for its significant and varied effects on different types of cells. However, its mechanism of action is not fully understood. It has been shown to exhibit antiproliferative effects and induce cell death in cancer cell lines, mainly through an apoptotic mechanism.²⁴ Isorhamnetin inhibits tyrosinase and functions as an anti-pigmentation agent.²⁵ Flavonoids are notable phytochemicals with therapeutic properties, found in fruits and vegetables as secondary metabolites. To date, approximately 9,000 types of flavonoids have been identified in natural foods. It has been reported that less than 8% of ingested quercetin reaches the kidneys and liver in rodents. In food, quercetin primarily exists in complex forms, combined with carbohydrates, phenolic acids, and alcohol. After ingestion, quercetin derivatives are mainly hydrolyzed in the gastrointestinal tract before being absorbed and utilized.²³ Quercetin is found in a conjugated form in human blood, primarily as a glycoside. Boulton et al.²⁶ found that when quercetin is bound to plasma proteins such as albumin, it accounts for 99.4% of its presence, which reduces its bioavailability to cells. Quercetin exhibits both neuroprotective and neurotoxic effects. Studies have shown that quercetin acts as a neuroprotectant in rat brains, especially when used in conjunction with fish oil. It demonstrates beneficial effects against neurodegenerative diseases, such as Alzheimer’s disease, where it inhibits acetylcholinesterase.²⁷

EBOV VP30 is a multifunctional protein that also acts as a vital transcriptional activator,¹⁴ which has led to the idea that it could be a potential target for an antiviral cure.¹⁵ VP30 contains a zinc-binding spot, phosphorylation spots, and an RNA-binding spot in its N-terminal segment. VP30 is highly phosphorylated, and this posttranslational modification regulates its activity during transcription.¹⁴ As a key transcriptional activator, VP30 is responsible for the transcription and replication of viral RNA. It establishes interactions with the NP and polymerase of the virus and initiates transcription of the viral mRNA that leads to viral protein synthesis. Inhibition of VP30 prevents the assembly of functional ribonucleoproteins and thus suppresses transcription of genes involved in EBOV replication.²⁸ Given its critical role in viral transcription, VP30 acts as a promising candidate for viral therapeutics. It is involved in the N-terminal oligomerization of VP30, which can be inhibited by generating peptides that induce inactive forms of VP30, thereby preventing viral replication.¹⁵ EBOV comprises nucleoproteins, structural proteins, non structural proteins, and glycoproteins. The glycoproteins consist of two vital subunits, glycoprotein subunit 1 (GP1) and glycoprotein subunit 2 (GP2), which play a significant role in interacting with host cells. The GP1 subunit features potential receptor-binding regions and is crucial for cell attachment, whereas the GP2 subunit facilitates the union of viral and cellular membranes. VP30 is a phosphoprotein essential for EBOV transcription. It has been shown that VP30 plays a role in the reinitiation of gene transcription and is regulated by modifications at its phosphorylation sites.²⁹

In this work, we investigate the potential of quercetin derivatives as inhibitors of the EBOV VP30 protein. The selection of quercetin analogues for this study was driven by the fact that flavonoids exhibit extensive antiviral effects against various RNA viruses as a result of their inhibition of viral enzymes, transcriptional factors, and replication-associated proteins.³⁰ Due to the presence of multiple hydroxyl moieties and aromatic rings in its structure, quercetin polyphenolic framework can bind to viral proteins via hydrogen bonding and π-stacking interactions through polar and charged amino acid residues.³¹ As a key transcription factor involved in the replication process and formation of hexamers, there is the presence of positive charges and polar amino acid residues in EBOV VP30, which can facilitate binding to the flavonoid framework.³² Additionally, previous experimental and computational studies on the antiviral effects of flavonoid-based compounds against Ebola-related viral targets support the rationale for assessing the quercetin derivatives described in this study as possible inhibitors of VP30 protein function.³³

We assess the binding interactions between these flavonoid derivatives and VP30, along with the stability of the resulting complexes via molecular docking and molecular dynamics (MD) simulations. additionally, adsorption distribution metabolism, excretion, and toxicity (ADMET) predictions were carried out to evaluate the pharmacokinetics and toxicity of the most active quercetin derivatives. Besides, density-functional theory (DFT) calculations were performed to facilitate an understanding of the essential information about the electronic structure, charge transfer, and the biological activity of the selected quercetin derivatives. The outcomes of this investigation underscore the promising ability of the selected compounds to interfere with VP30 hexamerization, a key step in viral reproduction, laying a foundation for future experimental studies aimed at developing effective therapeutics for EBOV infections.

Methodology

Molecular docking

The molecular docking method was applied to evaluate the interaction between quercetin derivatives (D1-D17) with the EBOV VP30 protein (PDB ID: 5T3T).³⁴ The protein with PDB ID: 5T3T was chosen because it corresponds to the VP30 protein of the EBOV, which acts as an important transcription factor involved in RNA transcription and replication. Since VP30 plays an important part in the expression of genes within the virus, inhibition of the protein would prevent the proliferation of the virus, thus making it a promising target for anti-virus therapy.³⁵ Firstly, the structural data of the VP30 protein (pdb format) were extracted from the Protein Data Bank.³⁶ Subsequently, ChemDraw³⁷ was used to draw the structures of the ligands (D1-D17), and the conformations were generated in three dimensions with Chemical Three-Dimensional (Chem3D, version 19.0) modelling software. Ligand structures were geometry-optimized with the molecular mechanics force field 2 (MM2 force field) and saved in Mol2 format for docking studies. The EBOV VP30 protein (receptor structure) was prepared before docking by the removal of crystallographic waters, co-crystallized ligands, and heteroatoms using molecular operating environment (MOE 2022).³⁸ Polar hydrogens were added to the structure, and the protein was energy minimised. Subsequently, the Cavity Detection-Guided Blind Docking (CB-Dock) web server³⁹ was utilized to conduct docking studies. The free CB-Dock web server uses AutoDock Vina⁴⁰ as the docking program. Initially, the pre-processed ligands and receptor molecules were uploaded in order to recognize the binding sites (based on the complementarity of shapes) by CB-Dock. Docking analysis was performed by using the cavity-detection blind docking method available in the CB-Dock2 server, which predicted the possible binding cavities on the protein surface and docked with the help of AutoDock Vina. Prepared protein structure and best ligands were uploaded to the server, and default docking parameters were used. CB-Dock automatically identified the top binding cavities according to their size and center coordinate, and automatically performed the docking in these predicted binding active sites. Finally, the binding poses along with their corresponding binding energies were acquired for further studies. The docking results were analyzed by using Discovery Studio Visualizer 2025.⁴¹ For all the ligands, docking poses were ranked based on binding energies calculated in kcal mol^-1. Docking conformation with minimum binding energy and better intermolecular forces was chosen as the best binding pose for subsequent analysis. Further investigations of the docked complexes were done to determine the possible interactions between the receptor and ligand molecules, including hydrophobic, hydrogen bonds, and electrostatic interactions. Interaction images were constructed to predict important amino acid residues in binding.

The chemical structures and nomenclature of the quercetin derivatives (D1-D17) selected for this study are summarized in Table 1.

Table 1
Codes, chemical names, and structures of quercetin derivatives (D1-D17) investigated in this study

Code	Quercetin derivative	Chemical structure
D1	2-(3,4-dihydroxyphenyl)-3,5-dihydroxy-7-(((2S,3S,4S,5S,6S)-3,4,5,6-tetrahydroxytetrahydro-2H-pyran-2-yl)oxy)-4H-chromen-4-one
D2	2-(3,4-dihydroxyphenyl)-5,7-dihydroxy-3-(((2S,3S,4R,5R,6S)-3,4,5-trihydroxy-6-(hydroxymethyl)tetrahydro-2H-pyran-2-yl)oxy)-4H-chromen-4-one
D3	2-(3,4-dihydroxyphenyl)-5,7-dihydroxy-3-(((1R,2S,3R,4S,5R)-2,3,4-trihydroxy-5-(1-hydroxyvinyl)cyclohexyl)oxy)-4H-chromen-4-one
D4	2-(3,4-dihydroxyphenyl)-5,7-dihydroxy-3-(((2S,3R,4S,5S,6R)-3,4,5-trihydroxy-6-(hydroxymethyl)tetrahydro-2H-pyran-2-yl)oxy)-4H-chromen-4-one
D5	3,5,7-trihydroxy-2-(4-hydroxy-3-methoxyphenyl)chroman-4-one
D6	3,5,7-trihydroxy-2-(4-hydroxyphenyl)chroman-4-one
D7	(2R,3R)-3,5,7-trihydroxy-2-(4-hydroxy-3-methoxyphenyl)chroman-4-one
D8	3,5,7-trihydroxy-2-(3-hydroxy-4-methoxyphenyl)chroman-4-one
D9	2-(3,4-dihydroxyphenyl)-3,5,7-trihydroxychroman-4-one
D10	2-(3,4-dihydroxyphenyl)-5,7-dihydroxy-3-methylchroman-4-one
D11	5,7-dihydroxy-2-(4-hydroxyphenyl)-3-methoxychroman-4-one
D12	2-(3,4-dihydroxyphenyl)-5,7-dihydroxy-3-methoxychroman-4-one
D13	(2S,3S)-2-(3,4-dimethoxyphenyl)-3,5,7-trimethoxychroman-4-one
D14	2-(3,4-dihydroxyphenyl)-5,7-dihydroxychroman-4-one
D15	7-hydroxy-2-(3-hydroxyphenyl)-3-methoxychroman-4-one
D16	(4S)-4-amino-5-(2-(3,4-dihydroxyphenyl)-3,6,8-trihydroxy-4-oxochroman-7-yl)-5-oxopentanoic acid
D17	(4R)-4-(((5,7-dihydroxy-2-(4-hydroxy-3-methoxyphenyl)-4-oxochroman-3-yl)oxy)amino)-5-oxopentanoic acid

MD simulations

MD simulation studies were carried out on the two highest-ranking docked complexes using derivatives D1 and D2, which were referred to as COM-1122 (VP30-D1 complex) and COM-1133 (VP30-D2 complex), respectively. The protein-ligand complex was built, and energy minimized using SYBYL-X 2.1.1 (Certara, SYBYL-X Suite, version 2.1.1, 2012).⁴² For the apo system, the crystal structure of VP30 was downloaded from the Protein Data Bank (PDB ID: 5T3T).³⁵ PyMOL 1.8.⁴³ was used for structural visualization and inspection. All the simulations were done using the Assisted Model Building with Energy Refinement 12 (AMBER12) software package,⁴⁴ and for the addition of hydrogen atoms as well as achieving neutralization, the LEaP module was used. In order to neutralize the systems before the solvent environment is applied, sodium (Na⁺) ions were introduced. The number of Na⁺ ions added in the apo system was 6, while in the ligand-bound complex system, 7 Na⁺ ions were introduced into each ligand-bound complex system. All systems were placed in a solvent box with a truncated octahedral geometry filled with transferable intermolecular potential with 3 points (TIP3P) water molecules,⁴⁵ where the minimal distance between the solute and the boundary is 10 Å. The simulation box size for the apo system was 82 × 82 × 82 Å, while for the ligand-complexed protein, the simulation box size was 84 × 84 × 84 Å. The number of water molecules present in the apo system was approximately 21,500, while that of the COM-1122 and COM-1133 complexes was 22,100 and 22,300, respectively. The van der Waals and direct space electrostatic nonbonded interactions were truncated at a value of 8 Å, which is in line with conventional AMBER simulations. The long range electrostatic interactions were accurately calculated via the Particle Mesh Ewald (PME) method.⁴⁶ All the simulations had been done with the ff12SB force field.⁴⁷ Temperature was regulated during equilibration using a Langevin thermostat with a collision frequency of 1 ps^-1, while production runs employed a Berendsen thermostat.⁴⁸ Similarly, all the simulations were sped up with a Compute Unified Device Architecture (CUDA) enabled Particle Mesh Ewald Molecular Dynamics (PMEMD) engine on NVIDIA Tesla K20 GPUs. Initially, energy minimization with the solute restrained was done to remove unfavorable steric contacts in the solvated systems. This step was made up of 2,000 conjugate gradient minimization steps followed by 500 steepest descent minimization steps. This step was then followed by another minimization of the whole system without any restraints, 1,000 steepest descents, and 19,000 conjugate gradient minimizations. After the energy minimization procedure was done, each system was heated from 0 to 310 K in increments and equilibrated for 400 ps at the NVT level before the production stage in the MD simulations. Then, pressure equilibration was achieved under NPT (constant number of particles, pressure, and temperature) conditions at 1 bar and 310 K for 1 ns using isotropic position scaling to stabilize the system density before the production phase of the MD simulations. The average density of the equilibrated systems stabilized near 1.0 g cm ³. MD simulations of both apo and protein-ligand complex systems were compared for 100 ns at 310 K. The coordinates of the trajectory were saved every 1 ps. Trajectory analyses were performed using the Processing TRAJectory (PTRAJ) and C++ Processing TRAJectory (CPPTRAJ) modules of AMBER. The changes in structure and movement patterns of proteins were monitored by calculating root mean square deviation (RMSD), root mean square fluctuation (RMSF), radius of gyration (R_g), potential energy (P.E.), and H-bonding analyses.⁴⁹

To validate the simulation results accurately, structural analyses were carried out on the equilibrated segments of the MD trajectories. RMSD, RMSF, and Rg were calculated as time-averaged parameters over the stable regions of the trajectories, whereas the estimation of molecular mechanics/Poisson-Boltzmann surface area (MM/PBSA) binding free energies was done using multiple frames randomly selected from the production run last 50 ns. This approach of trajectory sampling and averaging helps in reducing the statistical noise and, at the same time, offers a more trustworthy depiction of the ligand-protein interaction stability throughout the simulation.

MM/PBSA calculations

Binding free energy estimations were executed via the MM/PBSA technique. This technique can be used for analyzing the non-bonding interactions between ligands and receptors.⁵⁰ The MM/PBSA free energy estimation, combined with MD simulations, provides quantifiable calculations of VP30-ligand binding energies. The binding free energy (∆G_bind) is calculated by equation 1.

(1) $$\mathrm{\Delta }{G}_{\text{bind}}=G_{R+L}-(G_R+G_L)$$

where G_R+L, G_R, and G_L show the free energies of the VP30 complex, free VP30, and the ligand, respectively. In the MM/PBSA method, each free energy term of equation 1 is determined as described by equation 2.

(2) $$\mathrm{G}={\mathrm{E}}_{\mathrm{bond}}+{\mathrm{E}}_{\mathrm{vdw}}+{\mathrm{E}}_{\mathrm{elec}}+{\mathrm{G}}_{\mathrm{PB}}+{\mathrm{G}}_{\mathrm{SA}}-{\mathrm{TS}}_{\mathrm{S}}$$

The bond energy, consisting of bond, angle, and dihedral energy, is basically represented by E_bond (equation 2). The Van der Waals energy component is written as E_vdw, while the electrostatic energy is indicated by E_elec. The solvation energy of polar and nonpolar parts is indicated by G_PB and G_SA, respectively, whereas T denotes the absolute temperature, and S_S stands for solute entropy. Calculation of the TS_S term was not done in the MM/PBSA approach due to high computational expense and insignificant contribution to the analysis of the binding energy trend. The PBSA model in MMPBSA was applied for binding free energy calculation. Binding free energy computations were executed by using only the last 50 ns of MD trajectories from the equilibration phase.

Principal component analysis (PCA)

The PCA relies on the covariance matrix and measured values to visualize atomic displacements in proteins. Before PCA execution, solvents and ions were removed from the apo and complex structures via the PTRAJ unit. The structures were then positioned based on the stripped coordinates. PCA was conducted, and the covariance matrix was created for the first three key modules (PC1-PC3).⁵¹

ADMET predictions

Pharmacokinetics profile, drug likeness, toxicity, physicochemical properties, and pharmacodynamics profile of the prospective drug molecule form an important step in the drug discovery process.⁵² These analyses will help to determine the overall drug potential of the compound. Hence, the factors such as solubility in water, lipophilicity, Lipinski’s rule of five compliance, topological polar surface area (TPSA), blood-brain barrier (BBB) penetration, and toxicity are considered for determining drug likeness.⁵³ The ADMET properties of the selected quercetin derivatives were predicted using in silico methods. Critical descriptors such as physiochemical properties, aqueous solubility coefficient (Log S), lipophilicity (Log Po/w), pharmacokinetics, and drug-likeness according to the Lipinski criteria were obtained from the SwissADME server.⁵⁴ Each physicochemical and pharmacokinetic descriptor has been selected for its unique property related to drug likeness and does not overlap with any other descriptor. Evaluation of adherence to Lipinski’s rule of five was done based on numerical values, including (i) molecular weight (MW) ≤ 500 Da, (ii) lipophilicity (Log P) ≤ 5, (iii) hydrogen bond donor (HBD) ≤ 5, and (iv) hydrogen bond acceptor (HBA) ≤ 10. According to the rule, compounds violating more than one of these four criteria may show poor absorption or permeation properties.⁵⁵ These numerical parameters are normally accepted criteria for permeability and drug-like properties of orally administered drugs; nevertheless, these parameters are not considered reliable measures for pharmacological efficacy or bioavailability of molecules. Notably, molecules that violate one or more of the Lipinski parameters can still be highly effective drug candidates, and their potential should not be dismissed solely based on this rule.⁵⁶ Toxicity was then predicted using the pkCSM web server⁵⁷ and ProTox 3.0 web server.⁵⁸

DFT analysis

DFT is another quantum mechanical method that has encouraged much interest as a very powerful tool in drug discovery. Its key features are low computational cost and accuracy in the energetics of molecular and biological systems, enabling it for most of the medicinal purposes.⁵⁹ The 3D structures of the selected quercetin derivatives were optimized using Gaussian 05.⁶⁰ The HOMO (highest occupied molecular orbital) and LUMO (lowest unoccupied molecular orbital), molecular electrostatic potential (ESP) surfaces, Mulliken charge distributions, and natural bond orbital (NBO) visualization plots were generated and analyzed using GaussView 5.0.8.⁶¹ Single point energy computations were done by B3LYP density function, which is a combination of Becke’s three-parameter hybrid exchange function⁶² and the Lee-Yang-Parr correlation function.⁶³ All computations were done with the 6-311G (d,p) basis set. The choice of the B3LYP density function is based on the fact that it is a proven method of calculating DFT-derived geometric parameters with accuracy compared to experimental findings.⁶⁴ Further clarification of electronic properties was provided by the use of the frontier molecular orbital (FMO) theory, which was used to analyze the chemical stability, reactivity, chemical potential (μ), electronegativity (χ), electrophilicity index (ω), hardness (η), softness (S), and dipole moment (μ) of the selected quercetin derivatives. Additionally, the ESP surface of both compounds was drawn in order to determine their electrostatic charge distribution in three dimensions. A Mulliken population analysis was undertaken to determine how the electrons are distributed across the various atoms in the selected compounds to provide more understanding of the ligand-protein interactions. The NBO population analysis was also undertaken to determine the degree of electron delocalization, donor-acceptor interactions, charge transfer within the molecule, and general stability of the quercetin derivatives.

Results and Discussion

Docking was executed between the EBOV VP30 and the selected quercetin derivatives.

Structural features of EBOV VP30 CTD protein bound to NP (PDB ID: 5T3T)

The crystal structure of the EBOV VP30 C-terminal domain (CTD) in complex with NP (PDB ID: 5T3T) (Figure 1) was chosen for molecular docking analysis owing to its high quality and functional role in viral transcription. The crystal structure has been reported via X-ray diffraction analysis at a resolution of 2.20 Å. It is worth mentioning that the validation results for this structure were satisfactory, since the values of the R-work and R-free parameters equaled 0.214 and 0.250, respectively, while there were no Ramachandran outliers in the molecule. The macromolecule had a mass of 220.89 kDa and included 11,795 atoms with 1,930 deposited residues. There were five protein chains in the molecule under consideration (A, B, C, D, and E), representing the fusion protein of NP and the minor NP VP30; their lengths were the same and equal to 193 amino acid residues. Architecturally, the VP30 CTD protein adopted an α-helical structure having a defined NP-binding pocket that can stabilize ligand binding using hydrogen bonds and hydrophobic interactions, thus making it a potential candidate for conducting structure-based docking against the EBOV.⁶⁵

Figure 1
Crystal structure of the C-terminal domain (CTD) of the EBOV VP30 bound to NP (PDB ID: 5T3T). The protein structure has an α-helical fold with the binding site of VP30 and NP marked. The crystallographic structure of the protein has been reported at a resolution of 2.20 Å and comprises five amino acid chains (A-E), with a molecular mass of 220.89 kDa and 11,795 atoms.⁶⁵

Domain flexibility and functional hotspots

The sequence-based structural annotation of EBOV VP30 CTD (Figure 2) revealed principally an ordered domain with very low levels of disorder. So, VP30 CTD represents a structurally stable receptor with an accessible ligand-binding surface suitable for molecular docking. Hydropathy plots revealed alternating hydrophilic and hydrophobic regions within the structure. Such regions would be expected to generate accessible solvent interaction surfaces and ligand binding grooves. The high disordered binding propensity in combination with low levels of disorder, within a few of the most conserved regions of the protein, indicated functional hotspots of VP30 and NP interaction. All these features made VP30 CTD an attractive target for structure-based design of inhibitors through molecular docking.

Figure 2
Sequence-based structural profiling of the EBOV VP30 CTD. Profiles with relative hydropathy and inherent structural disorder and disorder-based binding activity of the site are shown. The domain seems generally ordered and stable with localized hydrophobic regions and binding sites that could be exploited as possible hot-spots by potential ligands to improve specificity of binding and/or facilitate VP30 and NP recognition. A disordered binding site does not rule out the ability of the protein to dock with drug molecules and perform structure-based drug design experiments.

Protein structure validation

Stereochemistry of EBOV VP30 Structure (Protein Data Bank ID: 5T3T) was assessed using Protein Structure Checking Program (PROCHECK, v.3.5) based on the Ramachandran plot (Figure 3). It showed that stereochemistry was excellent because 1,200 residues out of 1,223 non-glycine and non-proline (98.1%) were found in most-favored regions, while the other 23 residues (1.9%) occupied other allowed regions. There were no residues detected in generously allowed and disallowed regions, signifying the absence of unfavorable backbone conformations. These values are considerably greater than those specified by the PROCHECK quality criteria (> 90% residues in the most favored regions), indicating that the crystal structure has been determined reliably. Also, the residue-by-residue stereochemistry analysis did not indicate any significant local deviations in side-chain dihedral angles. On the whole, the great Ramachandran statistics and stereochemistry of the VP30 confirm that it has a very stable and energy-favorable conformation, making it suitable for further molecular docking and structure-based drug design studies.

Figure 3
PROCHECK Ramachandran plot of the structure of VP30 Protein from EBOV (PDB code: 5T3T). Overall, 98.1% of the residues lie in the most favorable regions, whereas 1.9% lie in additionally allowed regions, and no residue is found to lie in generous or disallowed regions, validating the high-quality stereochemistry and feasibility of performing molecular docking/MD studies on this protein.

Results of molecular docking

According to studies, the binding energies below -4.25 kcal mol^-1 suggest that the ligand-receptor binding is effective, the energy below -7.0 kcal mol^-1 indicates a strong binding, and the energy below -5.0 kcal mol^-1 indicates good ligand-receptor binding.⁶⁶

Moreover, the intermolecular interactions between ligands and the target protein provide important information for the stability and binding affinity of the target docked complexes. Among these interactions, hydrogen bonds play a key role in molecular recognition. Generally, shorter hydrogen-bond distances correspond to stronger interactions and better stabilization of the protein-ligand complex, provided that favorable donor-acceptor geometries are preserved. Based on donor-acceptor distance, hydrogen bonds are categorized into: strong hydrogen bonds (2.2 2.5 Å), which are dominantly covalent; moderate hydrogen bonds (2.5-3.2 Å), which are dominantly electrostatic; and weak hydrogen bonds (> 3.2 Å), which are dominantly electrostatic. Each hydrogen-bonding interaction was evaluated using these hydrogen bond categories, which were considered to contribute to the stability of the ligand.⁶⁷-⁶⁹ Besides hydrogen bonding, hydrophobic and aromatic interactions such as alkyl, π-alkyl, π-π stacked, π-σ, amide-π stacked, and π-π T-shaped contacts also play an important role in ligand-protein binding through favorable van der Waals and dispersion forces. These noncovalent interactions are usually found at varying intermolecular distances, depending on the nature of the interacting groups. They also play a role in binding affinity by promoting optimal ligand orientation and hydrophobic complementarity within the binding pocket.⁷⁰,⁷¹

Based on the molecular docking analysis performed on the interaction of D1-D17 with the EBOV VP30, all of the derivatives were found to create several stabilizing interactions with the receptor within its binding site (Figure S1 in the Supplementary Information (SI) section). Specifically, it was observed that the stabilizing interactions between the derivatives and VP30 are created mostly by means of hydrogen bonds and hydrophobic interactions with some key amino acids such as Arg179, Lys180, Phe181, Ser182, Lys183, Ser184, His215, Lys218, Gly219, and Ser254. The results revealed substantial differences in affinity and interaction patterns among the investigated derivatives, with binding energies ranging from -9.1 to -5.1 kcal mol^-1. The observed variations in binding energies were closely associated with the number of interactions formed, the nature of these interactions, the physicochemical characteristics of the participating amino acid residues, and the corresponding interaction distances.

Derivative D1

Among all derivatives, derivative D1 exhibited the highest affinity toward VP30, with a binding energy of -9.1 kcal mol^-1 and established a total of eleven intermolecular interactions. The complex was primarily stabilized by three hydrogen bonds with Arg179 (3.16, 3.18, and 4.42 Å), with the first two representing moderate hydrogen bonds that contributed substantially to binding. Additional hydrogen bonds with Ser184 (3.15 and 4.23 Å), Glu252 (4.92 Å), and Lys218 (5.60 Å) further reinforced complex stability. D1 also formed hydrophobic and aromatic interactions, including π-σ interactions with Arg179 (5.38 Å) and His215 (6.22 Å), π-alkyl interaction with Lys183 (6.25 Å), and a π-π T-shaped interaction with His215 (5.96 Å), which enhanced ligand orientation and hydrophobic complementarity within the binding pocket. The combination of multiple hydrogen bonds and aromatic contacts involving key residues such as Arg179, Ser184, Lys183, and His215 accounts for the strong affinity and stability of the D1-VP30 complex.

Derivative D2

The second-best affinity was observed for D2 with a binding energy of -8.7 kcal mol^-1. It established nine intermolecular interactions with VP30. The strongest interactions were hydrogen bonds with Gln212 (2.71 Å) and Gln233 (2.75 Å), which are relatively short and contribute considerably to complex stabilization. Additional hydrogen bonds with Ser187 (4.08 and 4.88 Å) and Ser254 (4.66 Å) provided weaker stabilizing effects. Hydrophobic stabilization was achieved through π-alkyl interactions with Lys183 and Lys218 (4.91 and 4.85 Å), while a C-H bond interaction with Ser182 (4.54 Å) contributed to ligand positioning. However, D2 also exhibited an unfavorable positive-positive interaction with Lys218 (4.61 Å), which may reduce overall binding stability. Although D2 formed two moderately strong hydrogen bonds, the lower number of interactions and the presence of electrostatic repulsion likely explain its slightly lesser affinity for the target protein compared with D1.

Derivative D3

D3, with a binding energy of -7.8 kcal mol^-1, established six interactions and was characterized by the dominant participation of Arg179. Two hydrogen bonds with Arg179 were observed at distances of 4.42 and 5.84 Å, accompanied by a π-alkyl interaction at 5.70 Å. A hydrogen bond with Glu252 (5.05 Å) and an alkyl interaction with Lys218 (3.90 Å) further stabilized the complex. The relatively short hydrophobic interaction with Lys218 contributed to the binding stability, whereas the longer hydrogen bond with Arg179 contributed mainly to ligand orientation. The repeated involvement of Arg179 highlights its importance as a key anchoring residue within the VP30 binding site.

Derivative D4

D4 exhibited a binding energy of -7.5 kcal mol^-1 and formed seven interactions involving multiple interaction types. Hydrogen bonds with Arg179 (4.23 Å) and Ser254 (4.45 Å) were complemented by a π-alkyl interaction with Arg179 (4.17 Å), a π-anion interaction with Arg179 (5.40 Å), and a π-cation interaction with Glu252 (6.60 Å). Furthermore, a π-σ interaction with Ser182 (4.05 Å) and an amide-π stacking interaction with Phe181 (7.63 Å) contributed to stabilization. Although the aromatic stacking distance is slightly longer than ideal, the combination of hydrophobic and electrostatic contacts resulted in a favorable docking score.

Derivative D5

D5 exhibited a binding energy of -6.3 kcal mol^-1 and formed a total of four interactions with the amino acid residues of VP30. Two hydrogen bonds were noticed with Lys180 (4.66 Å) and Arg179 (6.29 Å). Moreover, Arg179 also participated in two hydrophobic interactions: a π-alkyl (6.64 Å) and an alkyl (4.62 Å). However, the relatively long interaction distances, particularly the hydrogen bond at 6.29 Å and π-alkyl interaction at 6.64 Å, reduced their stabilizing contributions.

Derivative D6

D6 exhibited a binding energy of -6.1 kcal mol^-1 and was involved in five intermolecular interactions with key active-site residues Arg179, Phe181, and Ser182. Specifically, the ligand formed an alkyl interaction with Arg179 (4.19 Å), a relatively weak hydrogen bond with Arg179 (4.74 Å), an amide-π stacking interaction with Phe181 (8.63 Å), a π-σ interaction with Ser182 (4.22 Å), and a C-H bond interaction with Ser182 (4.26 Å). While the hydrophobic interactions occurred within favorable distances (4.19-4.26 Å), the amide-π stacking interaction with Phe181 occurred at 8.63 Å, indicating relatively weak aromatic stabilization.

Derivative D7

D7 demonstrated a binding energy of -6.6 kcal mol ¹ and formed six intermolecular interactions with key residues of VP30. These included a hydrogen bond with Arg179 (7.24 Å), π-σ interactions with Ser182 (4.26 Å) and Arg179 (4.29 Å), a C-H bond interaction with Ser182 (4.12 Å), an amide-π stacking interaction with Phe181 (8.61 Å), and a π-cation interaction with Arg179 (5.40 Å). The π-σ and C-H interactions occurred within favorable distances and contributed to ligand stabilization, whereas the hydrogen bond and amide-π stacking interaction were considerably longer than their optimal ranges, indicating weaker contributions to binding.

Derivative D8

D8 demonstrated a binding energy of -6.3 kcal mol^-1 and established six intermolecular interactions with Lys183, Ser184, Gly219, Ser254, and Val256. D8 formed a π-alkyl interaction with Lys183 (5.59 Å), a hydrogen bond with Ser184 (4.13 Å), a π-donor hydrogen bond with Ser254 (4.51 Å), a C-H bond interaction with Ser254 (5.30 Å), an alkyl interaction with Val256 (4.31 Å), and a C-H bond interaction with Gly219 (3.29 Å). Of these contacts, the C-H bond with Gly219 and the alkyl interaction with Val256 are within acceptable contact ranges, and likely make the largest contributions to ligand stability. The hydrogen bond with Ser184 and the π-donor hydrogen bond with Ser254 are too long for optimum hydrogen bonding, so these are likely weak contributors to binding. In addition, the π-alkyl interaction with Lys183 and the C-H bond with Ser254 are fairly long, so these likely contribute relatively little stabilization.

Derivative D9

D9 demonstrated a binding energy of -6.8 kcal mol^-1 and established seven intermolecular interactions, such as hydrogen bonds, C-H bonding, π-alkyl, and π-cation interactions, which support complex stability. The strongest interaction is a classical hydrogen bond with Ser184 at a distance of 2.71 Å, which is in the range of moderate hydrogen bonds and should contribute noticeably to the stability of the binding. Another hydrogen bond with Ser184 (4.23 Å) also exists; this bond is weaker than the first bond but still has some stabilizing contribution. Two more hydrogen bonds are found with Gln212 (5.55 and 5.66 Å); nevertheless, these are at long distances and therefore are likely to have a small contribution to the complex stabilization. In addition, a C-H bond interaction with Lys183 (5.72 Å) adds stabilization, although the energetic contribution is limited by the extended distance. Also, a π-alkyl interaction is observed between the ligand and Lys218 (4.94 Å), which is an acceptable distance for hydrophobic contacts and may contribute to the retention of the ligand in the binding pocket. Additionally, a π-cation interaction with Arg179 (6.21 Å) is observed; although π-cation interactions are generally favorable, this relatively long distance suggests a weaker electrostatic contribution.

Derivative D10

D10 exhibited a binding energy of -6.4 kcal mol ¹ and formed six interactions with His215, Phe181, Glu252, Lys218, and Arg179. The ligand established a C-H bond with His215 (5.93 Å), an amide-π stacking contact with Phe181 (8.57 Å), and a π-anion contact with Glu252 (7.05 Å). Furthermore, satisfactory hydrophobic contacts were observed through an alkyl interaction with Lys218 (3.94 Å) and a π-alkyl interaction with Arg179 (4.23 Å), which play a role in ligand stabilization inside the binding pocket. Nevertheless, an unfavorable donor-donor interaction with Arg179 (5.19 Å) causes electrostatic repulsion that may decrease binding efficiency. Overall, D10 binding is mainly governed by hydrophobic and π-mediated interactions, whereas the absence of conventional hydrogen bonds and the existence of the unfavorable donor-donor contact possibly restrict its total binding strength.

Derivative D11

D11 exhibited a binding energy of -6.2 kcal mol^-1 and established four interactions with the target protein. The ligand established an H-bond with His125 (5.55 Å) and a C-H bond with Asp231 (5.28 Å); nevertheless, these distances are substantially longer than those characteristically related to strong hydrogen bonds, signifying comparatively weak stabilizing effects. Other interactions comprise a π-alkyl contact with Ala119 (4.50 Å), which may contribute to a reasonable hydrophobic stabilization, and a π-donor hydrogen bond with Gln233 (6.55 Å), which could play a role in orienting the ligand within the binding pocket. Overall, the primary modes of binding of D11 are through weak hydrogen bonding and hydrophobic interactions. The large interaction distances and lack of strong conventional hydrogen bonds result in a less stable protein-ligand complex.

Derivative D12

D12 exhibited a binding energy of -5.9 kcal mol^-1 and formed only three hydrogen bonds with Ser234 (4.74 Å), Arg123 (5.83 Å), and Arg196 (4 Å). Of all these interactions, hydrogen bonding with Arg196 is the shortest and thus is expected to make the major contribution to stabilizing the complex, whereas the longer hydrogen bonds with Ser234 and especially Arg123 are expected to be weaker and mostly electrostatic. The results further indicate that hydrogen bonding alone is insufficient to achieve high affinity of the ligand for the target protein in the absence of complementary hydrophobic interactions.

Derivative D13

D13 demonstrated a binding energy of -5.3 kcal mol^-1 and established five interactions with the target protein. It formed two hydrogen-bonding interactions with Gln233 (5.67 and 5.77 Å). The two interactions are well outside the distances expected for strong hydrogen bonds, so only a weak electrostatic stabilization effect can be predicted from these contacts. However, other favorable interactions are formed by D13 with Tyr122 and Met233 via two alkyl interactions (at 4.41 and 5.20 Å, respectively), and Ala119 via a π-alkyl interaction (at 4.36 Å). These interactions with Tyr122 and Ala119 are at distances which might be considered reasonable for efficient hydrophobic contacts and probably play a role as main stabilizing interactions. However, the lack of any well-formed hydrogen bonds and few strong hydrophobic interactions leads to only a weak contribution to the stability of the protein-ligand complex and consequently to a comparatively low binding energy for the derivative, suggesting that D13 has a low affinity for the protein.

Derivative D14

D14 has a binding energy of -5.8 kcal mol^-1 and forms five interactions within the binding pocket of the target protein. It forms a π-σ interaction with Ala119 (4.25 Å), resulting in favorable hydrophobic stabilization of the complex. There is also a hydrogen bond formed with Gln233 (3.87 Å). This bond represents the strongest interaction for D14 and likely contributes most significantly to stabilization of the complex. Additional interactions include hydrogen bonds with Arg123 (6.60 Å), a π-donating hydrogen bond with Tyr122 (5.21 Å), and a C-H bond with Pro120 (5.78 Å). However, since these interactions occur over relatively long distances, the contributions are expected to be weak and result mainly from electrostatic and orientational forces, rather than strong binding interactions and stabilization.

Derivative D15

D15 exhibited a binding energy of -5.1 kcal mol^-1 and established three interactions within the active site of the target protein. D15 established a π-anion interaction with Asp231 (5.86 Å) and both π-donor hydrogen-bonding (5.61 Å) and conventional hydrogen-bonding (6.14 Å) interactions with Gln233. The low affinity of D15 may be attributed to the relatively small number of protein-ligand interactions and their longer interaction distances. These factors weaken the overall stabilizing forces between the D15 and the active-site residues, leading to reduced complex stability.

Derivative D16

D16 binds within the active-site cavity with a binding energy of -5.8 kcal mol^-1, making two hydrogen bonding interactions with Lys180 (4.07 and 4.76 Å), an amide-π stacking interaction with Gly219 (4.10 Å), a π-cation interaction with Lys218 (4.93 Å), and π-alkyl interactions with Val257 (6.12 Å) and Lys183 (6.74 Å). Results propose that D16 is accommodated inside the binding pocket via a combination of electrostatic and hydrophobic contacts; nevertheless, the unsatisfactory interaction distances may be responsible for its moderate affinity towards the target protein.

Derivative D17

D17 showed a binding energy of -5.7 kcal mol^-1 and made six contacts inside the binding pocket of the target protein. It established hydrogen bonds with Asp217 (4.03 Å) and Arg179 (3.40 and 4.35 Å), plus π-alkyl interactions with Arg179 (4.38 Å) and Lys218 (5.55 Å). The hydrogen bond linking D17 and Arg179 acts as the best interaction here and probably helps to stabilize the ligand quite well. Still, an unfavorable acceptor-acceptor interaction with Glu252 (5.52 Å) was observed too, which might cause electrostatic repulsion and lower total binding efficiency. These results show that D17 stays stable mostly via connections to Arg179, but the unfavorable contact and some longer-distance interactions restrict its binding strength.

The overall molecular docking results of derivatives (D1 D17) are summarized in Table 2, while comprehensive protein-ligand interaction profiles of all docked quercetin derivatives with EBOV VP30 are provided in Figure S1 (SI section).

Table 2
Binding energies, number and nature of interactions, types of amino acid residues, and bond distances for the molecular docking of quercetin derivatives (D1-D17) with the VP30 protein

Quercetin derivative	Binding energy / (kcal mol-1)	Number of interactions	Residue/nature of interaction	Distance of interaction / Å
D1	-9.1	11	Arg179/H-bonding Arg179/H-bonding Arg179/H-bonding Arg179/π-σ Lys183/π-alkyl His215/π-σ His215/π-π T-shaped Ser184/H-bonding Ser184/H-bonding Glu252/H-bonding Lys218/H-bonding	3.16 3.18 4.42 5.38 6.25 6.22 5.96 3.15 4.23 4.92 5.60
D2	-8.7	9	Lys183/π-alkyl Lys218/unfavorable positive-positive Lys218/π-alkyl Gln212/H-bonding Gln233/H-bonding Ser187/H-bonding Ser187/H-bonding Ser254/H-bonding Ser182/C-H bonding	4.91 4.61 4.85 2.71 2.75 4.88 4.08 4.66 4.54
D3	-7.8	6	Ser187/H-bonding Lys218/alkyl Glu252/H-bonding Arg179/H-bonding Arg179/H-bonding Arg179/π-alkyl	4.86 3.90 5.05 4.42 5.84 5.70
D4	-7.5	7	Phe181/amide-π stacking Ser182/π-σ Arg179/π-alkyl Arg179/H-bonding Arg179/π-anion Glu252/π-cation Ser254/H-bonding	7.63 4.05 4.17 4.23 5.40 6.60 4.45
D5	-6.3	4	Lys180/H-bonding Arg179/H-bonding Arg179/π-alkyl Arg179/alkyl	4.66 6.29 6.64 4.62
D6	-6.1	5	Arg179/alkyl Arg179/H-bonding Phe181/amide-π stacking Ser182/π-σ Ser182/C-H bonding	4.19 4.74 8.63 4.22 4.26
D7	-6.6	6	Arg179/H-bonding Ser182/π-σ Ser182/C-H bonding Phe181/amide-π stacking Arg179/π-σ Arg179/π-cation	7.24 4.26 4.12 8.61 4.29 5.40
D8	-6.3	6	Lys183/π-alkyl Ser184/H-bonding Ser254/π-donor H-bonding Ser254/C-H bonding Val256/alkyl Gly219/C-H bonding	5.59 4.13 4.51 5.30 4.31 3.29
D9	-6.8	7	Arg179/π-cation Gln212/H-bonding Gln212/H-bonding Lys183/C-H bonding Lys218/π-alkyl Ser184/H-bonding Ser184/H-bonding	6.21 5.66 5.55 5.72 4.94 2.71 4.23
D10	-6.4	6	His215/C-H bonding Phe181/amide-π stacking Glu252/π-anion Lys218/alkyl Arg179/π-alkyl Arg179/ unfavorable donor-donor	5.93 8.57 7.05 3.94 4.23 5.19
D11	-6.2	4	His125/H-bonding Asp231/C-H bonding Ala119/π-alkyl Gln233/π-donor hydrogen bonding	5.55 5.28 4.50 6.55
D12	-5.9	3	Ser234/H-bonding Arg123/H-bonding Arg196/H-bonding	4.74 5.83 4.00
D13	-5.3	5	Gln233/H-bonding Gln233/H-bonding Tyr122/alkyl Met233/alkyl Ala119/π-alkyl	5.67 5.77 4.41 5.20 4.36
D14	-5.8	5	Ala119/π-σ Arg123/H-bonding Tyr122/π-donor hydrogen bonding Gln233/H-bonding Pro120/C-H bonding	4.25 6.60 5.21 3.87 5.78
D15	-5.1	3	Asp231/π-anion Gln233/π-donor hydrogen bonding Gln233/H-bonding	5.86 5.61 6.14
D16	-5.8	6	Lys180/H-bonding Lys180/H-bonding Gly219/amide-π stacking Val257/π-alkyl Lys218/π-cation Lys183/π-alkyl	4.07 4.76 4.10 6.12 4.93 6.74
D17	-5.7	6	Asp217/H-bonding Glu252/unfavorable acceptor-acceptor Arg179/H-bonding Arg179/H-bonding Arg179/π-alkyl Lys218/π-alkyl	4.03 5.52 4.35 3.40 4.38 5.55

Altogether, derivatives D1 and D2 showed the most favorable binding profile among the seventeen quercetin derivatives examined against the VP30 protein, with the least binding energies and the highest number of bonding interactions with some key residues of VP30, such as Arg179, Lys183, Gln212, Lys218, Gly219, Ser184, Ser254, Glu252, Ala119, His215, and Val256. These results point toward stable ligand-protein complexes and highlight D1 and D2 as the most prospective derivatives in this series. D15, on the other hand, showed the lowest affinity as seen by its highest binding energy (-5.1 kcal mol ¹) with only three bonding interactions, suggesting its little stabilization inside the binding pocket of VP30 protein. The docking results show a distinct variation in binding efficacy among the derivatives; D1 turns out to be the most powerful candidate and D15 the least effective in terms of anticipated interaction with the VP30 protein.

D1 and D2 were selected for further computational studies because of their superior docking results. The dynamic stabilities of the VP30-ligand complexes were studied using MD simulations, and their pharmacokinetic suitability and electronic properties were evaluated by performing ADMET and DFT analyses, respectively. This stepwise computational workflow allowed for the detailed characterization of the most promising candidates identified in the initial docking screen.

MD simulation results

Structural parameters were calculated from the average values of the trajectories collected from the stabilized regions of the MD simulations, thus providing a reliable measure of the stability and dynamic characteristics of the VP30-ligand complexes.

The docked complexes of the two highest-ranked derivatives (D1 and D2) were selected for MD simulation. The MM/PBSA was determined for the VP30-flavonoid complexes. Table 3 shows that the MM/PBSA-computed binding free energy values for the COM-1122 and COM 1133 complexes are -22.1562 and -19.7180 kcal mol ¹, respectively. These values predict good hydrophobic and H-bonding interactions between the ligand and the protein in the complex systems, reflecting the binding interactions of the ligand within the dynamic site of the protein.

Table 3
MM/PBSA values for complexes COM-1122 and COM-1133

Complex	Electrostatic energy / (kcal mol-1)	Van der Waals energy / (kcal mol-1)	Non-polar energy / (kcal mol-1)	Polar energy / (kcal mol-1)	∆Gbind / (kcal mol-1)
COM-1122	-2.7788	-22.9897	-2.1656	5.7779	-22.1562
COM-1133	-1.0156	-21.3869	-0.1326	2.5519	-19.7180

RMSD results

RMSD calculation was conducted for the purpose of analyzing the stability of the apo VP30 protein structure and the ligand complex structures throughout the MD simulation process lasting 100 ns (Figure 4). The RMSD value of the apo VP30 protein (black line) displayed fluctuation at the beginning of the MD simulation process and attained the level of 4 Å between 0.5 and 0.7 ns, which suggested that the apo structure had higher flexibility without ligand binding. After equilibration, the apo system remained relatively stable with moderate conformational variations. The two complexes with bound ligands, COM-1122 and COM-1133, appeared to maintain stable trajectories after ca. 20 ns into the simulation. COM-1122 had large fluctuations up to 4.4 Å from 10 to 25 ns and then became relatively stable for the rest of the simulations, except for a small fluctuation at 89 ns. Conversely, COM 1133 showed smoother fluctuations in terms of RMSD with fewer fluctuations during the entire simulation process. Despite fluctuations that occurred in the latter part of the simulations, these fluctuations probably represent the conformational adaptations of the ligand receptor complex towards a better fit. Generally, it can be seen from the RMSD graphs that ligand binding played a role in stabilizing the VP30 structure through the reduction of its excessive flexibility when compared with the apo state. Stable trajectories obtained along with good MM/PBSA and Rg values indicate stability of the complexes formed during the simulation.

Figure 4
Root mean square deviation (RMSD) graph of Apo form and the complexes COM-1122 and COM-1133 during MD simulation.

RMSF results

Figure 5 shows that the RMSF value for the Apo form was 4.6 Å for residues 100 to 105, with maximum fluctuations observed between 130 and 250 Å. In contrast, COM-1122 exhibited an RMSD value of 3 Å for residue 230. Furthermore, large variations of up to 2.8 Å were noticed for residues 225 to 230. On the other hand, COM 1133 showed somewhat more discrepancy, with a maximum fluctuation of 3.5 Å noted for residue 220. Overall, the complexes exhibited obvious flexibility, probably due to conformational discrepancies. These results recommend that the conformational variations are owing to the binding of ligands to the VP30 protein.

Figure 5
Root mean square fluctuation (RMSF) profiles of the Apo state and the complexes COM-1122 and COM-1133 during simulation.

R_g results

The R_g is a crucial indicator for the assessment of protein stability, as it measures the firmness of the protein in both Apo and complex shapes. Figure 6 demonstrates the R_g data, signifying the folding and unfolding of the protein throughout simulation. The Apo form shows significant fluctuations between 70 and 85 ns but stabilizes from 85 to 100 ns, with an average R_g value of 18.2 Å. COM 1122 exhibits fluctuations up to 18.6 Å at 79.9 ns, with an average R_g value of 18.3 Å. In contrast, COM 1133 does not exhibit any significant fluctuations and remains stable between 25 and 59 ns, with an average R_g value of 18 Å. Overall, the relatively stable R_g values of the complexes signify firm folding states, which can be attributed to the ligand-protein binding.

Figure 6
Plot of radius of gyration (R_g) values of Apo form and the complexes COM-1122 and COM-1133 during MD simulation.

PCA results

PCA is a mathematical procedure based on covariance matrices used to analyse atomic displacements and loop dynamics in proteins.⁵¹ PCA was employed to detect conformational changes in the complex systems during a 100 ns simulation. During this period, shifts in motion were observed in the two major elements of the systems.

The PCA of the two principal components for the Apo form and the complex COM-1122 is represented by scattered and clustered dots. Figure 7 displays that the Apo form exhibits discrete movement for both PC1 and PC2, as indicated by the color-coded points. In contrast, the complex COM-1122 demonstrates more clustered and compact behavior, with PC1 and PC2 ranging from -80 to 50. PC1 is slightly dispersed between 20 and 40. The shaded area in the figure highlights the drug-binding region, which ranges from -50 to 60. The complex is compact and exhibits minimal dispersion, indicating stability and strong ligand-protein interactions.

Figure 7
Principal component analysis (PCA) of two protein systems: (a) Apo form and (b) COM-1122 form. Blue and red colours represent different motions in the PCA, highlighting the correlation between the dispersed motion of the Apo form and the compact motion of the complex form.

Overall, it can be concluded from the MD simulations that complexes with maximum interactions showed increased stability due to drug binding. The simulation data endorsed the docking results, with stable RMSD values. A comparative study of the RMSD values points out the flexibility and differences in the dynamic behavior of VP30 throughout the simulation, thus identifying the stability of the complexes. This result reinforced the docking data, as the complexes exhibited appropriate docking energies and binding free energies. RMSF analysis was employed to estimate the flexibility of the Apo state and the complexes, demonstrating conformational variations attributable to ligand binding to VP30. The Apo state exhibited the highest variations, whereas both complexes showed prominent variations as well.

The complexes showed comparatively stable structural behavior during the entire course of the simulation and also exhibited balance in their R_g values. Mainly, the moderately stable R_g values of the complexes specified steady folding forms, possibly emerging from the ligand binding. The designated inhibitors display probable efficacy for in vitro investigations versus the EBOV VP30.

Plants are invaluable sources of natural products that frequently offer greater efficiency with negligible side effects to cure many disorders. Docking investigations indicated that flavonoids established several non-bonded interactions within the active spots of EBOV VP30 proteins and showed suitable binding free energies. The conformational changes were confirmed through MD simulation data, together with MM/PBSA, RMSD, H-bonding, RMSF, P.E., R_g, and PCA, which revealed the resilience and strength of the VP30-flavonoid complexes. The results demonstrated the ligands binding behavior in the protein active sites and the protein-complex interactions during simulation. The MD simulation results supported the molecular docking study findings.

ADMET results

Due to their excellent binding energies, stability in molecular dynamics, and positive MM/PBSA binding energies, quercetin derivatives D1 and D2 were chosen for further ADMET studies to determine whether or not they have suitable properties as potential EBOV VP30 inhibitors.

ADMET predictions for quercetin derivatives D1 and D2 give useful information regarding their efficacy as potential drug candidates against the EBOV VP30 protein. This is because the VP30 protein is crucial in activating transcription of genes in the EBOV,⁷² therefore making it necessary to have drugs capable of effectively binding to the protein while having good pharmacokinetics and minimal toxicity. In this regard, the ADMET properties suggest that the two quercetin derivatives D1 and D2 have ideal safety and metabolism profiles, but certain limitations exist with respect to oral bioavailability, polarity, and toxicity.

Physicochemical characteristics of D1 and D2 make it evident that both molecules exhibit structural characteristics that could potentially render them useful as antiviral agents based on their flavonoid core. Molecular weight values of 450.35 Da for D1 and 464.38 Da for D2 fall into the acceptable range for small molecules as antivirals, although being quite close to the upper limit imposed by Lipinski’s rule of five. Medium fraction of sp³ hybridized carbons (0.25 and 0.29), and a low number of rotatable bonds (3 and 4) indicate relatively stiff molecules, which can be advantageous in terms of stabilizing ligand-protein interactions essential for blocking of viral transcription progressions.

Both derivatives possessed high hydrogen bonding ability both as donors and acceptors (HBD = 8; HBA = 12), together with an increased polar surface area (TPSA = 210.51 Ų). The above properties indicate the potential of the two derivatives to form hydrogen bonds with amino acids of the VP30 functional or regulatory sites, leading to complex formation and potential antiviral effects. However, the increase in polarity and violation of Lipinski’s rules could simultaneously decrease passive cell membrane permeability and bioavailability via the oral route. Nevertheless, numerous antiviral flavonoids display biological activity via extensive hydrogen bonding and electrostatic interactions with viral proteins,⁷³,⁷⁴ thus allowing the above-mentioned physicochemical properties to be suitable for developing antiviral drugs.

The results of lipophilicity are consistent with the above finding. The values of XLOGP3 were low (-0.09 and 0.36 for D1 and D2, respectively) for both derivatives, while the consensus logP values were also negative (-0.49 and -0.27 for D1 and D2, respectively), indicating hydrophilicity. The low to moderate lipophilicity would help to enhance the aqueous solubility while minimizing the potential toxic effects, although it might affect membrane penetration. The moderately predicted aqueous solubility of derivatives D1 and D2 (Log S values of 2.75 and 3.04, respectively) is beneficial because it ensures easier dissolution and transport within biological systems, which is necessary for antiviral drugs that require systemic distribution.

The pharmacokinetic predictions showed that both derivatives had moderate intestinal absorption, but D2 (44.93%) had better absorption than D1 (34.15%). Such a feature can be explained by the higher lipophilic properties and higher molecular flexibility of derivative D2. However, the very low values of BBB permeation (logBB < -2) and low brain/central nervous system (CNS) permeability (logPS ca. -4.7) confirm that neither derivative will be able to pass through the BBB. When developing antiviral compounds for EBOV VP30, having poor CNS penetration is actually advantageous when taking antivirals systemically, since this will reduce chances of any neurological side effect and increase the potency of the drug in peripheral tissues where the virus replicates. Also, it is predicted that there is no evidence of each derivative being a P-glycoprotein substrate, which would lead to decreased active efflux and give a better retention within infected cells.

One of the key benefits for both derivatives from a pharmacological perspective is the absence of inhibition against any of the important cytochrome P450 enzymes (CYP1A2, CYP2C19, CYP2C9, CYP2D6, and CYP3A4). This would be a very positive trait for potential antiviral drugs, since it significantly decreases the possibility of any metabolic interactions between drugs, especially in cases where a combination therapy of antivirals is required for the treatment of serious viral illnesses.

The estimated bioavailability score of 0.17 for both molecules implies a low level of oral bioavailability because of the high polarity and strong hydrogen-bonding ability of both molecules. However, the limitation could possibly be circumvented using strategies like nano-encapsulation and prodrug formation, methods employed for polyphenolic antiviral drugs.⁷⁵,⁷⁶ Further, both molecules produced a single Pan-Assay INterference compoundS (PAINS) hit for the presence of the benzen-1,2-diol group (catechol_A). The presence of such a group is quite common among flavonoid molecules and could cause interference in the assay because of its potential redox reactions or non-specific binding properties. Although this alert does not nullify the derivatives, it suggests that the experimental justification should be judiciously interpreted.

Clearance predictions indicated that D2 is predicted to have marginally higher total clearance (0.626 log mL min ¹ kg ¹) compared to D1 (0.504 log mL min^-1 kg^-1), indicating a slightly faster rate of elimination. This will probably give D2 lower overall systemic accumulation. No evidence of renal OCT2 substrate properties also minimizes the chances of renal drug interaction.

Toxicological predictions revealed an overall good safety profile. Both derivatives presented predictions in the high range of LD₅₀ values (5000 mg kg^-1), and they are both categorized as class 5 toxicity, which corresponds to a pretty low general toxicity. Besides hepatotoxicity, neurotoxicity, cardiotoxicity, cytotoxicity, and mutagenicity appear inactive, endorsing their potential safety as drug candidates. Such traits are obligatory for antiviral drugs that may require repeated dosing during viral infection treatment.

However, there are certain toxicity issues that should be considered. The likelihood that both compounds possess nephrotoxicity is relatively high. However, the probability of nephrotoxicity in compound D2 is higher (0.76) compared to that of compound D1 (0.68). Compound D1 is likely to cause the inhibition of the hERG II channel, which could lead to electrophysiological disorders of the heart, while compound D2 does not have such a tendency. Clinical toxicity can become an issue in the case of using compound D2 because it poses some level of threat (probability is 0.51). Nutritional toxicity was also predicted for both derivatives, potentially associated with the interference with nutrient metabolism or gastrointestinal tolerance. The overall ADMET results of derivatives D1 and D2 are summarized in Table 4.

Table 4
In silico prediction of ADMET parameters for quercetin derivatives D1 and D2

entry	Property	Quercetin derivative D1	Quercetin derivative D2
1	MW / Da	450.35	464.38
2	fraction Csp3	0.25	0.29
3	number of rotatable bonds	3	4
4	HBA	12	12
5	HBD	8	8
6	TPSA / Å2	210.51	210.51
7	XLOGP3 (Log Po/w)	-0.09	0.36
8	consensus (Log Po/w)	-0.49	-0.27
9	Log S (calculated with the ESOL model)	-2.75	-3.04
10	intestinal absorption in humans (absorbed) / %	34.151	44.93
11	BBB permeability (Log BB)	-2.372	-2.243
12	CNS permeability (Log PS)	-4.714	-4.752
13	P-gp substrate	no	no
14	CYP1A2 inhibitor	no	no
15	CYP2C19 inhibitor	no	no
16	CYP2C9 inhibitor	no	no
17	CYP2D6 inhibitor	no	no
18	CYP3A4 inhibitor	no	no
19	Lipinski’s rules (violations)	2 violations: HBA > 10, HBD > 5	2 violations: HBA > 10, HBD > 5
20	bioavailability score	0.17	0.17
21	PAINS # alerts	1 alert: catechol_A	1 alert: catechol_A
22	predicted median lethal dose (LD50) / (mg kg-1)	5000	5000
23	max. tolerated dose in humans (log (mg kg-1 day-1))	0.809	0.913
24	hERG I inhibitor	no	no
25	hERG II inhibitor	yes	no
26	total clearance (log (mL min-1 kg-1))	0.504	0.626
27	renal oct2 substrate	no	no
28	predicted toxicity class	5	5
29	hepatotoxicity	inactive	inactive
30	neurotoxicity	inactive	inactive
31	nephrotoxicity	active (probability = 0.68)	active (probability = 0.76)
32	cardiotoxicity	inactive	inactive
33	cytotoxicity	Inactive	Inactive
34	mutagenicity	inactive	inactive
35	clinical toxicity	inactive	active (probability = 0.51)
36	nutritional toxicity	active (probability = 0.64)	active (probability = 0.55)
37	skin sensitization	no	no

Interference compounds

The bioavailability radar charts of quercetin derivatives D1 and D2 are given in Figure 8 for an instant assessment of drug-likeness. The pink area designates the optimum range for properties including lipophilicity (XLOGP3: -0.7 to +5.0), molecular weight (150-500 Da), topological polar surface area (TPSA: 0-140 Ų), aqueous solubility coefficient (logS ≤ 6), fraction of sp³ carbons (C sp³ ≥ 0.25), and a maximum of nine rotatable bonds.

Figure 8
The bioavailability radar of quercetin derivatives (a) D1 and (b) D2 via the SwissADME predictor. The two compounds remain mostly within the favorable physicochemical window, indicating balanced properties. LIPO: lipophilicity, measured by XLOGP3; SIZE: molecular weight; POLAR: polarity, indicated by TPSA (topological polar surface Area); INSOLU: insolubility in water, assessed by the logS parameter; INSATU: Instauration, the fraction of sp3 hybridized carbon atoms; and FLEX: flexibility, quantified by the number of rotatable bonds.

Figure 9 explains the Boiled Egg plots for quercetin derivatives D1 and D2. The white part shows where we expect good human intestinal absorption, and the yellow part highlights possible CNS penetration. Drug candidates estimated to be absorbed via non-oral routes are located in the gray region of the plot.⁷⁷ The analysis shows that D1 and D2 have limited gastrointestinal absorption and negligible BBB penetration; however, these characteristics do not detract from their potential use as anti-Ebola drugs. Since the action of VP30 occurs in peripheral tissues instead of CNS tissue, the lack of BBB permeability may represent a benefit in reducing neurological side effects of D1 and D2. Additionally, the moderate predicted intestinal absorption and strong affinity for the VP30 target make D1 and D2 attractive lead candidates; their limited permeability is primarily due to their high polarity and significant ability to form hydrogen bonds, both of which create strong interactions with the VP30 binding site. Consequently, D1 and D2 should be regarded as potential lead structures for further optimization through medicinal chemistry and enhanced formulation strategies to increase bioavailability.

Figure 9
The Boiled Egg model predicted for quercetin derivatives (a) D1 and (b) D2. The models depict the compounds GI absorption and Blood-Brain Barrie (BBB) penetration. The white part of the egg shows compounds that are possibly to be passively absorbed in the GI tract, whereas the yellow (yolk) part indicates BBB permeability. The compounds positioned outside the egg are predicted to be poorly absorbed and scarcely penetrate the brain.

The toxicity radar charts (Figure 10) quickly demonstrate a graphical comparison of the toxicity profiles of quercetin derivatives D1 and D2, focusing on the confidence of positive toxicity results for both in comparison to the mean of their respective classes of compounds.

Figure 10
The toxicity radar charts of quercetin derivatives (a) D1 and (b) D2, which are proposed to quickly illustrate the confidence of positive toxicity results compared to the average of their class.

Figure 11
The network charts of quercetin derivative (a) D1 and (b) D2 are proposed to swiftly demonstrate the connection between their molecular characteristics and anticipated toxicity classes.

Overall, the ADMET predictions show that quercetin derivatives D1 and D2 have several features that favor their development as effective EBOV VP30 inhibitors. The two derivatives have good metabolic interaction patterns, a non-toxic acute profile, moderate absorption properties, and strong intermolecular interactions with the viral proteins. The derivative D2 is pharmacokinetically superior compared to D1 due to better absorption in the gut and no hERG II inhibition. However, the predicted nephrotoxicity, low oral bioavailability, and clinical toxicity of the two derivatives need thorough evaluation. The whole ADMET results indicate that the two derivatives have the potential to act as anti-EBOV VP30 antiviral molecules and warrant further justification through in vitro antiviral analyses and in vivo pharmacological investigations.

DFT results

The frontier molecular orbitals (FMO), the HOMO, and the LUMO, are the simplest yet essential characteristics of the effective electron donors and acceptors, respectively. The HOMO is where the electrons are more likely to interact with the chemical environment as an effective donor; the LUMO is the orbital into which the electrons should be transferred. The HOMO-LUMO gap (∆E) determines the effective charge transfer into the molecule and is also strongly related to the optical polarizability, chemical reactivity, and stability.⁷⁸,⁷⁹

The computed values of FMO energies and global reactivity indicators offered important information concerning the stability, electron transfer capacity, and binding tendency of the quercetin derivatives D1 and D2 with the viral target. Compound D1 demonstrated HOMO and LUMO energies of -3.1663 and -2.3276 eV, respectively, while D2 possessed slightly lower values of -3.5769 and -2.6895 eV (Figure 12). The lower value for HOMO energy in D2 implies increased electronic stability and resistance to electron donation as compared to D1. In contrast, the higher HOMO energy level of D1 suggests that it possesses better donating properties, allowing it to bond well with the electron-deficient residues of the amino acids present in the VP30 binding site. The HOMO-LUMO energy gap (∆E) determines the stability and reactivity of molecules. D1 showed a lower energy gap (0.8387 eV) compared to D2 (0.8874 eV) (Figure 12). The computed HOMO-LUMO energy gaps for D1 and D2 correspond to a low energy gap, which is typically related to high molecular activity, electronic polarizability, and charge-transfer properties. The presence of relatively low energy gaps would allow for a greater probability of electronic rearrangement during molecular interactions with biological macromolecules, resulting in greater molecular recognition and intermolecular binding actions. Accordingly, the small ∆E values calculated for D1 and D2 might aid in the promising interaction with the VP30 protein as observed in the molecular docking study results.⁸⁰

Figure 12
HOMO-LUMO energy gap for (a) derivatives D1 and (b) D2 computed via DFT technique at the B3LYP/6-311G (d,p) level.

Subsequently, DFT analysis was also employed to calculate the molecular descriptors of the chemical reactivity of compounds D1 and D2. These descriptors include the chemical hardness (η), softness (S), potential (μ), electrophilicity index (ω), electronegativity (χ), and dipole moment (μ). Among these indices, hardness is the most crucial one, which is an index of the chemical stability against fluctuations in electron distribution and charge transfer. Based on FMO theory, the hardness of a molecule is determined as the energy difference between the HOMO and LUMO, and is quantitatively expressed as η = (E_LUMO - E_HOMO)/2. Whereas, the softness of a molecule can be easily defined by the inverse value of hardness, i.e., S = 1/η. A small energy gap between them means low hardness, low stability, and higher chemical reactivity of a compound.⁸¹ The electronegativity (χ) is represented by the χ = - (E_HOMO + E_LUMO)/2 equation. It is a measure of an atom/molecule tendency to attract electrons. The negative value of electronegativity is referred to as chemical potential (μ) and quantified using the formula μ = (E_HOMO + E_LUMO)/2. The electrophilicity index (ω) is Parr’s concept, which reflects the electron-attracting ability of a molecule and is calculated as ω = μ²/2η.⁸²

The computed chemical hardness (η) and softness (S) also reinforce this hypothesis. D1 had a smaller hardness (0.4194 eV) but larger softness (1.1920 eV^-1) than D2, which had hardness and softness equal to 0.4437 eV and 1.1269 eV^-1, respectively. The softer compounds tend to be more polarizable and chemically flexible, allowing them to form favorable non-covalent interactions like hydrogen bonding, π-π interactions, and electrostatic interactions with the VP30 amino acids. Consequently, the higher softness of D1 may imply a higher level of electronic flexibility in interacting with VP30.

Also, electronegativity (χ) and chemical potential (µ) are significant parameters for controlling electron transfer in ligand-protein interactions. D2 possessed higher electronegativity (3.1332 eV) and a more negative chemical potential (-3.1332 eV) compared with D1 (2.7470 and -2.7470 eV, respectively). The results suggest that D2 has a higher affinity for electrons and can act as a better electron acceptor in the interaction with VP30 residues. Such a behavior might contribute to the stabilization of ligand-protein complexes through charge-transfer interactions.

Electrophilicity index (ω) gives further insight into the ability of these compounds to act as electron acceptors. Compound D2 has shown higher electrophilicity (11.0679 eV) compared to that of compound D1 (8.996 eV). This clearly indicates that the compound D2 has much better potential to accept the electron density provided by the surrounding protein. Increased electrophilicity may contribute to enhanced intermolecular interactions, which can then ultimately lead to enhanced biological activity.

A significant difference between the two derivatives was noted in the dipole moment values. It was observed that the dipole moment of D2 was significantly higher (6.4998 D) than that of D1 (0.7184 D). Due to higher polarity in D2, it can be anticipated that its intermolecular interactions will be better facilitated via hydrogen bonding as well as electrostatic interactions with polar amino acids present in VP30.

Overall, the DFT results show that these quercetin derivatives exhibit good electronic parameters conducive to good antiviral activity for EBOV VP30. The lower HOMO-LUMO gap and higher level of chemical softness of D1 could account for the stronger docking affinities and stable MD behavior by enabling effective charge transfer interactions inside the VP30 binding site. In contrast, D2 shows high electrophilicity, electronegativity, and dipole moment value indicative of strong electron-accepting ability and stronger intermolecular interactions with the amino acid residue of VP30. Taken together, each of these compounds holds the potential to be a promising VP30 inhibitor, with D1 exhibiting good reactivity and D2 showing adherence to strong intermolecular electrostatic interactions. The overall data of DFT analysis for derivatives D1 and D2 are presented in Table 5.

Table 5
Density functional theory (DFT)-derived chemical reactivity descriptors of quercetin derivatives D1 and D2

entry	Computational method (functional/basis set)	B3LYP/6-311G (d, p) quercetin derivative D1	B3LYP/6-311G (d, p) quercetin derivative D2
1	EHOMO / eV	-3.1663	-3.5769
2	ELUMO / eV	-2.3276	-2.6895
3	ELUMO - EHOMO gap / eV	0.8387	0.8874
4	chemical hardness (η) / eV	0.4194	0.4437
5	chemical softness (S) / eV-1	1.1920	1.1269
6	electronegativity (χ) / eV-1	2.7470	3.1332
7	chemical potential (µ) / eV-1	-2.7470	-3.1332
8	electrophilicity index (ω) / eV-1	8.996	11.0679
9	dipole moment (µ) / D	0.7184	6.4998

DFT-based electronic descriptors reveal important details about the nature of the docking process exhibited by both D1 and D2 derivatives. Specifically, smaller HOMO-LUMO energy differences mean that there is more electronic polarizability and charge transfer, which usually implies better molecular recognition and stronger non-covalent interactions between ligands and proteins. For example, the smaller energy gap of D1 implies that this derivative is more reactive and adaptable to the electronic environment of the VP30 binding pocket, which explains its better docking score (-9.1 kcal mol^-1) than D2 (-8.7 kcal mol^-1). Likewise, the computed softness values mean that both derivatives are capable of electronic redistribution during complex formation and, therefore, can form hydrogen bonds and hydrophobic interactions with important residues in VP30, including Arg179, His215, Lys183, Ser184, Ser187, and Glu252.

Additionally, the electrophilicities indicate that both derivatives have electron-acquiring capabilities, which may enhance the stabilization of donor-acceptor interactions inside the active site. HOMO dispersal on the aromatic framework of the flavonoid and the hydroxyl-rich regions suggests sites that potentially donate electrons to electron-deficient amino acid residues, and LUMO shows probable spots where electron density would be acquired from protein sources. These electronic features are in agreement with the wide-ranging hydrogen-bonding interactions and strong binding energies noticed in the docking investigations. Consequently, the DFT results offer a theoretical description for the robust VP30-binding affinities of D1 and D2 and support their choice as the most promising quercetin derivatives identified in this study.

ESP surface analysis

Additional information regarding the interaction between the quercetin derivatives D1 and D2 and the VP30 receptor pocket can be gained from the ESP maps. The ESP surfaces for D1 and D2 are shown in Figure 13, mapped on the surfaces of electron density for both derivatives. As both derivatives have different extremum values of electrostatic potentials, separate ESP range scales were utilized for each derivative in order to depict all the values of potential. These scales are also presented in Figure 13. Red/orange zones denote the electron-rich areas, while blue/cyan zones denote electron-deficient areas. The highly negative electrostatic potential zones, located on the carbonyl groups, phenolic hydroxyls, and glycosidic oxygen atoms, designate these locations as the likely sites of nucleophile attack and hydrogen bonding. Conversely, positive regions were generally localized on those hydroxyl hydrogen atoms where individual molecules would be expected to act as hydrogen donors. These observations are consistent with the docking study results, where the ligands, D1 and D2, created hydrogen bonds with positively polarized or charged amino acids, such as Arg179, Ser184, Glu252, Gln212, Ser187, and Ser254. Therefore, the electrostatic compatibility of the negatively charged oxygen atoms of the ligand and the positively polarized residues of VP30 plays a vital role in complex formation.

Figure 13
Molecular electrostatic potential (ESP) maps of (a) D1 and (b) D2 derivatives of quercetin optimized with DFT. Negative charge areas are shown in red/orange, and positive charge regions are depicted in blue/cyan. The maps were created via distinct electrostatic potential ranges conforming to the minimum and maximum values of each molecule. The individual scale bars are provided along with each surface and should be measured when comparing the electrostatic potential distributions.

For D1, the wider spread of negative electrostatic potential indicates that electron density is more delocalized over both the flavonoid scaffold and sugar moiety, contributing to synchronous intermolecular contacts across the binding pocket of VP30. This distributed charge density correlates with the multiple hydrogen-bond interactions and a higher binding energy obtained for D1. On the other hand, D2 has more confined regions of electron density near catechol and carbonyl groups. This localized charge enrichment could serve to increase local electrostatic attractions and directed hydrogen bonds, which may stabilize a robust interaction network with residues Gln212, Lys218, Ser187, Ser254, and Ser182. Hence, the ESP findings show that both derivatives hold numerous reactive sites capable of contributing to hydrogen bonding, dipole-dipole interactions, and charge-transfer processes vital for suppressing the VP30 transcriptional activator.

Moreover, the joint presence of electron-donating O-functional groups and electron-deficient aromatic carbon atoms gives rise to a polarized molecular surface that promotes dual electrostatic and hydrophobic-type interactions. This balance is especially beneficial for the situation of ligand recognition since it can interact with polar residues and, at the same time, leverage aromatic residues (e.g., His215) through complementary non covalent interactions. The favorable binding energies and MM/PBSA binding free energies obtained for D1 and D2 might be attributed to this electronic arrangement.

Mulliken charge distribution analysis

The Mulliken charge analysis describes quantitatively the electronic polarization that is responsible for the reactive nature of the derivatives D1 and D2. The presence of negative charge density on the oxygen atoms is confirmed in both molecules. The atoms in D1 possessing considerable negative charges include O3 (-0.444), O5 (-0.434), O14 (-0.407), and O20 (-0.415). On the other hand, some oxygen atoms, which carry noticeable negative charges, are O7 (-0.402), O27 (-0.421), O28 (-0.409), O30 (-0.421), O31 (-0.404), and O32 (-0.407) in the case of D2. The negatively charged atoms are likely to form favorable interactions with protonated and positively charged amino acids within the binding site of VP30, especially with the Lys and Arg residues, thus improving the electrostatic stabilization of the ligand-protein complexes.

The positive charge of carbon atoms dispersedly allocated along the aromatic skeleton is a sign of intramolecular charge separation and π-electron polarization. This polarization renders these regions more prone to electrophilic interaction and also facilitates π-related contacts with aromatic amino acids. This results in two different complementary reactive sites competing for reaction during the assembly process: electron-rich oxygen atoms and electron-deficient carbon centers, which thus co-occur to allow simultaneous hydrogen bonding, electrostatic attraction, π-cation interactions, and hydrophobic contacts. Such electronic diversity is in agreement with the extensive interaction profiles that are identified during molecular docking.

A comparison (color-coded) of the two derivatives indicates different types of electronic properties. Derivative D1 possesses better charge delocalization, indicating its higher electronic flexibility and higher capacity for adaptability inside the VP30 cavity. This can be attributed to the presence of a lower HOMO-LUMO gap and higher molecular softness of D1 that favor charge transfer and reactivity. On the other hand, derivative D2 displays a higher degree of charge localization and possesses a larger dipole moment, signifying its stronger directional electrostatic interactions with polar side chains of amino acids. The larger electrophilicity and electronegativity of D2 also confirm its potential for participating in charge transfer interactions with the target. Hence, D1 has better electronic adaptability and interaction capabilities, while D2 has stronger localized electrostatic recognition. Overall, these two complementary electronic features describe the robust binding affinities and stable intermolecular interactions observed for both derivatives with the VP30 protein.

The total results of the Mulliken charge distribution investigation for derivative D1 are provided in Figure 14 and Table 6, whereas the corresponding data for derivative D2 are given in Figure 15 and Table 7.

Table 6
Mulliken atomic charges of carbon, oxygen, and hydrogen atoms in the optimized quercetin derivative D1 calculated at the B3LYP/6 311G(d,p) level of theory

Atom	Mulliken charge / e	Atom	Mulliken charge / e	Atom	Mulliken charge / e
C1	+0.207	O18	-0.389	H35	+0.167
C2	+0.182	C19	+0.148	H36	+0.262
O3	-0.444	O20	-0.415	H37	+0.262
C4	+0.307	C21	-0.122	H38	+0.100
O5	-0.434	C22	-0.029	H39	+0.110
C6	-0.274	C23	+0.374	H40	-0.082
C7	+0.363	C24	+0.089	H41	+0.051
O8	-0.235	C25	+0.093	H42	+0.019
C9	-0.077	C26	+0.007	H43	+0.261
C10	+0.109	C27	+0.050	H44	+0.115
O11	-0.294	O28	-0.201	H45	+0.012
C12	-0.135	O29	-0.382	H46	-0.060
C13	+0.265	O30	-0.328	H47	+0.092
O14	-0.407	O31	-0.065	H48	+0.252
C15	-0.150	O32	-0.064	H49	+0.150
C16	-0.040	H33	+0.279	H50	+0.024
C17	+0.155	H34	+0.122

Table 7
Mulliken atomic charges of carbon, oxygen, and hydrogen atoms in the optimized quercetin derivative D2 calculated at the B3LYP/6 311G(d,p) level of theory

Atom	Mulliken charge / e	Atom	Mulliken charge / e	Atom	Mulliken charge / e
C1	-0.116	C19	+0.141	H37	+0.256
C2	+0.217	C20	-0.133	H38	+0.119
C3	-0.160	C21	-0.045	H39	+0.115
C4	+0.201	C22	+0.020	H40	+0.104
C5	-0.094	C23	+0.128	H41	+0.206
C6	+0.169	C24	-0.149	H42	+0.136
O7	-0.402	C25	+0.077	H43	+0.128
C8	+0.187	O26	-0.315	H44	+0.017
C9	+0.122	O27	-0.421	H45	+0.096
C10	+0.378	O28	-0.409	H46	+0.260
O11	-0.390	C29	+0.064	H47	+0.260
O12	-0.361	O30	-0.421	H48	+0.137
O13	-0.207	O31	-0.404	H49	+0.073
O14	-0.388	O32	-0.407	H50	+0.235
C15	+0.152	O33	-0.386	H51	+0.211
C16	+0.380	H34	+0.123	H52	+0.174
C17	-0.007	H35	+0.120	H53	+0.218
C18	+0.145	H36	+0.260

Figure 14
Mulliken charge distribution and colored electrostatic representation of quercetin derivatives D1. Red-colored regions symbolize electron-rich atoms with negative charge localization, black regions denote almost neutral charge distribution, and green-colored regions represent electron-poor atoms with positive charge delocalization for potential intermolecular interaction with the EBOV VP30 protein.

Figure 15
Mulliken charge distribution and colored electrostatic representation of quercetin derivatives D2. Red-colored regions symbolize electron-rich atoms with negative charge localization, black regions denote almost neutral charge distribution, and green-colored regions represent electron-poor atoms with positive charge delocalization for potential intermolecular interaction with the EBOV VP30 protein.

NBO population analysis

The NBO population analysis of the quercetin analogues D1 and D2 also offered important information regarding the electronic structure, delocalization of charges, and possible interactions with the VP30 protein of the EBOV. The NBO charge maps demonstrate various electrophilic and nucleophilic areas that can play an essential role in antiviral action due to the presence of several polar or charged amino acids in the VP30 protein structure.

The electronic nature of D1 (Figure 16, Table 8) is characterized by a significant delocalization of electrons along the conjugated flavonoid backbone. The oxygen atoms, especially O3, O5, O8, O11, O14, O18, O20, O28, O29, O30, O31, and O32, have a high negative charge (ca. -0.462 to -0.772). As a result, they can be regarded as electron donors that can form hydrogen bonds with positively charged groups or hydrogen bond donors in the vicinity of the VP30 binding site, such as Lys, Arg, Ser, or Glu amino acids. The aromatic ring is characterized by alternating positive and negative charges. This indicates that the ring contains an effective π-electron delocalization system that can potentially stabilize ligand-protein interactions through π-π stacking and hydrophobic contacts. Color-mapping of the molecule D1 supports the fact that there are electron-rich areas present mostly in the functional groups, such as hydroxyl and carbonyl groups, while electron-deficient areas are spread across certain aromatic carbons. Such polarization can increase the possibility of electrostatic complementarity between D1 and VP30. The presence of multiple hydroxyl groups can lead to the formation of multiple hydrogen bonds, which is useful for stable inhibition of viral transcription proteins.

Table 8
NBO atomic charges showing the electron-density distribution over carbon, oxygen, and hydrogen atoms in the optimized quercetin derivative D1 at the B3LYP/6-311G(d,p) level of theory

Atom	Mulliken charge / e	Atom	Mulliken charge / e	Atom	Mulliken charge / e
C1	+0.321	O18	-0.686	H35	+0.263
C2	+0.206	C19	+0.299	H36	+0.473
O3	-0.702	O20	-0.705	H37	+0.471
C4	+0.462	C21	-0.283	H38	+0.201
O5	-0.670	C22	-0.175	H39	+0.214
C6	-0.194	C23	+0.292	H40	+0.145
C7	+0.200	C24	+0.228	H41	+0.167
O8	-0.464	C25	-0.133	H42	+0.257
C9	+0.146	C26	-0.013	H43	+0.509
C10	+0.302	C27	+0.084	H44	+0.430
O11	-0.603	O28	-0.462	H45	+0.412
C12	-0.296	O29	-0.657	H46	+0.408
C13	+0.398	O30	-0.772	H47	+0.282
O14	-0.506	O31	-0.554	H48	+0.497
C15	-0.096	O32	-0.710	H49	+0.224
C16	-0.245	H33	+0.494	H50	+0.013
C17	+0.305	H34	+0.226

Figure 16
Natural bond orbital (NBO) population analysis and charge distribution of quercetin derivatives D1. Red areas denote high electron density or negative charges, whereas green areas denote low electron density or positive charges.

On the other hand, the structure of derivative D2 (Figure 17, Table 9) is characterized by a very polar electronic arrangement. The oxygen atoms, namely O7, O11, O12, O13, O14, O26, O27, O28, O30, O31, O32, and O33, have a negative charge density (from -0.332 to -0.739), suggesting their tendency towards nucleophilicity. In addition, derivative D2 can be predicted to have high hydrogen-bond acceptor properties and higher electrostatic interaction capabilities with the VP30 active site residues. The charge separation in D2 seems more localized rather than widely delocalized, which may enhance site-specific binding energy. The D2 aromatic fused ring system also enhances considerable conjugation and potential π-π interactions with aromatic amino acids (e.g., His215) of VP30. The presence of powerful electron-rich oxygen atoms and a moderate amount of electron-deficient carbon atoms indicates that not only can D2 form hydrogen bonds and dipole-dipole interactions, but also hydrophobic interactions are feasible. These multiple interactions are thought to encourage synergistic binding to antiviral targets.

Table 9
NBO atomic charges showing the electron-density distribution over carbon, oxygen, and hydrogen atoms in the optimized quercetin derivative D2 at the B3LYP/6-311G(d,p) level of theory

Atom	Mulliken charge / e	Atom	Mulliken charge / e	Atom	Mulliken charge / e
C1	+0.318	C19	+0.262	H37	+0.469
C2	+0.375	C20	-0.232	H38	+0.359
C3	-0.231	C21	-0.174	H39	+0.222
C4	+0.392	C22	+0.099	H40	+0.206
C5	+0.324	C23	-0.015	H41	+0.263
C6	+0.365	C24	+0.076	H42	+0.211
O7	-0.519	C25	-0.080	H43	+0.321
C8	+0.337	O26	-0.391	H44	+0.170
C9	+0.171	O27	-0.714	H45	+0.262
C10	+0.453	O28	-0.609	H46	+0.471
O11	-0.684	C29	-0.103	H47	+0.473
O12	-0.652	O30	-0.735	H48	+0.167
O13	-0.332	O31	-0.708	H49	+0.129
O14	-0.539	O32	-0.739	H50	+0.449
C15	-0.106	O33	-0.734	H51	+0.448
C16	+0.421	H34	+0.222	H52	+0.421
C17	-0.193	H35	+0.222	H53	+0.470
C18	+0.268	H36	+0.469

Figure 17
Natural bond orbital (NBO) population analysis and charge distribution of quercetin derivatives D2. Red areas denote high electron density or negative charges, whereas green areas denote low electron density or positive charges.

Comparatively, D1 seems to exhibit a greater extent of global electronic delocalization and conformational flexibility to adapt to the VP30 binding cavity. Conversely, D2 seems to exhibit a greater extent of the localized charge polarization over the molecular structure and a potentially greater hydrogen-bonding efficiency. Regarding the antiviral activity, D2 seems to exhibit a greater extent of binding specificity due to the concentrated electronegative regions, while D1 can exhibit a more extensive interaction coverage with its distributed electron density.

Overall, the NBO population analysis indicates that each of the two quercetin derivatives is electronically suited to inhibit the VP30 transcriptional activator protein of Ebola. The highly electronegative oxygens, conjugated aromatic ring systems, and a balanced electron distribution of these derivatives would result in good binding energy for the VP30 transcriptional activator protein. D2 appears to have a slight edge, due to a more centralized negative charge density and a higher chance of establishing electrostatic and hydrogen bonds with the key amino acid residues of VP30.

To sum up, our results in this study justify experimental investigations to test quercetin derivatives and other related natural compounds as molecules with potential antiviral activity against the EBOV. Growing interest has been observed in the search for small natural compounds as ligands to target viral replication machinery proteins as a therapy against the deadly EBOV, which has been a target of several non-structural viral proteins to be used as drug targets. While in the ongoing efforts towards the development of an effective anti-EBOV vaccine, a perfect therapeutic approach is fairly limited, our focus in this study was set upon the protein VP30, which served as the molecular target after a promising drug binding pocket was discovered. Two quercetin derivatives emerged as potential ligands toward VP30, binding with a relatively good affinity in the active site due to several interactions they formed in the pocket that stabilize the binding interaction, a characteristic of the compounds which was ascribable to the hydroxyl and carbonyl groups capable of stabilizing the interaction by hydrogen-bond or electrostatic interactions. Moreover, the target compounds displayed adequate binding energy and stable interaction patterns within the active pocket of VP30. These results further endorsed the potential of the studied compounds as possible antiviral lead scaffolds. Thus, these investigations provide a sound theoretical foundation upon which in vitro and in vivo experiments can be designed to discover new and effective therapies for EBOV and other viral diseases.

While the current study has given a thorough computational insight into the interactions of quercetin derivatives with EBOV VP30, it is only through actual experiments in the laboratory, like biochemical or antiviral assays, that we can validate the inhibitory potential of these compounds against the virus, as our predictions suggest. Such experiments will also help us assess their viability as potential antiviral drugs.

Structure-activity relationship (SAR) analysis

An initial SAR analysis of the quercetin derivatives (D1 D17) and their molecular docking study against the VP30 protein of the EBOV was performed. It was found that the antiviral activity depended significantly on the number and orientation of hydroxyl substituents, glycosylation, polarity, and electronic effects across the flavonoid moiety.

Derivatives D1 (-9.1 kcal mol^-1) and D2 (-8.7 kcal mol ¹), which showed the most potent activity, had several hydroxyl groups and glycosidic functionalities that were able to interact strongly through hydrogen bonds with amino acid residues like Lys218, Lys283, Glu252, Lys218, Ser184, Arg179, Gln212 and Gly219. As the increase in HBD and HBA groups increased the electrostatic complementarity to the positively charged binding site of VP30, it increased the binding strength of the molecules. However, derivatives having fewer numbers of hydroxyl and more methoxylated ones, such as D13-D17, showed relatively weaker docking energy score values (-5.1 to -5.8 kcal mol^-1).

Data of the electronic descriptors calculated by the DFT technique further proved the SAR results. The smaller the energy gap, the greater the reactivity and affinity. In D1, the smallest energy gap (0.8387 eV) was observed, alongside the highest softness (1.1920 eV^-1). This shows better electronic flexibility and charge transfer within the VP30 binding region. Likewise, higher electrophilicity and dipole moments in D2 favored the formation of stronger electrostatics and intermolecular forces in VP30. Thus, soft electronic structure, high polarity, and charge distribution favor antiviral interaction.

A correlation was also observed between molecular polarity and the efficiency of docking. Derivatives with higher TPSA and those with several oxygen-containing functionalities are expected to show good interactions through hydrogen bonding and van der Waals forces. Yet, inconsistent molecular polarity with high hydrogen bond-forming ability led to the violation of Lipinski’s rules, thus lowering the expected oral bioavailability. This means that a balance should be maintained between binding capability and pharmacokinetics while optimizing leads.

The results obtained from SAR studies indicate that the antiviral effect of quercetin analogs against the Ebola VP30 protein is mainly controlled by:

(i) number of hydroxyl moieties and their spatial arrangement; (ii) glycoside substitution with hydrogen bonding capabilities; (iii) low HOMO-LUMO energy difference along with greater molecular softness; (iv) higher electrophilicity and charge transfer capability, and (v) molecular polarity that facilitates stable protein-ligand complexes.

These insights will serve as a foundation for the rational design of novel inhibitors against Ebola VP30 with potent antiviral properties and superior pharmacokinetics in the future.

Limitations of in silico analyses

Despite the positive ADMET predictions and docking outcomes for the quercetin derivatives D1 and D2, which indicate favorable pharmacokinetics and lower toxicities, it is essential to understand the limitations associated with the in silico approach. This method relies heavily on predictive models and databases, which are useful in establishing preliminary drug likeness, but do not account for the complexity of biological systems. Important considerations like metabolism, biodistribution, off-target activity, and biological dynamics require experimental validation.⁸³,⁸⁴ In addition, although DFT analysis provides useful information about the electronic structure and reactivity of molecules, the results obtained from such calculations are essentially theoretical. In the case of the derivatives D1 and D2, DFT calculations predicted favorable values for certain parameters related to chemical reactivity, which can explain their antiviral effects against Ebola.⁸⁵

In summary, despite the convenience and cost-effectiveness of in silico investigations, it is crucial to interpret their outcomes with caution and validate them experimentally to support the pharmacological properties of the target compounds.

Conclusions

In this study, we used an integrated computational approach encompassing molecular docking, MD simulations, ADMET predictions, and DFT calculations to assess seventeen quercetin derivatives as potential inhibitors of the EBOV VP30 protein. Molecular docking predicted D1 and D2 to be the best candidates with promising binding energies of -9.1 and -8.7 kcal mol ¹, respectively. The two derivatives formed several interactions that stabilize the binding with the key residues in the active site.

The 100 ns MD simulations demonstrated that the VP30-ligand complexes remained stable over time. The binding free energies for D1 and D2 were calculated to be -22.1562 and -19.7180 kcal mol^-1, respectively, via MM/PBSA, while RMSD, RMSF, Rg, and PCA show that structural stability and sustained interaction of ligand-protein occurred during MD simulation.

ADMET analysis predicted acceptable intestinal absorption (34.15% for D1, 44.93% for D2), low acute toxicity (LD₅₀ = 5000 mg kg^-1), and no expected inhibition of the key CYPs; therefore, D1 and D2 had a reasonable safety and metabolism profile even with poor BBB permeability and low oral bioavailability. Additionally, DFT calculations revealed small HOMO-LUMO energy gaps (0.8387 eV for D1 and 0.8874 eV for D2), indicating good electronic reactivity and charge-transfer properties that would favor a strong interaction with the VP30 binding site.

Altogether, the molecular docking and MD simulations, ADMET predictions, and DFT calculations all indicate that D1 and D2 are the most promising quercetin derivatives in targeting the EBOV VP30 protein. These findings provide a promising theoretical basis for subsequent experimental validation and development of quercetin-based anti-Ebola therapeutics.

Supplementary Information

Supplementary Information (docked binding interactions of derivatives D1-D17) is available free of charge at http://jbcs.sbq.org.br as a PDF file.

Supplementary PDF

Supplementary PDF

Acknowledgments

The authors wish to thank Princess Nourah bint Abdulrahman University Researchers Supporting Project (number PNURSP2026R33), Princess Nourah bint Abdulrahman University, Riyadh, Saudi Arabia, for the financial support.

ChatGPT (OpenAI, GPT-5.3) was utilized during the preparation of the graphical abstract for image refinement, layout optimization, and visual presentation suggestions. All scientific content, interpretation, and final design decisions were reviewed and approved by the authors.

Data Availability Statement

All data supporting the findings of this study are available within the article and its Supplementary Information.

References

• ¹ Report of a WHO/International Study Team; Bull. World Health Organ. 1976, 56, 247. [Link] accessed in July 2026
  » Link
• ² Takada, A.; Kawaoka, Y.; Trends Microbiol. 2001, 9, 506. [Crossref]
  » Crossref
• ³ Zhao, Y.; Ren, J.; Harlos, K.; Jones, D. M.; Zeltina, A.; Bowden, T. A.; Padilla-Parra, S.; Fry, E. E.; Stuart, D. I.; Nature 2016, 535, 169 [Crossref]
  » Crossref
• ⁴ Halfmann, P.; Neumann, G.; Feldmann, H.; Kawaoka, Y.; EBioMedicine 2014, 1, 2. [Crossref]
  » Crossref
• ⁵ Moghadam, S. R. J.; Omidi, N.; Bayrami, S.; Moghadam, S. J.; SeyedAlinaghi, S.; Asian Pac. J. Trop. Biomed. 2015, 5, 260. [Crossref]
  » Crossref
• ⁶ Cenciarelli, O.; Pietropaoli, S.; Malizia, A.; Carestia, M.; D’Amico, F.; Sassolini, A.; Di Giovanni, D.; Rea, S.; Gabbarini, V.; Tamburrini, A.; Int. J. Microbiol. 2015, 2015, 769121. [Crossref]
  » Crossref
• ⁷ Wilson, J. A.; Bray, M.; Bakken, R.; Hart, M. K.; Virology 2001, 286, 384. [Crossref]
  » Crossref
• ⁸ Hayman, D. T.; Emmerich, P.; Yu, M.; Wang, L.-F.; Suu-Ire, R.; Fooks, A. R.; Cunningham, A. A.; Wood, J. L.; PLoS One 2010, 5, e11978. [Crossref]
  » Crossref
• ⁹ Goeijenbier, M.; Van Kampen, J.; Reusken, C.; Koopmans, M.; Van Gorp, E.; Neth. J. Med. 2014, 72, 442. [Link] accessed in July 2026
  » Link
• ¹⁰ Basler, C. F.; Wang, X.; Mühlberger, E.; Volchkov, V.; Paragas, J.; Klenk, H.-D.; García-Sastre, A.; Palese, P.; Proc. Natl. Acad. Sci. U. S. A. 2000, 97, 12289. [Crossref]
  » Crossref
• ¹¹ Lindenbach, B. D.; Rice, C. M.; J. Virol. 1999, 73, 4611. [Crossref]
  » Crossref
• ¹² Huang, Y.; Xu, L.; Sun, Y.; Nabel, G. J.; Mol. Cell 2002, 10, 307. [Crossref]
  » Crossref
• ¹³ Kirchdoerfer, R. N.; Moyer, C. L.; Abelson, D. M.; Saphire, E. O.; PLoS Pathog. 2016, 12, e1005937. [Crossref]
  » Crossref
• ¹⁴ Modrof, J.; Becker, S.; Mühlberger, E.; J. Virol. 2003, 77, 3334. [Crossref]
  » Crossref
• ¹⁵ Hartlieb, B.; Muziol, T.; Weissenhorn, W.; Becker, S.; Proc. Natl. Acad. Sci. U. S. A. 2007, 104, 624. [Crossref]
  » Crossref
• ¹⁶ Jun, S.-R.; Leuze, M. R.; Nookaew, I.; Uberbacher, E. C.; Land, M.; Zhang, Q.; Wanchai, V.; Chai, J.; Nielsen, M.; Trolle, T.; FEMS Microbiol. Rev. 2015, 39, 764. [Crossref]
  » Crossref
• ¹⁷ Saeidnia, S.; Abdollahi, M.; Daru 2014, 22, 70. [Crossref]
  » Crossref
• ¹⁸ Ahmad, N.; Rehman, A. U.; Badshah, L. S.; Ullah, A.; Mohammad, A.; Khan, K.; J. Mol. Struct. 2020, 1203, 127428. [Crossref]
  » Crossref
• ¹⁹ Takahashi, K.; Halfmann, P.; Oyama, M.; Kozuka-Hata, H.; Noda, T.; Kawaoka, Y.; J. Virol. 2013, 87, 8862. [Crossref]
  » Crossref
• ²⁰ Nishino, H.; Nagao, M.; Fujiki, H.; Sugimura, T.; Cancer Lett. 1983, 21, 1. [Crossref]
  » Crossref
• ²¹ Nam, J.-S.; Sharma, A. R.; Nguyen, L. T.; Chakraborty, C.; Sharma, G.; Lee, S.-S.; Molecules 2016, 21, E108. [Crossref]
  » Crossref
• ²² Ferreira, L. G.; Dos Santos, R. N.; Oliva, G.; Andricopulo, A. D.; Molecules 2015, 20, 13384. [Crossref]
  » Crossref
• ²³ Materska, M.; Pol. J. Food Nutr. Sci. 2008, 58, 407. [Link] accessed in July 2026
  » Link
• ²⁴ Sotnikova, R.; Nosalova, V.; Navarova, J.; Interdiscip. Toxicol. 2013, 6, 9. [Crossref]
  » Crossref
• ²⁵ Si, Y.-X.; Wang, Z-J.; Park, D.; Jeong, H. O.; Ye, S.; Chung, H. Y.; Yang, J.-M.; Yin, S.-J.; Qian, G.-Y.; Biosci. Biotechnol. Biochem. 2012, 76, 1091. [Crossref]
  » Crossref
• ²⁶ Boulton, D. W.; Walle, U. K.; Walle, T.; J. Pharm. Pharmacol. 1998, 50, 243. [Crossref]
  » Crossref
• ²⁷ Maalik, A.; Khan, FA.; Mumtaz, A.; Mehmood, A.; Azhar, S.; Atif, M.; Karim, S.; Altaf, Y.; Tariq, I.; Trop. J. Pharm. Res. 2014, 13, 1561. [Crossref]
  » Crossref
• ²⁸ Ma, Y. H.; Hong, X.; Wu, F.; Xu, X. F.; Li, R.; Zhong, J.; Zhou, Y. Q.; Liu, S. W.; Zhan, J.; Xu, W.; Acta Pharmacol. Sin. 2023, 44, 1487. [Crossref]
  » Crossref
• ²⁹ Martínez, M. J.; Biedenkopf, N.; Volchkova, V.; Hartlieb, B.; Alazard-Dany, N.; Reynard, O.; Becker, S.; Volchkov, V.; J. Virol. 2008, 82, 12569. [Crossref]
  » Crossref
• ³⁰ Di Petrillo, A.; Orrù, G.; Fais, A.; Fantini, M. C.; Phytother. Res. 2022, 36, 266. [Crossref]
  » Crossref
• ³¹ Boo, H. J.; Yoon, D.; Choi, Y.; Kim, Y.; Cha, J. S.; Yoo, J.; Biomolecules 2025, 15, 313. [Crossref]
  » Crossref
• ³² Raj, U.; Varadwaj, P. K.; Interdiscip. Sci. 2016, 8, 132. [Crossref]
  » Crossref
• ³³ Fanunza, E.; Iampietro, M.; Distinto, S.; Corona, A.; Quartu, M.; Maccioni, E.; Horvat, B.; Tramontano, E.; Antimicrob. Agents Chemother. 2020, 64, e00530-20. [Crossref]
  » Crossref
• ³⁴ Meng, X. Y.; Zhang, H. X.; Mezei, M.; Cui, M.; Curr. Comput.-Aided Drug Des. 2011, 7, 146. [Crossref]
  » Crossref
• ³⁵ Kirchdoerfer, R. N.; Moyer, C. L.; Abelson, D. M.; Saphire, E. O.; PLoS Pathog 2016, 12, e1005937. [Crossref]
  » Crossref
• ³⁶ Protein Data Bank; Research Collaboratory for Structural Bioinformatics (RCSB): Piscataway, NJ, USA. [Link] accessed in May 2026
  » Link
• ³⁷ PerkinElmer Informatics; ChemDraw, version 19.0; PerkinElmer Informatics, Waltham, MA, USA, 2019.
• ³⁸ Labute, P.; Molecular Operating Environment (MOE); 2022.02 Chemical Computing Group ULC, 910-1010. Sherbrooke St. W., Montreal, QC H3A 2R7, Canada, 2023.
• ³⁹ Cao, Y.; CB-Dock Web Server; Sichuan University, Chengdu, China. [Link] accessed in July 2026
  » Link
• ⁴⁰ Trott, O.; Olson, A. J.; AutoDock Vina, version 1.2.0; The Scripps Research Institute, La Jolla, CA, USA, 2021.
• ⁴¹ Biovia, D. S.; Discovery Studio Visualizer, version 2025; Dassault Systèmes: San Diego, CA, USA, 2025. [Link] accessed in July 2026
  » Link
• ⁴² Certara, SYBYL-X Suite, version 2.1.1; Certara Inc., Princeton, NJ, USA, 2012; Certara; SYBYL-X Suite, version 2.1.1, Certara, Princeton, NJ, 2013.
• ⁴³ Schrödinger, LLC; The PyMOL Molecular Graphics System, version 1.8; Schrödinger, LLC, NY, USA, 2015.
• ⁴⁴ Case, D. A. et al.; AMBER 12; University of California, San Francisco, 2012.
• ⁴⁵ Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L.; J. Chem. Phys. 1983, 79, 926. [Crossref]
  » Crossref
• ⁴⁶ Darden, T.; York, D.; Pedersen, L.; J. Chem. Phys. 1993, 98, 10089. [Crossref]
  » Crossref
• ⁴⁷ Perez, A.; MacCallum, J. L.; Brini, E.; Simmerling, C.; Dill, K. A.; J. Chem. Theory Comput. 2015, 11, 4770. [Crossref]
  » Crossref
• ⁴⁸ Berendsen, H. J.; Postma, J. V.; Van Gunsteren, W. F.; DiNola, A. R.; Haak, J. R.; J. Chem. Phys. 1984, 81, 3684. [Crossref]
  » Crossref
• ⁴⁹ Khan, M. T.; Khan, A.; Rehman, A. U.; Wang, Y.; Akhtar, K.; Malik, S. I.; Wei, D.-Q.; Sci. Rep. 2019, 9, 1. [Crossref]
  » Crossref
• ⁵⁰ Martínez, A. M.; Kak, A. C.; IEEE Trans. Pattern Anal. Mach. Intell. 2001, 23, 228. [Crossref]
  » Crossref
• ⁵¹ Abdi, H.; Williams, L. J.; Wiley Interdiscip. Rev.: Comput. Stat. 2010, 2, 433. [Crossref]
  » Crossref
• ⁵² van de Waterbeemd, H.; Gifford, E.; Nat. Rev. Drug Discovery 2003, 2, 192. [Crossref]
  » Crossref
• ⁵³ Guan, L.; Yang, H.; Cai, Y.; Sun, L.; Di, P.; Li, W.; Liu, G.; Tang, Y.; MedChemComm 2018, 10, 148. [Crossref]
  » Crossref
• ⁵⁴ Swiss Institute of Bioinformatics; SwissADME, Swiss Institute of Bioinformatics, Lausanne, Switzerland. [Link] accessed in July 2026
  » Link
• ⁵⁵ Lipinski, C. A.; Lombardo, F.; Dominy, B. W.; Feeney, P. J.; Adv. Drug Delivery Rev. 2001, 46, 3. [Crossref]
  » Crossref
• ⁵⁶ Lohit, N.; Singh, A. K.; Kumar, A.; Singh, H.; Yadav, J. P.; Singh, K.; Kumar, P.; Lett. Drug Des. Discovery 2024, 21, 1334. [Crossref]
  » Crossref
• ⁵⁷ Pires, D. E. V.; pkCSM Web Server, University of Cambridge, Cambridge, UK. [Link] accessed in July 2026
  » Link
• ⁵⁸ Preissner, R.; ProTox 3.0 Web Server; Charité - University Medicine Berlin, Berlin, Germany, 2024. [Link] accessed in July 2026
  » Link
• ⁵⁹ Tandon, H.; Chakraborty, T.; Suhag, V.; Res. Med. Eng. Sci. 2019, 7, 791. [Link] accessed in July 2026
  » Link
• ⁶⁰ Frisch, M. J. et al.; Gaussian 05, revision C.02; Gaussian, Inc.: Wallingford, CT, USA, 2004.
• ⁶¹ Dennington, R. D.; Keith, T. A.; Millam, J. M.; GaussView 5.0.8, Gaussian, 2008.
• ⁶² Becke, A. D.; Phys. Rev. A 1988, 38, 3098. [Crossref]
  » Crossref
• ⁶³ Lee, C.; Yang, W.; Parr, R. G.; Phys. Rev. B 1988, 37, 785. [Crossref]
  » Crossref
• ⁶⁴ Ayalew, M. E.; J. Biophys. Chem. 2022, 13, 29. [Crossref]
  » Crossref
• ⁶⁵ Kirchdoerfer, R. N.; Moyer, C. L.; Abelson, D. M.; Saphire, E. O.; PLoS Pathog. 2016, 12, e1005937. [Crossref]
  » Crossref
• ⁶⁶ Sahu, N.; Madan, S.; Walia, R.; Tyagi, R.; Fantoukh, O. I.; Hawwal, M. F.; Akhtar, A.; Almarabi, I.; Alam, P.; Saxena, S.; Saudi Pharm J 2023, 31, 101788. [Crossref]
  » Crossref
• ⁶⁷ Jeffrey, G. A.; An Introduction to Hydrogen Bonding, vol. 32; Oxford University Press: New York, USA, 1997.
• ⁶⁸ Sarkhel, S.; Desiraju, G. R.;. Proteins: Struct., Funct., Bioinf. 2004, 54, 247. [Crossref]
  » Crossref
• ⁶⁹ Panigrahi, S. K.; Desiraju, G. R.; Proteins: Struct., Funct., Bioinf 2007, 67, 128. [Crossref]
  » Crossref
• ⁷⁰ Bissantz, C.; Kuhn, B.; Stahl, M.; J. Med. Chem 2010, 53, 5061. [Crossref]
  » Crossref
• ⁷¹ Salentin, S.; Schreiber, S.; Haupt, V. J.; Adasme, M. F.; Schroeder, M.; Nucleic Acids Res 2015, 43, W443. [Crossref]
  » Crossref
• ⁷² Biedenkopf, N.; Schlereth, J.; Grünweller, A.; Becker, S.; Hartmann, R. K.; J. Virol. 2016, 90, 7481. [Crossref]
  » Crossref
• ⁷³ Aguiar-Pech, J.; Borges-Argáez, R.; Puerta-Guardo, H.; Pathogens 2025, 14, 1156. [Crossref]
  » Crossref
• ⁷⁴ Torabfam, M.; Celebi Torabfam, G.; Osonga, F.; Dias, C.; Sadik, O.; Sci. Rep. 2025, 15, 43140. [Crossref]
  » Crossref
• ⁷⁵ Yang, B.; Dong, Y.; Wang, F.; Zhang, Y.; Molecules 2020, 25, 4613. [Crossref]
  » Crossref
• ⁷⁶ El Monfalouti, H.; Kartah, B. E. In Biochemistry; Gouvinhas, I.; Carmona-Ribeiro, A. M.; Barros, A. N., eds.; IntechOpen: London, UK, 2024, ch. 5. [Crossref]
  » Crossref
• ⁷⁷ Daina, A.; Zoete, V.; ChemMedChem 2016, 11, 1117. [Crossref]
  » Crossref
• ⁷⁸ Sertbakan, T. R.; CBU J. Sci. 2017, 13, 851. [Crossref]
  » Crossref
• ⁷⁹ Sathyanarayanmoorthi, V.; Karunathan, R.; Kannappan, V.; J. Chem. 2013, 2013, 258519. [Crossref]
  » Crossref
• ⁸⁰ Murugavel, S.; Manikandan, N.; Lakshmanan, D.; Naveen, K.; Perumal, P. T.; J. Chil. Chem. Soc. 2015, 60, 3015. [Crossref]
  » Crossref
• ⁸¹ Hoque, M. J.; Ahsan, A.; Hossain, M. B.; Biomed. J. Sci. Tech. Res. 2018, 9, 7360. [Crossref]
  » Crossref
• ⁸² Kumar, S.; Saini, V.; Maurya, I. K.; Sindhu, J.; Kumari, M.; Kataria, R.; Kumar, V.; PLoS One 2018, 13, e0196016. [Crossref]
  » Crossref
• ⁸³ Sacan, A.; Ekins, S.; Kortagere, S.; Methods Mol. Biol. 2012, 910, 87. [Crossref]
  » Crossref
• ⁸⁴ Shah, A.; Jain, M. In Computer Aided Drug Design (CADD): From Ligand-Based Methods to Structure-Based Approaches; Rudrapal, M.; Egbuna, C., eds.; Elsevier: Amsterdam, Netherlands, 2022, ch. 9. [Crossref]
  » Crossref
• ⁸⁵ Guan, H.; Sun, H.; Zhao, X.; Int. J. Mol. Sci. 2025, 26, 3262. [Crossref]
  » Crossref

Edited by

• Editor handled this article:
  Paulo Augusto Netz (Associate)

Publication Dates

• Publication in this collection
  17 Aug 2026
• Date of issue
  2026

History

• Received
  26 May 2026
• Published
  28 July 2026