key: cord-0831836-5puk5vd7
authors: Ibrahim, Mahmoud A.A.; Mohamed, Eslam A.R.; Abdelrahman, Alaa H.M.; Allemailem, Khaled S.; Moustafa, Mahmoud F.; Shawky, Ahmed M.; Mahzari, Ali; Hakami, Abdulrahim Refdan; Abdeljawaad, Khlood A.A.; Atia, Mohamed A.M.
title: Rutin and flavone analogs as prospective SARS-CoV-2 main protease inhibitors: In silico drug discovery study
date: 2021-03-20
journal: J Mol Graph Model
DOI: 10.1016/j.jmgm.2021.107904
sha: 20343f5cc978d4012e36b3754ecf4f4c62058ade
doc_id: 831836
cord_uid: 5puk5vd7

Coronavirus disease 2019 (COVID-19) is a new pandemic characterized by quick spreading and illness of the respiratory system. To date, there is no specific therapy for Severe Acute Respiratory Syndrome coronavirus 2 (SARS-CoV-2). Flavonoids, especially rutin, have attracted considerable interest as a prospective SARS-CoV-2 main protease (M(pro)) inhibitor. In this study, a database containing 2017 flavone analogs was prepared and screened against SARS-CoV-2 M(pro) using the molecular docking technique. According to the results, 371 flavone analogs exhibited good potency towards M(pro) with docking scores less than −9.0 kcal/mol. Molecular dynamics (MD) simulations, followed by molecular mechanics-generalized Born surface area (MM/GBSA) binding energy calculations, were performed for the top potent analogs in complex with M(pro). Compared to rutin, PubChem-129-716-607 and PubChem-885-071-27 showed better binding affinities against SARS-CoV-2 M(pro) over 150 ns MD course with ΔG(binding) values of −69.0 and −68.1 kcal/mol, respectively. Structural and energetic analyses demonstrated high stability of the identified analogs inside the SARS-CoV-2 M(pro) active site over 150 ns MD simulations. The oral bioavailabilities of probable SARS-CoV-2 M(pro) inhibitors were underpinned using drug-likeness parameters. A comparison of the binding affinities demonstrated that the MM/GBSA binding energies of the identified flavone analogs were approximately three and two times less than those of lopinavir and baicalein, respectively. In conclusion, PubChem-129-716-607 and PubChem-885-071-27 are promising anti-COVID-19 drug candidates that warrant further clinical investigations.

The human severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) is a universal menace to world health and economy [1, 2] . Eventually, the World Health Organization (WHO) officially labeled the COVID-19 outbreak as a pandemic on 11 th March 2020 [3] . Based on WHO's proclamation, more than 118 million people around the world have been infected with COVID-19, including some 2.6 million deaths as of 12 th March 2021 [4] .

Unfortunately, the numbers of infected cases are on the increase [5] . Until now, no specific drugs are available to eliminate COVID-19 infection [6] . More recently, remdesivir and dexamethasone were declared to have positive influences in clinical trials [6] [7] [8] . Ongoing clinical trials are being conducted to assess the efficacy of several drugs to treat COVID-19 disease [9] .

Recently, WHO announced 200 vaccine products at various stages of clinical development and issued an Emergency Use Listing (EUL) for Pfizer and AstraZeneca/Oxford COVID-19 vaccines [10] . Although the developed vaccines are currently being distributed and administered, a mass vaccination program will take a long time. Therefore, continued efforts are needed to develop safe and effective anti-COVID-19 drugs.

SARS-CoV-2 M pro (also called 3CL pro ) is one of the most charming targets for preventing viral replication [11] . Very recently, the first crystal structure of SARS-CoV-2 M pro in complex with a peptidomimetic inhibitor (N3) was resolved to enable the rational design of specific inhibitory compounds towards SARS-CoV-2 M pro [12] .

Utilizing in silico drug discovery techniques, several attempts have been made to repurpose known pharmaceutical drugs as potential therapeutic candidates for the treatment of COVID-19 [13] [14] [15] . Virtual screening combined with molecular dynamics (MD) simulations of J o u r n a l P r e -p r o o f chemical libraries against SARS-CoV-2 targets was conducted for the sake of hunting the most promising inhibitors that could block the viral replication [16] .

