Virtual screening of natural products inspired in-house library to discover potential lead molecules against the SARS-CoV-2 main protease

Abstract SARS-CoV-2, a new coronavirus emerged in 2019, causing a global healthcare epidemic. Although a variety of drug targets have been identified as potential antiviral therapies, and effective candidate against SARS-CoV-2 remains elusive. One of the most promising targets for combating COVID-19 is SARS-CoV-2 Main protease (Mpro, a protein responsible for viral replication. In this work, an in-house curated library was thoroughly evaluated for druggability against Mpro. We identified four ligands (FG, Q5, P5, and PJ4) as potential inhibitors based on docking scores, predicted binding energies (MMGBSA), in silico ADME, and RMSD trajectory analysis. Among the selected ligands, FG, a natural product from Andrographis nallamalayana, exhibited the highest binding energy of −10.31 kcal/mol close to the docking score of clinical candidates Boceprevir and GC376. Other ligands (P5, natural product from cardiospermum halicacabum and two synthetic molecules Q5 and PJ4) have shown comparable docking scores ranging −7.65 kcal/mol to −7.18 kcal/mol. Interestingly, we found all four top ligands had Pi bond interaction with the main amino acid residues HIS41 and CYS145 (catalytic dyad), H-bonding interactions with GLU166, ARG188, and GLN189, and hydrophobic interactions with MET49 and MET165 in the binding site of Mpro. According to the ADME analysis, Q5 and P5 are within the acceptable range of drug likeliness, compared to FG and PJ4. The interaction stability of the lead molecules with viral protease was verified using replicated MD simulations. Thus, the present study opens up the opportunity of developing drug candidates targeting SARS-CoV-2 main protease (Mpro) to mitigate the disease. Communicated by Ramaswamy H. Sarma


Introduction
Since the first case of 2019-nCOV infection reported in December 2019, from the local seafood market of Hubei Province of Wuhan, China; it has infected a total of 276,436,619 confirmed cases, including deaths 5,374,744 as of 22 December 2021 worldwide, causing a pandemic (https://covid19.who.int/; Zhu et al., 2020). 2019-nCOV (later named SARS-CoV-2) is a zoonotic, beta coronavirus of bat origin (https://talk.ictvonline.org/ information/w/news/1300/page; . It spreads via airborne respiratory droplets or direct contact with the virus-contaminated surface (Fu et al., 2020).
M pro catalyzes the cleavage of viral peptides, which is required for viral replication and transcription ( Figure 1a). Before processing the overlapping polyproteins pp1a & pp1ab, M pro autocleaves itself between nsp4 & nsp6 (Lee et al., 2020). M pro monomer is basically inactive, but the homodimer is the primary active species with both protomers almost orthogonally aligned to each other. M pro protomer contains three domains, Domain I includes 8-101 amino acid residues, domain II contains 102-184 amino acid residues, and domain III is made up of amino acid residues 201-303 ( Figure 1b). The cleft between domain I & II, termed as substrate binding site, contains CYS-HIS catalytic dyad along with other important amino acid residues, thought to have vital role in proteolytic activity ( Figure  1c). Two protomers bind to each other via N-terminal. Dimerization of M pro could be aborted when substrate binds to the active site of protomeric M pro . Eventually it deactivates M pro catalyzed proteolytic cleavage of polyproteins, which make protomeric M pro as an attractive drug target against SARS-CoV-2 Goyal & Goyal, 2020;Grottesi et al., 2020). Although a number of natural and synthetic products were identified against M pro , but there is still need for more potent, structurally diverse and selective inhibitor's identification for arresting the viral replication.
Traditional drug discovery is a prolong process with the involvement of massive cost, thus computational drug discovery program has been evolved as an alternative route to predict potent lead molecules ( Romano & Tatonetti, 2019;Singh, Bhardwaj, Das, et al., 2021) with a very short rage of time. Although, in silico drug discovery may not truly be accurate, however, indispensable method in catastrophic pandemic conditions like COVID-19. Natural products have played a crucial role in treating and preventing human diseases. Most of the clinical drugs are either naturally occurring or structurally inspired by the natural products (Lahlou, 2013;2007).
In continuation of our efforts in drug discovery exploiting isolation of natural products (Acharya et al., 2018;Basu et al., 2021;Goel et al., 2021;Majumder et al., 2019;Nandi et al., 2012) and synthetic, semi-synthetic modification (Bhattacharjee et al., 2020;Chatterjee et al., 2014Chatterjee et al., , 2017Chatterjee et al., , 2019Dey et al., 2007;Gupta et al., 2014;Kaldhi et al., 2019;Putta et al., 2020), we hereby report the identification of new scaffolds which may possibly abort the viral replication based on their affinity towards SARS-CoV-2 M pro on structure-based virtual screening technique by using the Schr€ odinger suite. A total of 150 molecules have been shortlisted from our in-house library of synthetic and natural products based on the structural resembles of previously reported antiviral agents (Li et al., 2018;Yu et al., 2017;Zemtsova et al., 2011;Zhang et al., 2015). To evaluate the efficacy of the selected library, USFDA approved various Protease inhibitors for viral replication such as: HIV-1 protease inhibitors and clinical candidate for covid-19 M pro Lopinavir, Ritonavir (Cao et al., 2020;Markowitz et al., 1995), Serine protease inhibitor Boceprevir and Cysteine protease inhibitor GC3766 (Hu et al., 2021;Ma et al., 2020;Pavan et al., 2021) were used as standard reference in the docking study. In this study, we have used the reported crystal structure of SARS-CoV-2 M pro (PDB ID: 7BRP) co-crystallized with clinical candidate Boceperevir for virtual screening of the curated in-house library of molecules for identification of leads based on the determination of binding affinity by molecular docking studies, followed by rigorous analysis for binding free energy calculations (MMGBSA). In addition, ADME prediction for top scoring ligands and MD simulation studies in replicates were carried out to understand the trajectory and interaction profile in order to confirm the potential leads against SARS-CoV-2 M pro .

Protein and ligands preparation
We employed molecular docking studies with M pro considering the crystal structure of the receptor (PDB-ID: 7BRP) with resolution 1.80 Å deposited in the Protein Data Bank (www.rcsb.org; Fu et al., 2020). Protein was prepared using the protein preparation wizard tool in Schr€ odinger release 2021-1 (Glide, Schr€ odinger, LLC, New York, NY, 2021). Through the process, the missing side chains and loops were incorporated using Prime. Hydrogen atoms were added after deleting any original ones for amino acid residues and the ligand. The protonation and tautomeric states of ASP, GLU, ARG, LYS and HIS were adjusted to match a pH of 7.4 using ProPka. Possible orientations of ASN and GLN residues were generated. Active site water molecules beyond 5.0 Å from the ligand were deleted, followed by the optimization and minimization of protein conformation (Thakur et al., 2021). Then, the protein-ligand complex was subjected to geometry refinement using an OPLS4 force field restrained minimization with the convergence of heavy atoms to an RMSD of 0.3 Å (Paul et al., 2021).

Receptor grid generation and molecular docking
Binding site analysis was performed using the sitemap module in Schrodinger. A glide grid module was utilized to generate a grid within a 5.0-Å radius of co-crystalized ligand. All the ligand molecules were prepared using the LigPrep wizard tool (LigPrep, Schr€ odinger, LLC, New York, NY, 2021),which can generate several structures from each input structure with various ionization states, tautomer's, stereo chemical characteristics, and ring conformations. The OPLS4 force field was used for the optimization, which produced the low energy conformations of the ligands. A glide module (Glide, Schr€ odinger, LLC, New York, Figure 1. (a) Role of M pro in viral replication and probable mechanism of inhibition (inhibition process is highlighted in the box); (b) Structural details of M pro (nsp5) protomer (Domain I, II and III are designated by orange, green and blue color respectively, cyan color protein fragment indicates the active binding pocket); (c) 3 D representation of M pro ligand binding site (HIS41 and CYS145 catalytic dyad amino acid residues are highlighted in blue color). NY, 2021) was utilized to perform the docking experiment at the active site. Molecular mechanics-generalized born surface area (MM-GBSA) method in Prime was used for rescoring the docked pose of ligand (Prime, Schrodinger LLC, New York, NY, 2011) (Genheden & Ryde, 2015;Zhang et al., 2017). These poses were taken as inputs for the energy minimization of the protein-ligand complexes (E complex ), the free protein (E protein ), and the free ligands (E ligand ). The binding free energy DG bind was determined according to the following equation:

Molecular dynamic simulation
Here we used the SPC water model because it is computationally cost-effective. Initially, the solvated system was minimized for 100 ps (Mark & Nilsson 2001) and the complete system was annealed and equilibrates using ensembles. After minimization, a production run was carried out for 100 ns to investigate the structural modification of the complex. Default relaxation protocols were used. M-shake algorithm was used to keep hydrogen constrained, and the Particle Mesh Ewald method was used to calculate long-range electrostatic interaction (Harvey & De Fabritiis, 2009;Kr€ autler et al., 2001). After the simulation, the trajectory obtained was loaded in the Simulation Interaction Diagram (SID) module of Desmond to analyze the results. Each complex was subjected to specific parameters such as Root Mean Square Deviation (RMSD), Root Mean Square Fluctuation (RMSF), and Protein Secondary Structure Element (SSE), conformational modification of ligands and intermolecular interactions to analyze the level of structural changes (Jena et al., 2021). MD simulations were performed in triplicate for 100 ns using periodic boundary conditions (Khan et al., 2020).

Results and discussion
Protein and ligand preparation and molecular docking studies Boceprevir (Boc) and GC376, FDA approved clinical candidates inhibiting viral proteases, were considered as a standard reference molecule in present study Sacco et al., 2020). Crystal structure of SARS-CoV-2 M pro was adopted as an origin point for computer-aided drug design due to its functional relevance in the survival of CoV. Our docking protocol was validated by re-docking of Boceprevir into active pocket of prepared protein (PDB ID: 7BRP), has shown null deviation (0.00 Å RMSD) from the initial co-crystalize pose reported in literature (Fu et al., 2020;Jain, 2008;Mandal et al., 2016). Along with the in-house library (a set of 150 molecules), four clinical antiviral candidates (Boceprevir,GC376,Lopinavir and Ritonavir) were then docked into the validated active pocket of M pro . The best conformational pose of the each ligand has chosen on the basis of the best docking score and all the best docking score of each ligands were tabulated in decreasing order with the docking cutoff score À3.50 Kcal/mole (Table S1). Furthermore, the interactions of top scorer ligands (5 selected ligands, FG, Q5, P5, PJRJ_37A and PJ4, see Figure 2) with the amino acid residues of Mpro were carefully investigated using the Maestro (Schr€ odinger, LLC, New York, NY, 2021) and Discovery studio visualizer (BIOVIA DS Visualizer 4.5, San Diego: Dassault Syst emes, 2020).
The binding site of M pro is surrounded by hydrophilic as well as hydrophobic amino acid residues, buried inside the key residues HIE41, GLY143, CYS145, HIE163, GLU166, ASP187, ARG188, and GLN192 Fu et al., 2020;Grottesi et al., 2020;Pavan et al., 2021;Rizvi et al., 2021). After exploring the molecular interactions of 5 isolated natural products from in house library, it was perceived that the Flavone glucoside (FG) (5,2 0 ,6 0 -trihydroxy-7-methoxyflavone 2 0 -O-b-D-glucopyranoside), a natural product isolated from Andrographis nallamalayana, manifests hydrogen bonds (H-bond) with the key residues GLU166 (2.21 Å), and ARG188 (2.53 Å) with another residue SER144 (2.68 Å) located at the edge of domain II. The H-bond interactions in FG with binding pocket resembles with Boceprevir (GLU166, 1.88 Å) and GC376 (GLU166, 2.88 Å and SER144, 2.85 Å). The other interactions such as carbon-hydrogen bond (C-H bond), Pi-alkyl bond (ASP187, MET165, and ASN142), Polar and Hydrophobic interactions with the binding pocket residues was observed for the ligand FG. The similar interactions, other than H-bonding, were also found in Boceprevir and GC376 (see Figure 3a,b and Table 1). Interestingly, a strong p-alkyl (3.23 Å) interaction was observed between aryl methoxy group of FG and electron deficient imidazole ring of catalytic dyad histidine (HIS41) residue. However, another catalytic dyad residue, CYS145, only involved in polar interaction through the sugar ring of FG, without any significant bond formation. Since, overall interaction of ligand FG with key amino acid residues of binding pocket were found similar respect to Boceprevir and GC376, they exhibited comparable docking score of À10.31, À11.01 and À9.44 kcal/mole respectively. We presumed that, the slight inferior docking score of FG ligand with respect to Boceprevir may be due to longer distance of H-bond interactions. In contrary, another flavone ligand, 3-methoxy derivative of Luteolin (P5) naturally product from cardiospermum halicacabum, has shown significant alteration of binding nature within the active pocket of M pro . Ligand P5, shows Hbonding with TYR54 (2.75 Å), ASN142 (2.71 Å) and GLY143 (2.47 Å) (Figure 3e,f) along with other interactions involving chromone ring of flavone moiety to form Pi-sulfur bond with residues like MET49 (4.92 Å) and CYS145 (3.53 Å). Interestingly, the same chromone ring is also involved in Pialkyl bonding with MET49 (4.22 Å). In comparison to FG, P5 has shown low docking score of À7.60 kcal/mole, presumably due to the lacks strong polar bonding with the pocket residues in absence of sugar ring (Cherrak et al., 2020). Notably, in previous in silico study of Luteolin have executed similar interactions such as Pi-alkyl and Pi-sulfur with MET49, CYS145 in the active site of M pro (Gogoi et al., 2020).
After exploring binding architecture of natural products, we try to elucidate binding pattern of synthetic molecules in our library within active pocket of M pro . Among the 145 synthetic molecules, ligand Q5 (quinolone scaffold based derivate), PJRJ_37A (catechol dimer) and PJ4 (indole based derivative) have offered the best docking scores À7.75, À7.60, and À7.18 kcal/mole respectively which is significantly lower than FG and reference molecules Boceprevir and GC376 (Table S1). Interestingly, although, those 4 molecules have shown quite similar docking scores, their binding nature with the amino acid residues in the active pocket are quite diverse. In deep analysis has shown that, all 3 synthetic ligands interact with catalytic dyad amino acid residues (HIS41 and CYS145); Q5 and PJ4 have Pi-Pi interaction with both the HIS41 and CYS145, whereas PJRJ_37A only interacts with HIS41 through Pi-bond and no existing interaction was encountered for CYS145 (Figure 3c,d,g,h). Other than catalytic dyad, PJRJ_37A executed H-bond interactions with THR190 distancing 1.87 Å. Whereas, PJ4 able to form two distinct H-bonds with GLU166 (2.57 Å) and ARG188 (1.97 Å) similar to FG. Notably, no H-bond interaction has observed in Q5. We speculate that free hydroxyl or amino functional group of ligand could play a pivotal role in the formation of H-bond within the active pocket of M pro . Additionally, strong Pi-interactions (Pi-Pi, Pi-Alkyl and Pi-S) were found in all four ligands with active pocket amino acids such as LEU27 (Q5), MET49 (Q5, PJRJ_37A, PJ4), MET165 (PJRJ_37A, PJ4) and LEU167 (Q5). Remarkably, the molecule PJ4 was found to be less docking score compare to standard reference molecules, but while analyzing the binding affinities and interaction it attracts the interest of authors to take this for further validation process along with other high scoring molecules. Even though the docking score of Ritonavir, Lopinavir (Protease inhibitors and clinical candidates for SARS-CoV-2 M pro ) are found less interacting in the binding pocket of the M pro based upon docking score (-5.92 and À5.86 kcal/mole) as well as interaction profile with pocket residues. Formed interactions with the binding pocket are found to be weak due to large distance between the ligand and residues when compared to the other docked ligands (Table S1).
Docking of molecules in protein's binding site is a robust approach to elucidate the correct binding pose among several predicted poses of a compound. However, the approach still lacks the ability to correctly rank the affinities of the small molecules to the target protein Choudhary et al., 2020;Sahlgren et al., 2017;Wan et al., 2020). In order to have a better ranking of the ligands, and determination of predictive binding energies, MMGBSA calculations were performed (Yadav et al., 2021). The values obtained after calculation were approximate free energies of binding with a more negative value indicating a stronger binding, also depicts the predictive binding energies of ligand in the binding site of the viral protease. FG exhibited highest negative binding energy of approximately À39.83 kcal/mole with in the binding site and Q15, P_5, PJ_RJ_37A, PJ4 (synthetic derivatives) manifested predictive binding energy of À39.36, À35.25, À24.60, À24.60 kcal/mole in the binding site of M pro (PDB: 7BRP) correspondingly. These scoring results by molecular docking and predictive binding energy studies suggest that the in-house library of natural products and synthetic molecules may bind and inhibit SARS-Cov-2 Main protease enzyme.
The summary of the Hit molecules along with reference molecules is presented in Table 1. While analyzing all the hits were arranged in their descending order of Docking scores and then important amino acid residues were identified from the binding pocket region based on the interaction profile of the ligands. Domain I and Domain II forms the binding pocket of M pro protein with connection of a loop. However based on the ligands interaction profile, we have identified LEU27, MET49, PRO52 and TYR54 are the four amino acid residues from domain I plays a critical role to bind the selected ligands. Whereas, Domain II is found to have H-Bond, Pi-Pi along with hydrophobic interactions. MET165 (Hydrophobic) and GLU166 (H-bond) from domain II were found to be important residues forming interaction with almost all in study ligands. Interestingly, ALA191 from the loop region was also found to be important to form hydrophobic interaction with ligands. Table 1 indicates the M pro active site's interactive amino acid mainly in terms of Hbonding, Hydrophobic interaction and Pi-Pi interacting residues along with docking scores (GScore), pose confirmer score (Emodel) and binding free energy (MMGBSA score) in kcal/mole. Based on a comprehensive study of molecular docking data as well as binding free energy calculations, five compounds were chosen and subjected to dynamic simulations to evaluate their stability.