Natural products have been the focus of several investigations in the search for discovering anti-COVID-19 drug candidates [17, 18] . Among natural products, flavonoids have attracted considerable interest as potential SARS-CoV-2 inhibitors [19] [20] [21] .

Rutin is a flavonol glycoside molecule extracted from many plants, including buckwheat, tobacco, forsythia, hydrangea, viola, etc. Rutin exemplifies a significant component of many antiviral medicines [22] . Rutin is also known to demonstrate good affinity against avian influenza virus [23] , HSV-1 [24] , HSV-2 [25] , and parainfluenza-3 virus [26] . Several studies have documented that rutin has high efficacy as a potent inhibitor to treat COVID-19 infection [27, 28] . The potency of rutin against the SARS-CoV-2 M pro was underpinned and attributed to the hydroxyl groups within the sugar group involved in its structure, demonstrating several noncovalent interactions with the heteroatoms of amino acids of the M pro 's active site [21, 29, 30] . Therefore, the current study was set out to evaluate binding affinities of flavone analogs as SARS-CoV-2 M pro inhibitors using in silico structure-based drug discovery techniques. A database of 2017 flavone analogs was retrieved, prepared and screened virtually against SARS-CoV-2 M pro . Based on the predicted docking scores, the most promising hits were then subjected to molecular dynamics (MD) simulations. The stabilities and affinities of the top potent analogs were further investigated over a 150 ns MD course. Prediction of drug-likeness was carried out to disclose the bioavailabilities of the identified drug candidates. The current data highlight the importance of the identified potential anti-SARS-CoV-2 M pro inhibitors, which can be tested in vitro and in vivo, to fight the global threat of COVID-19.

J o u r n a l P r e -p r o o f

The resolved experimental 3D crystal structure of SARS-CoV-2 main protease (M pro ;

PDB code: 6LU7 [12] ) in complex with peptidomimetic inhibitor (N3) was chosen for all molecular docking and molecular dynamics simulations. Heteroatoms, water molecules and ions were omitted. H++ server was utilized to assign the protonation state of SARS-CoV-2 M pro . As well, all missing hydrogen atoms were appropriately added [31] . The pKa for SARS-CoV-2 M pro residues was investigated under the physical conditions of salinity = 0.15, internal dielectric = 10, pH = 6.5 and external dielectric = 80.

Prior to molecular docking calculations, a set of 2017 flavone analogs was retrieved in SDF format from the PubChem database (https://pubchem.ncbi.nlm.nih.gov). The investigated flavone analogs were selected based on the chemical skeleton of rutin, as shown in Figure 1 . The 3D chemical structures of the retrieved analogs, in addition to rutin (PubChem-528-080-5), were generated with the help of Omega2 software [32, 33] . For each compound, a conformational search was first carried out and all possible conformers within the energy window value of 10 kcal/mol were generated. The lowest energy conformer was selected and subsequently minimized using an MMFF94S force field with the help of accessible SZYBKI software [34, 35] . A schematic representation of the utilized in silico techniques and the filtration process of the database is depicted in Figure 1 . 

AutoDock4.2.6 software was employed to carry out all molecular docking calculations [36] . AutoDock protocol was followed to prepare the pdbqt file of SARS-CoV-2 M pro [37] . The genetic-algorithm number (GA) and the maximum number of energy evaluations (eval) were J o u r n a l P r e -p r o o f adjusted to 250 and 25,000,000, respectively. All other options were kept at their default parameters. The dimensions of the docking grid were defined to enclose the active site of SARS-CoV-2 M pro (60 Å × 60 Å × 60 Å). As well, the grid spacing value was set to 0.375 Å. The grid center's coordinates were positioned at −13.069, 9.740, 68.490 (XYZ assignments, respectively).

Atomic partial charges of the investigated analogs were assigned using the Gasteiger method [38] . The probable binding modes for each studied analog were processed by built-in clustering analysis (2.0 Å RMSD tolerance) and the conformation with the lowest energy in the largest cluster was picked out as representative.

All molecular dynamics (MD) simulations were executed for the most promising flavone analogs complexed with SARS-CoV-2 M pro using AMBER16 software [39] . General AMBER force field (GAFF2) [40] and AMBER force field 14SB [41] were utilized to describe the flavone analogs and M pro , respectively. In the current study, implicit and explicit MD simulations were conducted. In implicit MD simulations, the AM1-BCC method was employed to calculate atomic partial charges of the analogs [42] . There were no cutoff or periodic boundary conditions applied for nonbonded interactions (more precisely, a cutoff value of 999 Å was utilized). Solvent effect was considered using igb=1 implicit solvent model [43] . The docked analog-M pro complexes were initially energy minimized for 500 steps and gradually heated from 0 K to 300 K over 10 ps under NVT condition, utilizing Langevin thermostat. Eventually, 250 ps in addition to 1,000 ps production stages were carried out and snapshots were recorded every 1 ps, giving 250 and 1,000 snapshots, respectively. The CPU version of pmemd (pmemd.MPI) in AMBER16 was used to perform all implicit MD simulations.

In explicit MD simulations, the restrained electrostatic potential (RESP) approach at the HF/6-31G* level was applied to calculate atomic partial charges of the analogs with the assistance of Gaussian09 software [44, 45] . The analog-M pro complexes were solvated in a cubic water box with a minimum distance to the box edge of 15 Å using the TIP3P water model with periodic boundary conditions [46] . Combined steepest descent and conjugate gradient method was employed to perform energy minimization for 5000 steps on the solvated analog-M pro complexes. The investigated systems were then gently annealed from 0 K to 300 K over 50 ps with a weak restraint of 10 kcal mol −1 Å −1 on the M pro . Besides, the systems were adequately equilibrated for 1 ns, and production stages were carried out under the NPT ensemble over simulation times of 5 n, 10 ns, 50 and 150 ns. Snapshots were collected every 10 ps, giving 500, 1000, 5000, and 150000 snapshots, respectively. Long-range electrostatic forces and energies were calculated with Particle Mesh Ewald (PME) method and the Lennard-Jones interactions were estimated with a 12 Å cutoff [47] . Langevin thermostat with a gamma_ln collision frequency set to 1.0 was applied to conserve the temperature at 298 K. A Berendsen barostat was employed for the pressure control with a pressure relaxation time of 2 ps [48] . SHAKE option to constrain all bonds involving hydrogen atoms was used with a time step of 2 fs [49] . Moreover, coordinates and energy values were gathered every 10 ps over the production stage for binding energy calculations and post-dynamics analyses. The GPU version of pmemd (pmemd.cuda) in AMBER16 was employed to execute all explicit MD simulations. All calculations, including molecular docking, molecular dynamics, and quantum mechanics, were performed on the CompChem GPU/CPU cluster (hpc.compchem.net). Molecular graphics were carried out using BIOVIA DS Visualize 2020 [50].

J o u r n a l P r e -p r o o f

The binding free energies of the most potent flavone analogs in complex with SARS-CoV-2 M pro were evaluated with the help of molecular mechanical-generalized Born surface area (MM/GBSA) approach [51] . In the current study, the modified GB model suggested by Onufriev (igb = 2) was utilized to determine the polar solvation energy. Uncorrelated snapshots were taken every 10 ps throughout the MD course, and the MM/GBSA (ΔG binding ) energy calculated as follows:

In silico drug-likeness prediction of the most potent flavone analogs was executed using J o u r n a l P r e -p r o o f

Virtual screening is a comprehensive approach at the initial stage of the drug discovery pipeline that permits discovering promising bioactive inhibitors at a great throughput [53] .

Because of its quickness and cost-efficiency, it is a recommended technique to recognize potential drug candidates against the universally spreading SARS-CoV-2 virus, mainly when time is of the essence.

Based on rutin's chemical structure, the PubChem database was explored, and compounds with similar substructures were identified, giving a total number of 2017 flavone analogs (see Table S1 . It can be seen from the data in Table S1 that 371 analogs displayed docking scores less than −9.0 kcal/mol. Evaluated docking scores, 2D chemical structures, number of conformations in the largest cluster and binding features for nine most potent analogs with SARS-CoV-2 M pro are summarized in Table 1 . Corresponding data for rutin is also considered for the purpose of comparison. As well, the 2D representations of interactions of those nine potent analogs with the essential amino acid residues of SARS-CoV-2 M pro are illustrated in Figure S1 . It is worth mentioning that those nine potent analogs were picked out based on the estimated 50 ns MD/MM/GBSA binding energy calculations described in latter sections.

J o u r n a l P r e -p r o o f From docking features in Table 1 and Figure S1 , it is apparent that all identified flavone analogs demonstrated similar binding modes, exhibiting a substantial hydrogen bond with GLU166 and other hydrogen bonds with several amino acid residues inside M pro 's binding pocket. Other interactions were also spotted, such as van der Waals, hydrophobic and pi-based interactions between the investigated analogs and M pro , giving a docking score less than −9.0 kcal/mol ( Figure S1 ). Table S2 .

The data in Table S2 interestingly displayed 227 analogs with binding energies (ΔG binding ) less than −40.0 kcal/mol. These identified analogs complexed with SARS-CoV-2 M pro were further inspected over 1,000 ps MD simulations to shortlist prospective potent analogs.

MM/GBSA binding energies estimated for the 227 analogs over 1,000 ps MD simulations are listed in Table S3. As Table S3 shows, about half of the screened compounds (≈46%) exhibited considerable binding energies (ΔG binding less than −50.0 kcal/mol). Those potent 105 analogs were further subjected to 5 ns MD simulations in explicit water solvent towards more reliable binding affinities of analog-M pro complexes. The corresponding 5 ns MD/MM/GBSA binding energies calculations were evaluated and presented in Table S4 . It can be seen from data in Table   S4 that only 33 compounds showed considerable binding energies (ΔG binding less than −55.0 kcal/mol).

Furthermore, 10 ns MD simulations in explicit water solvent, followed by MM/GBSA binding energy calculations, were performed on the top 33 potent analogs complexed with M pro (Table S5) . It is apparent from Table S5 that nine potent analogs showed MM/GBSA binding energies (ΔG binding ) less than −55.0 kcal/mol. Finally, those nine potent analogs were picked out and submitted to molecular dynamics (MD) simulations over 50 ns, and the corresponding binding energies were calculated ( Table 2 and Figure 3 ). As shown in Table 2 In order to unveil the main driving force in the binding of the identified flavone analogs with SARS-CoV2 M pro , decomposition of the MM/GBSA binding energies was executed ( Table 3) .

As can be seen from Table 3 

Structural and energetic analyses were executed over 150 ns MD simulations to demonstrate the analog's stability inside the SARS-CoV-2 M pro active site. The structural and energetic analyses involved binding energy per frame, hydrogen bond length, root-mean-square deviation (RMSD), and center-of-mass (CoM) distance.

The correlation between the binding energy per frame and time for PubChem-129-716- 

Hydrogen bond analysis was performed on the assembled trajectories over the 150 ns MD course. The results are listed in Table 4 The hydrogen bonds are investigated by the acceptor-donor atom distance of less than 3.5 Å and acceptor-H-donor angle of higher than 120°. b Only hydrogen bonds with occupancy higher than 50% were listed.

In order to obtain a more in-depth insight into the stability of inhibitor-M pro throughout the 

To estimate the stability of the identified analogs over the MD simulations, root-mean- 

In order to assess the activity and bioavailability of the discovered flavone analogs,

Lipinski's rules were utilized using Molinspiration cheminformatics online software Table 5 .

J o u r n a l P r e -p r o o f A striking feature of the figures in Table 5 is the milogP being less than five for the three 

Lopinavir 

As can be observed from data given in Table 6 and Figure (Table 6 ). It is apparent from Table 6 

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

The continuing 2019-nCoV epidemic threat of novel coronaviruses to global health -The latest 2019 novel coronavirus outbreak in Wuhan, China

A new coronavirus associated with human respiratory disease in China

WHO Director-General's Opening Remarks at the Media Briefing on COVID-19-11th

WHO Coronavirus Disease (COVID-19) Dashboard. Available online

Clinical characteristics and imaging manifestations of the 2019 novel coronavirus disease (COVID-19):A multi-center study in Wenzhou city

Coronavirus puts drug repurposing on the fast track

Effect of Remdesivir vs Standard Care on Clinical Status at 11 Days in Patients With Moderate

Dexamethasone in Hospitalized Patients with Covid

A real-time dashboard of clinical trials for COVID-19

World Health Organization, DRAFT landscape of COVID-19 candidate vaccines

Crystal structure of SARS-CoV-2 main protease provides a basis for design of improved alpha-ketoamide inhibitors

Structure of M(pro) from SARS-CoV-2 and discovery of its inhibitors

In-silico drug repurposing and molecular dynamics puzzled out potential SARS-CoV-2 main protease inhibitors

Drug Repurposing for Candidate SARS-CoV-2 Main Protease Inhibitors by a Novel In Silico Method

In Silico Evaluation of Prospective Anti-COVID-19 Drug Candidates as Potential

SARS-CoV-2 Main Protease Inhibitors

Putative Inhibitors of SARS-CoV-2 Main Protease from A Library of Marine Natural Products: A Virtual Screening and Molecular Modeling Study

Natural-like products as potential SARS-CoV-2 M(pro) inhibitors: in-silico drug discovery

In silico drug discovery of major metabolites from spices as SARS-CoV-2 main protease inhibitors

Flavonoids: promising natural compounds against viral infections

Potential bioactive glycosylated flavonoids as SARS-CoV-2 main protease inhibitors: A molecular docking and simulation studies

Flavonoids with inhibitory activity against SARS-CoV-2 3CLpro

The Pharmacological Potential of Rutin

Structural insights into the South African HIV-1 subtype C protease: impact of hinge region dynamics and flap flexibility in drug resistance

In vitro anti-HIV and -HSV activity and safety of sodium rutin sulfate as a microbicide candidate

Antiherpetic activities of flavonoids against herpes simplex virus type 1 (HSV-1) and type 2 (HSV-2) in vitro

Antibacterial, antifungal, and antiviral activities of some flavonoids

In-Silico approach for identification of effective and stable inhibitors for COVID-19 main protease (M(pro)) from flavonoid based phytochemical constituents of Calendula officinalis

In-silico investigation of phytochemicals from Asparagus racemosus as plausible antiviral agent in COVID-19

Possible SARS-coronavirus 2 inhibitor revealed by simulated molecular docking to viral main protease and host toll-like receptor

Structure-based lead optimization of herbal medicine rutin for inhibiting SARS-CoV-2's main protease

H++: a server for estimating pKas and adding missing hydrogens to macromolecules

Conformer generation with OMEGA: algorithm and validation using high quality structures from the Protein Databank and Cambridge Structural Database

MMFF94s option for energy minimization studies

AutoDock4 and AutoDockTools4: Automated docking with selective receptor flexibility

Computational protein-ligand docking and virtual drug screening with the AutoDock suite

Iterative partial equalization of orbital electronegativity -a rapid access to atomic charges

Development and testing of a general amber force field

Simmerling, ff14SB: improving the accuracy of protein side chain and backbone parameters from ff99SB

Fast, efficient generation of high-quality atomic charges

AM1-BCC model: II. Parameterization and validation

Implicit solvent models

A well-behaved electrostatic potential based method using charge restraints for deriving atomic charges -the RESP model

Gaussian 09

Comparison of Simple Potential Functions for Simulating Liquid Water

Particle mesh Ewald: AnN⋅log(N) method for Ewald sums in large systems

Moleculardynamics with coupling to an external bath

Settle -an Analytical Version of the Shake and Rattle Algorithm for Rigid Water Models

Combined molecular mechanical and continuum solvent approach (MM-PBSA/GBSA) to predict ligand binding

Rate-limited steps of human oral absorption and QSAR studies

Integration of virtual and high-throughput screening

Role of molecular dynamics and related methods in drug discovery

Molecular dynamics simulations in drug design

Re-assessing the rule of 5, two decades on

Lack of antiviral activity of darunavir against SARS-CoV-2