In silico ADME properties of selected ligands
The ADME (Absorption, Distribution, Metabolism, and Excretion) properties predicted via QikProp (Table S2). The ADME analysis determines the physicochemical properties and biological functions as well as drug likeness of the compound. This is important in terms of assessing the efficacy of drug molecules. The descriptors such as molecular weight, drug-likeness, dipole moment, solvent accessible surface area, hydrogen bond donor, and acceptor traits, octanol/water coefficient, aqueous solubility, binding to human serum albumin, number of likely metabolic reactions, brain/blood partition coefficient and human oral absorption were predicted for top nine scoring ligands obtained from docking studies. The values obtained for all understudy ligands were found to be in recommended range. Selected ligands exhibited good pharmacokinetic descriptors in terms of molecular weight drug-likeness, and human oral absorption. However, the ADME analysis of FG suggests that further derivatization in molecular structure can be done to improve the molecular descriptors such as acceptor hydrogen bond, number of likely reactions, and %human oral absorption. This will eventually lead to an improvement in the drug-likeness score (#star) of the FG. The ADME analysis of PJ4 Indole derivative shows low aqueous solubility and less binding with Human serum Albumin which might be due to non-polar nature of the ligand, further derivatization in molecular structure may lead to improve the predicted descriptors that will eventually lead to an improvement in the drug-likeness score (#star). Quinoline derivatives (Q5, Q14 and Q16) and naturally isolated P_5 exhibited very good pharmacokinetic descriptors in terms of ADME analysis (QikProp, version 4.4, Schr€ odinger, LLC, New York, NY, 2015).

Molecular dynamic (MD) simulation
The molecular dynamics (MD) simulation was performed to find out the stability, confirmation and intermolecular interaction of the ligand molecules with M pro protein (Choudhary et al., 2020). The time-dependent modification of the complexes was calculated over 100 ns using Desmond package. The MD simulation was performed under the thermodynamical conditions (applied volume, pressure and temperature) Khan et al., 2020). Based on the docking results, re-scoring by prime MM-GBSA and visual inspection, five protein ligand docking complexes (FG, Q5, PJRJ_37A, P5, and PJ4) were selected for MD simulations.

RMSD and RMSF
The protein backbone deviation (N, Ca, C) of protein was calculated from the triplicated average RMSD value for the MD simulation (Figure 4a) analysis. The complex structures were highly fluctuated (1 to 2.7 Å) up to equilibration. After equilibration (10 ns), the system gets stabilized, and this trend continued up to 100 ns. The backbone RMSD of the proteinligand docking complex (PDB ID: 7BRP) increased gradually until the 2.10 Å then gets stable till 100 ns. The average RMSD of triplicated simulation is taken for mean average of RMSD calculation and analysis Low RMSD during the simulation indicates the stable complex formation, FG shows good stability as this complex is equilibrated at 5 ns and remains stable throughout the simulation with the least conformational changes with an averaged RMSD of 1.939 Å. After 5 ns, the RMSD of the complexes varied $1 to 2.75 Å. notably, the RMSD of FG, Q5, P5 and PJ4 are low when compared with PJRJ_37A. It indicates that the molecules were highly stable during the MD simulation (Figure 4a). The RMSD analyzes confirm that the showed structural stability during the MD simulation. As suggested via protein backbone RMSD, ligand RMSD was also found stable throughout the simulation with minimal fluctuation. RMSD profile of protein and ligands are shown in Figures S1-S6.
Further, the flexibility of the complexes was analyzed, while the ligands present in the active site of the protein.
The RMSF was used to investigate the fluctuation of the complexes in the function of time. The N-terminal ($8 Å) has high fluctuation compared with C terminal ($6 Å). The FG, Q5, P5 and PJ4 ligands with Mpro protein complexes were highly stable in the catalytic region except for PJRJ_37A, due to the weak intermolecular interaction with the protein (Figure 4b). The complexes exhibit low fluctuation in beta sheets region (80-130); more fluctuations are observed with interactive amino acids such as HIS41, MET49, GLU166, ASP187, ARG188, THR190 and GLN192 is due to the intermolecular interactions of the ligand molecules. Moreover, the intermolecular interactions and secondary structure elements (alpha helices and beta strands) make the protein molecule is slightly rigid. RMSF profile of all understudy ligands are shown in (Figures S1-S6).

Intermolecular interaction
The atomic level information is essential to predict the binding mode of FG, Q5, PJRJ_37A, P5 and PJ4 in the binding site of M pro protein. For binding mode analysis, the intermolecular interactions such as hydrogen bond, hydrophobic contact, ionic interaction and salt bridge were analyzed over 100 ns MD simulation studies. The stability of the protein can be accounted for two main reasons; (1) The coordinates of M pro were extracted from the ligand-bound conformation and (2) in the docked conformation; FG, P5 and PJ4 were tightly bound to M pro active site. Similarly, the complexes are also found to be stable during the simulations. Interestingly other interactions present in Q5 (lacking conventional Hbond) were also found to be very strong as the ligand found to be bound and the complex was stable throughout the simulation study. The interaction diagram showing the detailed interactions of FG, Q5, PJRJ_37A, P5 and PJ4 with M pro is depicted in Figure 5. The H-bond interactions calculated by considering the donor acceptor bonds (D-H-A), shown in Figure 5a (FG), 5b (Q5) and 5c (PJ4), alongside with ionic and hydrophobic interactions are also shown. The percentage occupation of each interaction is also shown in this figure. From Figure 5a, it is evident that the Hydroxy group of 5 th position in flavone moiety is involved intermolecular   H-bonding with flavone ketone with 59% of occupancy along with H-bond with electropositive residue ARG188and polar residue GLN192 (32-45% occupancy). Hydroxy group present in the o-position of phenyl ring forms a strong H-bond with GLU166 comprising 62% occupancy. Other H-bonds which are present more than 30% of total simulation period are ASP187, GLN189 and THR190. Hydrophobic interaction with the non-polar residue MET165 is formed due to methoxy moiety present in the 7 th position of flavone moiety and this interaction remained more than 50% of occupancy along with various Salt bridge interaction with residues of domain II. Another flavonoid class of molecule P5 shows various water bridges along with strong hydrophobic interaction with pocket residues more than 40% of total simulation period (Figure 5g,h). This molecule shows H-bonding with acidic residue ASP187 (56% occupancy) directly as well as through water bridge (34% occupancy) with the hydroxyl group of flavone moiety present at 7 th position. Hydroxy group present at 5 th position forms intermolecular H-bond with ketone present at 4 th position with 60% occupancy. 46% occupancy is observed with H-bond between Hydroxy present at para position of phenyl group present with polar residue THR26. An interesting and strong Pi-Pi interaction is observed with catalytic residue HIS41 with 56% occupancy in final confirmation and 70% throughout the simulation period. Presence of methoxy group at 7 th position is responsible for the hydrophobic interaction with MET165 in FG and sugar moiety is responsible for fitting the molecule into active site with water bridges which lacks in P5 might be the reason for slight higher RMSD and reduction of the binding free energy. Synthetic derivative of Quinoline Q5 however lacks conventional H-bonding with pocket residues but formation of weak H-bonds during the simulation time was observed with HIS163, GLU166, GLN189 and GLN192 for a short period of around 17% occupancy (Figure 5c,d). However, presence of strong Hdrophobic interactions observed with the pocket residues due to highly aromatic ring systems present in ligand makes the complex stable throughout the simulation period. Indole derivative PJ4 found to be stable in the complex and interestingly the interactions observed with the catalytic dyad remained with 70% of total simulation period with HIS41 as Pi-Pi stack with the phenyl ring present at 2 nd position of indole and 35% with the CYS165 as hydrophobic interaction with the phenyl ring of the second indole ring (Figure 5e,f). Ketone in the conjugated chain was found to possess very strong H-bond with the acidic GLU166 (75% occupancy). Electropositive residue ARG188 was found to have 85% occupancy with the secondary amine of indole ring. Hydrophobic interaction of MET165 was observed with phenyl ring for more than 30% of simulation time along with water bridges 20%. Other synthetic molecules such as Catechol dimer was found to be less stable in the complex when compare to the other ligand complexes discussed above. Catechol dimer PJRJ_37A was found to possess Hbond with polar residues of pocket THR190 and GLN192 with 39% and 41% occupancy in final confirmation (see Figure S4) whereas more than 50% contact in overall simulation period. Direct H-bond with electronegative residue ASP187 (48% occupancy) with Hydroxy group present in meta-position. Less functionality and small structure might be the reason for PJRJ_37A to show less interactive profile with pocket residues but it did not destabilize the protein in larger extent during the simulation period. Thus polarity and hydrophobicity plays major role in the binding and stability with M pro and binding to the catalytic site may play a major role in M pro inhibition. The stable conformation obtained from MD simulations could be further utilized for screening large chemical libraries. Protein-ligand complex structures showed a slight deviation, which indicates that these ligands are stable inside the binding pocket and hydrogen bond and other interactions are durable and hold inside the SARS-CoV-2 main protease.

Conclusion
The new findings indicate that computational approaches can be used to find inhibitors of SARS-CoV-2, which is still causing a global pandemic. In this work, 150 natural and synthesized molecules were checked for the target of the Mpro of SARS-CoV-2. These compounds were screened against SARS-CoV-2 M pro with four potential repurposed drugs (Boceprevir,GC376,Ritonavir,and Lopinavir). Our study identified four compounds from an in-house chemical library of natural and synthetic compounds based on docking score, predicted binding energies, Predicted ADME score and averaged RMSD trajectories. SARS-CoV-2 main protease inhibitors are anticipated to be ligands FG, P5 (natural flavone class of products) and Q5, PJ4 (synthetic ligands). Even though the scores of molecular docking are not so high when compared to Standard reference molecules understudy except for the ligand FG which is comparable with the reference, but the interaction profiles of the top ligands attracted the interest of authors to take the ligands for further study as leads. We presume that logical structural modifications of the leads selected may give the potent catalytic inhibitors for the M pro to abort the readily mutating and replicating SARS-CoV-2 virus. To enhance the ADME profile of molecules FG and better binding to the site of M pro further derivatization in structure is suggested and in vitro investigations are required for further confirmation of the findings.