key: cord-295145-ry4e2izd
authors: Das, Pratik; Majumder, Ranabir; Mandal, Mahitosh; Basak, Piyali
title: 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
date: 2020-07-24
journal: Journal of biomolecular structure & dynamics
DOI: 10.1080/07391102.2020.1796799
sha: 
doc_id: 295145
cord_uid: ry4e2izd

The recent outbreak of the coronavirus disease COVID-19 is putting the world towards a great threat. A recent study revealed COVID-19 main protease (M(pro)) is responsible for the proteolytic mutation of this virus and is essential for its life cycle. Thus inhibition of this protease will eventually lead to the destruction of this virus. In-Silico Molecular docking was performed with the Native ligand and the 15 flavonoid based phytochemicals of Calendula officinals to check their binding affinity towards the COVID-19 main protease. Finally, the top 3 compounds with the highest affinity have been chosen for molecular dynamics simulation to analyses their dynamic properties and conformational flexibility or stability. In-Silico Docking showed that major phytochemicals of Calendula officinals i.e. rutin, isorhamnetin-3-O-β-D, calendoflaside, narcissin, calendulaglycoside B, calenduloside, calendoflavoside have better binding energy than the native ligand (inhibitor N3). MD simulation of 100 ns revealed that all the protease-ligand docked complexes are overall stable as compare to M(pro)-native ligand (inhibitor N3) complex. Overall, rutin and caledoflaside showed better stability, compactness, and flexibility. Our in silico (Virtual molecular docking and Molecular dynamics simulation) studies pointed out that flavonoid based phytochemicals of calendula (rutin, isorhamnetin-3-O-β-D, calendoflaside) may be highly effective for inhibiting M(pro) which is the main protease for SARS-CoV-2 causing the deadly disease COVID-19. Rutin is already used as a drug and the other two compounds can be made available for future use. Thus the study points a way to combat COVID-19 by the use of major flavonoid based phytochemicals of Calendula officinals. Communicated by Ramaswamy H. Sarma

The earliest history of human Corona Viruses can be traced back to 1965 when two scientists Tyrrell and Bynoe discovered that they could passage virus stain named B814 (Tyrrell & Bynoe, 1966) . With the onset of 2003, at least five new types of human coronaviruses came into identification, causing huge harm to mankind. Among all of them, the severe acute respiratory syndrome (SARS) coronavirus could be identified as the most harmful till the end of 2019 causing significant morbidity and mortality (Khan et al., 2020) . Some new groups of coronaviruses like the New Haven coronavirus (HCoV-NH Ã ) and NL & NL63 were also identified globally during the year 2004-2005 (Kahn & McIntosh, 2005) .

The end of 2019 marked the emergence of a deadly and highly contagious new coronavirus, the 2019-nCoV which now officially termed as the SARS-CoV-2 or the "Novel Coronavirus" (Mittal et al., 2020) . There has been a rapid outbreak of a pneumonia-like illness mainly originating from Wuhan City, capital of Hubei province, China which was later coined as COVID-19(Al-Khafaji et al., 2020) . Community transmission was soon confirmed in Guangdong Province in China . The World Health Organisation (WHO) has declared "Public Health Emergency of International Concern" on 30 January 2020 and with the outspread of this disease in over 100 countries they declared this as pandemic for the first time in 21 st century and urged all countries around the globe to consider it as a global public health emergency. With Lockdown in most of the countries, the most feasible and possible way to carry out research and serve the purpose of helping mankind is the in silico approach which gives us an overview and hope for new drug discovery effective against such deadly virus.

The Genome or genetic material of the virus composed of a single RNA strand which is identified as the largest among all the RNA viruses reported so far (Elfiky, 2020) . Once the virus infects the cell, its genome act as a messenger RNA (mRNA) and initiates the synthesis of two long polyproteins (pp1a, pp1ab). These synthesized long polyproteins helps the virus to replicate new viruses (Boopathi et al., 2020) . These polyproteins in terms make several RNAs, many more structural proteins that together contribute to the construction of a new virus and two proteases (Sarma et al., 2020) . These protease hence formed are further used for cutting the polyproteins into all essential functional slices (Cui et al., 2019; St. John et al., 2015) . The SARS-CoV-2 is recognized as b-coronavirus and it contains about 29.9 kb genome sequence . The angiotensin-converting enzyme 2 (ACE2) receptors are the major human receptors that are used by this virus to infect the cells (Lu et al., 2020; Nejadi Babadaei et al., 2020) . It is on January 7, 2020, that Chinese scientists attempted to isolate SARS-CoV-2 and finally sequenced the genome SARS-CoV2. He also pointed out a potential drug target protein for the inhibition of SARS-CoV-2 replication ( Figure 1 ). The COVID-19 Main protease (M pro ) is the major protease that plays a key role in the proteolytic mutation of the virus ( Figure S1 ) (Jin et al., 2020a) . Thus blocking the viral polypeptide cleavage would provide us a new direction towards effective treatment of this deadly disease (Gupta et al., 2020) . It is urgent to identify effective and active antiviral compounds that can block the infection at the individual stage and stop its spread worldwide.

India is always considered as the global hub for herbal medicines. More than 5000 natural plants and herbs are used for medicinal purposes (Agarwal & Raju, 2006) worldwide and the majority of them are from India. The Calendula officinals is a noted medicinal plant and used all over the world including Europe, Africa, and India (Bhatnagar & Sastri, 1960) . Calendula officinals or commonly known as marigold contain diverse phytochemicals or bio-active components hence they are used in various therapies starting wound healing to ulcers (Muley et al., 2009) . It is also considered to be a highly antiseptic and anti-bacterial agent and often used as antiseptic cream (Chakraborthy, 2010) . Calendula officinals is also used in homeopathy to treat wounds and pain. Different in vitro studies have also pointed out that Calendula extracts have antiviral, antimicrobial, and anti-genotoxic properties (AshwlayanVD & Verma, 2018) . Considering all the biological activities of Calendula officinals we have tried to screen the bioactive compounds from Calendula officinals against COVID-19 M pro .

The main goal of our in silico studies (Virtual molecular docking and Molecular dynamics simulation) is to identify bioactive compounds from Calendula officinals which can inhibit COVID-19 M pro effectively and hence providing a wayout to combat such deadly disease. This study perhaps first time giving an idea of how calendula can be used in blocking the actions of a deadly virus which is putting the existence of human beings towards a big threat. Further study is indeed required to validate its action on SARS-CoV-2 which is responsible for the spread of the deadly disease "COVID-19", but an initial footstep is what this work would like to provide. A medicine can be made easily available with the continuation of the work.

The crystal structure of a COVID-19 main protease (M pro ), vital for virus replication has been download from the protein data bank (PDB) website (http://www.rcsb.org/pdb) (PDB ID:6LU7: Resolution 2.16 Å) (Jin et al., 2020b) . Major Phytochemicals (majorly flavonoids) of Calendula officinals identified earlier (Muley et al., 2009) were acquired from https://pubchem.ncbi.nlm.nih.gov/ as SDF format. Later on, the SDF format of ligand was converted into mol2 format ( Figure S2 ) using Discovery Studio Visualizer 2017 R2. Prior to the docking, receptor (Covid-19 M pro ) and ligands have been prepared by the Dock Prep tool of UCSF Chimera. Briefly, the receptor and ligands have been optimized firstly by adding hydrogen to the system. Next charges have been assigned, to associate atoms with partial charges and other force field parameters. AMBER ff14SB force field had been selected for all standard residue and Gasteiger force field for all other residues.

The active site of the Covid-19 M pro Protein (PDB ID: 6LU7) has been predicted from different receptor cavities available. Site 1 with coordinates (X¼ À13.669, Y ¼ 15.134, Z ¼ 73.529) has been chosen as the active site as the native inhibitor ligand binds to that site.

Virtual molecular docking and analysis have been performed according to the method described previously (Majumder et al., 2019) . All the small molecules have been docked through autodock vina to active site 1 (X=-13.669, Y ¼ 15.134, Z ¼ 73.529) of Covid-19 M pro . Re-docking has been also performed to validate our virtual molecular docking protocol.

MD simulation has been performed to analyze the dynamic properties and conformational flexibility or stability of apo-M pro (COVID-19 M pro ) and selected docked ligand into the target site. GROMACS. 2019.4 package has been used to run and analyze 100 ns MD simulation. GROMOS96 43a1 in single point charge (SPC) water models and PRODRG server have been used to generate protease and ligand force field and parameter files respectively (Berendsen et al., 1995) , (van Aalten et al., 1996) . Each system has been solvated with water molecule (apo-M pro : 19388, M pro -rutin: 19341, M pro -isorhamnetin-3-O-ß-D: 19374, M pro -calendoflaside: 19356) in separated dodecahedron boxes (x ¼ 9.631, y ¼ 9.631, z ¼ 6.810) with at least 10 Å distance from the edge of the box. Each system has been neutralized by adding 4 NA þ counter ions. The energy minimization of all four systems has been performed through running the steepest descent minimization algorithm with 5000 steps to achieve stable system with maximum force < 1000 kJ mol À1 nm À1 . Prior to the running of real dynamics, the solvent and ions of all the systems have been equilibrated by NVT and NPT ensemble. First, all the systems have been run under an NVT ensemble (constant under number of particles, volume, and temperature) where the systems have been gradually heated from 0 to 300 K and equilibrated for 100 ps. Then all the systems have been run under NPT ensemble (constant number of particles, pressure, and temperature) and equilibrated at 300 K for 100 ps. Here, the backbone C a atoms have been restrained and all the solvent molecules (water and ions) have been allowed to move freely to achieve the solvent equilibrium in the system. The h-bonds of the system have been constrained. The long range electrostatic interactions have been achieved through the particle mesh Ewald method. The systems temperature (300 K) has been controlled through the V-rescale method. Parrinello-Rahman method has been applied to equilibrate and set the pressure (1 bar), density, and total energy of systems. After that, four well-equilibrate systems have been carried out for production run without any restrain for 100 ns with a time step of 2 fs, and after every 10 ps coordinates of the structure have been saved. After the completion of the production run, MD simulation trajectories have been used for various dynamics analysis such as root mean square deviation (RMSD), root mean square fluctuation (RMSF), radius of gyration (Rg), number of hydrogen bonds, columbic interaction energy, principle component analysis (PCA), and secondary structural analysis by different inbuilt scripts of GROMACS. Also, the MM-PBSA binding free energy of receptor-ligand docked complexes have been estimated by using the Molecular Mechanics/Poisson-Boltzmann Surface Area (MMPBSA) method which utilized the MD simulation trajectories and is one of the best-used methods in this kind of analysis (Kumari et al., 2014) .

The drug-likeness property of the finally selected phytochemicals has been predicted to check their pharmacological property and their likeness to be a standard drug. Drug Likness tool (DruLiTo) is an open-access software by the Department of Pharmacoinformatics, NIPER (INDIA) was used for the analysis.

Initially, re-docking of the endogenous or native co-crystal ligand has been performed for the validation of the whole docking procedure ensuring its reproducibility. The re-docking of the native ligand has shown that it binds to the same site of the Covid-19 main protease (M pro ) as the co-crystal ligand binds ( Figure 3 ) and the root mean square deviation (RMSD) value between the co-crystal ligand and re-docked native ligand position is 1.162 Å. RMSD where RMSD values below 2 Å is considered to be successful. Re-docking results depicted that all the major interaction between amino acid residues and co-crystal ligand resembles to the interaction between amino acid residues and the re-docked ligand. The binding energy of the re-docked ligand at the active site of M pro is À7.7 Kcal/mol. The major amino acid interaction includes: Phe140, His164, Asn142, Gly143, Cys145, His163, Met165, Glu, 166, Gln189, Arg188, His41, Thr26, Thr25, Ser144, Leu141, Met49, Ser146, Asp187, Leu27 ( Figure 2 ). Among this 19 amino acid residue 15 matches with that of the co-crystal ligand. Thus we can consider the re-docking to be successful.

A total of 15 major phytochemicals of Calendula officinals have been docked to M pro (COVID-19 main protease). The receptor-ligand interactions have indicated that all the compounds used for molecular docking have a substantial binding affinity towards the receptor as compared to the native ligand. The detailed list is depicted in Table 1 .

The compounds having binding energy greater than the native ligand (-7.7 Kcal/mol) have been taken for further analysis.

Six compounds have been taken into consideration for further analysis as reflected in Table 1 . Rutin, Isorhamnetin-3-O-b-D, Calendoflaside, Narcissin, Calendulaglycoside B, Calenduloside have shown better binding energy than that of the native ligand (Table 1 ). The detailed interactions with the different amino acid residues of M pro have been depicted in Table 2 and also the common interactions of native ligand observed with other ligands are highlighted in Bold front.

Rutin has shown the best binding affinity towards the receptor with a binding energy of À8.8 Kcal/mol and Kd value of 0.336 lM. The interaction is depicted in Figure 3 (a). Docking study has exhibited that Rutin interacts with 16 amino acid residue among which 15 residues (Ser144, His163, Asn142, Cys145, Gly143, His41, Phe140, Thr25, Thr26, Thr190, Arg188, Met165, Glu166, His164, Leu141, Gln189) matches with that to the inhibitor N3. Thus it can be predicted that Rutin binds to the entire amino acid residue needed for proper inhibition of receptor protein (M pro ). Rutin is not only one of the phytochemicals for marigold but many natural products including oranges, black tea, asparagus, buckwheat, onions, green tea, figs, and most citrus fruit. Thus sources of rutin are quite easily available.

Isorhamnetin-3-O-b-D another important phytochemical of Calendula officinals has shown a binding energy score of À8.7 Kcal/mol and Kd value of 0.398 lM which in terms of the native ligand is quite high. The interaction has been depicted in Figure 3 (b). Docking study has been exhibited that Isorhamnetin-3-O-b-D also interacted with 16 amino acid residue among which 13 residues (Cys145, Gly143, Asn142, Ser144, His163, Phe140, Gln189, Asp187, Arg188, Met165, His41, Thr26, Met49) matches with that to the inhibitor N3. Hence we can predict that Isorhamnetin-3-Ob-D also binds to all major amino acid residue which helps in inhibiting the receptor protein.

Calendoflaside has shown a significant binding affinity towards the target receptor with a binding energy of À8.5 Kcal/mol and Kd value of 0.558 lM. The detailed interaction has been depicted in Figure 3 (c). Docking study has shown that Calendoflaside also interacted with 16 amino acid residue among which 15 (Arg188, Asp187, Met165, His163, Ser144, Glu166, Phe140, Leu141, Cys145, Gly143, Asn142, Leu27, Met49, Gln189, His41) coincides with that of the native ligand which gives us a clear idea that Calendoflaside also binds to major amino acid residue responsible for inhibition of COVID-19 main protease (M pro ).

Narcissin a major flavonoid present in many natural products. Docking study with Narcissin has shown that it also exhibits a comparatively high binding affinity towards the receptor with a binding energy of À8.4 Kcal/mol and Kd value of 0.661 lM. A thorough interaction has been depicted in Figure 3 (d). Narcissin has interacted with 13 major amino acid residues among which 11 (Leu141, Asn142, Cys145, Gly143, Leu27, Met49, Asp187, Met165, His41, Gln189, Arg188) amino acid matches with that of the native ligand. Thus it is very clear from the result that Narcissin can also be considered to be a good choice for inhibitor as it also interacts with major amino acid residues which the native ligand interacts with.

Calendulaglycoside B has also exhibited a considerably good binding efficiency towards the receptor with a binding energy of-8.2 Kcal/mol and Kd value of 0.928 lM. The detailed interaction has been depicted in Figure 3 (e). This phytochemical interacts with 16 amino acid residues and among them, 14 amino acid residues interaction (Phe140, Ser144, His163, Glu166, Gln189, Arg188, Asp187, Leu141, His41, Met165, Gly143, Cys145, Asn142, His164) are similar to that of the native ligand. Thus we can infer that Calendulaglycoside B also interacts with major amino acid residues which are similar to the native ligand interaction.

Calenduloside on the other has a good binding affinity towards the target receptor with a binding energy of-7.9 Kcal/mol and Kd value of 1.542 lM which is greater than that of native ligand but not as good as the other ligands. The detailed interaction has been exhibited in Figure 3 (f). Calenduloside interacts with 15 amino acids among which 11 amino acid residues (Thr25, Asn142, Cys145, His164, Gln189, Glu166, Met165, Gly143, Leu27, Thr26, Met49) match that of the native ligand. Thus we can infer that although calenduloside interacts with major amino acid residues similar to that of the native ligand but among the 6 chosen ligand calenduloside is the lowest in terms of efficiency.

Apo-form of COVID-19 main protease (apo-M pro ), top three (Rutin, Isorhamnetin-3-O-b-D, Calendoflaside) docked ligands with higher binding affinity and the crystal structure of COVID-19 M pro with inhibitor N3 have been selected to find out their system stability, flexibility, and other dynamic properties through 100 ns MD simulation.

The Root mean square deviation (RMSD) designates the dynamic stability of all the systems. It also calculates the has showed a comparatively higher RMSD value from 20 ns onwards. The other selected ligands i.e. rutin and calendoflaside showed a relatively less RMSD than that of the native ligand inhibitor N3 and they are almost the same as that of the apo-M pro and all have got almost stable after 75 ns. The M pro -calendoflaside complex has showed a RMSD of slightly 0.4 nm. The crystal ligand i.e. the inhibitor N3 and the protease complex have showed the RMSD above 0.4 nm. The RMSD of isorhamnetin-3-O-b-D has been stable initially around 10 ns and maintained a constant RMSD of around 0.22 nm. Suddenly, it has risen to the RMSD of 0.4 nm during 10-20 ns. After 20 ns, it has shown a stable RMSD of 0.58 nm which indicates that the system got equilibrated after 20 ns concerning the initial structure of the protease-ligand complex. All the complexes are seen to be stable during the 100 ns simulation. But two complexes (M pro -rutin, M pro -calendoflaside) are found to be more stable than the original native co-crystal ligand: inhibitor N3. 

The compactness of the receptor-ligand complexes during the molecular dynamic simulation is best described by the Rg factor. The measurement of the distance between the center of mass of the receptor atoms with its terminal in a specified time frame is the major significance of this study. A compact receptor structure indicates less variation in the gyration value while an expanded or unstable form shows higher Rg value. The study has depicted that the variation in Rg values concerning different time points for the Apo-form of COVID-19 main protease (apo-M pro ) and all the proteaseligands complexes (Figure 4(b) ). 

The root mean square fluctuation (RMSF) value represents the mobility and flexibility of a structure. The RMSF of the amino acid residues of COVID-19 main protease apo-form (apo-M pro ) and protease (M pro ) complexes with rutin, isorhamnetin-3-O-b-D, calendoflaside, and inhibitor N3 have been analyzed (Figure 4(c) ). The average RMSF value for apo-M pro , M pro -inhibitor, M pro -rutin, M pro -isorhamnetin-3-O-b-D, and M pro -calendoflaside are 0.15, 0.19, 0.15, 0.19, and 0.13 nm respectively. It has been observed that apo-M pro , M pro -rutin, and M pro -calendoflaside structure have similar kind of low RMSF value. Whereas crystal structure (M proinhibitor) and M pro -isorhamnetin-3-O-b-D have shown higher RMSF value as compared to others. Among all the predicted ligand calendoflaside and rutin have shown very good RMSF value when they bind with M pro . All the ligands exhibited a higher RMSF peak around 0.7 nm at the C-terminal protease region. From the RMSF value analysis, we can conclude that all of our predicted ligands showed less conformational changes during binding and can act as a stable complex.

The hydrogen bond analysis is a very important factor while considering receptor-ligand stability as these bonds are transient interaction and are responsible for providing stability to the receptor-ligand complex. These are highly specific bonds. In this study, we have calculated all the hydrogen bonds for all the complexes including that of the native ligand inhibitor N3. The number of hydrogen bonds at different time points has been calculated and depicted in Figure 5 the native ligand inhibitor N3. The M PRO -calendoflaside complex showed the maximum number of hydrogen bonds and thus it could be concluded that this complex is the most stable among all the others. Thus calendoflaside is the best choice for ligand considering hydrogen bond into account while the other ligands including the inhibitor N3 also showed good stability towards the complex.

We have analyzed interaction energy to quantify the strength of the interaction of protease-ligand complexes. It can be useful to calculate non-bonded interaction energy between protease and ligands. In Figure 5(b) , the short-range columbic interaction energy of all the ligands with protease has been depicted. In this 100 ns MD simulation study, we 

The binding capacity of the ligand towards the receptor is quantitatively estimated by the binding free energy analysis. Binding free energy is the summation of all non-bonded interaction energies. The binding free energy of the interaction between M PRO (COVID-19 main protease) and docked compounds/ligands has been estimated using the G_MMPBSA tool (Kumari et al., 2014) . Our binding energy analysis from 100 ns MD simulation trajectories have shown that calendoflaside, isorhamnetin-3-O-b-D, and rutin have a higher binding affinity towards COVID-19 M pro inhibition site as compared to the native co-crystal ligand inhibitor N3 as reflected in Figure 5 (c) and Figure 6 , Table 3 ). Results indicate, van der Waals, electrostatic, and SASA energy negatively contribute to the total interaction energy while only polar solvation energy positively contributes to the total free binding energy. Considering the negative contribution, the contribution of van der Waals interaction is much more than that of the electrostatic interactions for all the cases. The contribution of SASA energy is relatively less as compared to the total binding energy. The high negative value of van der Waals energy also points to the massive hydrophobic interaction between the ligands and the COVID-19 main protease (Verma et al., 2016) . The DG binding energy of all the ligands binding to the COVID À19 main protease with respect to different time points has been calculated ( Figure 5(c) ). During the 100 ns MD simulations, there are few fluctuations in the binding energies at different time points. Close analysis shows that inhibitor N3 has more binding energies than rutin and isorhamnetin-3-O-b-D at some time points like 20 ns, 42 ns, and 78 ns but shows lower binding energy than calendoflaside throughout the 100 ns of MD simulation. M PRO -calendoflaside complex has shown very higher energy initially during 10 ns but with time it decreases but kept an average of À250 to À350 kJ/mol. The overall energies although have been found to higher for all the predicted ligands.

When considering structure-based drug designing, the residue-based binding energy is considered to be of utmost importance as based on the residue we can improve the lead compound. Thus the key residues were predicted played an important role in ligand binding. We have analysed the contribution of individual amino acid residues in the binding energy of M PRO -ligand complexes (Figure 5(c) ). From the graph, we have observed that HIS41, MET49, GLN127, SER139, GLU166, GLN189, HIS172 amino acid residues at the inhibition site play a key role in the M pro -ligand binding. Most of these amino acid residues are hydrophilic which indicates that in the receptor-ligand stable interaction non-bonded hydrogen bonds play a major role. 

To analyse the collective motion of all complexes as well as apo-form of COVID-19 main protease (apo-M pro ), the PCA analysis based on C-a atoms have been performed. It has been observed that the first few eigenvectors or the principal components (PCs) of the structures play an important role and describe overall motions of the entire system. Therefore, the first 50 eigenvectors of apo-M pro and its complexes have been considered and plotted (Figure 7(a) ). The first five eigenvectors contribute for 59.24%, 70.93%, 64.32%, 75.10%, 59.78% of the total motion found for the 100 ns MD simulation trajectory of apo-M pro , M PRO -inhibitor N3, M PRO -rutin, M PRO -isorhamnetin-3-0b-D, and M PRO -calendoflaside respectively. It has been found that M PRO -isorhamnetin-3-0-b-D and M PRO -inhibitor N3 have shown very high correlated motions as compared to other complexes as well as apo-M pro . The M PRO -rutin and M PRO -calendoflaside complexes have shown almost similar kind of correlated motions as compared to the apo-M pro . Such data suggest that rutin and calendoflaside have formed very stable complexes with M PRO (COVID-19 main protease) and can be considered as a lead compound.

We have found earlier that the first five eigenvectors contribute the majority portion of the total dynamics of the whole system. Therefore, we have plotted only the first two eigenvectors against each other where each dot represents correlated motions (Figure 7(b) ). The well stable clustered dots signify the more stable structure and low clustered dots represent the weaker stable structure. It has been observed that apo-M PRO , M PRO -rutin, and M PRO -calendoflaside complexes have shown very stable clustered dots as compared to others.

Covariance plots have been generated to understand how residue atoms move relative to each other (Figure 8 ). Motions can be positively correlated (red region), negatively correlated (blue region), or uncorrelated (white region). The deeper color signifies a higher positive correlation or negative correlation. We have found that negatively correlated motions are dominant in all the dynamic systems except M PRO -isorhamnetin-3-0-b-D complex (Figure 8(d) ). It has been observed that correlated residual motions are slightly higher in apo-M PRO as compared to other dynamic systems like M PRO -rutin, M PRO -inhibitor N3, and M PRO -calendoflaside. It has revealed that the binding of ligands induces a slightly different type of motion in M PRO -ligand complexes as compared to apo-M PRO . Here, M PRO -calendoflaside complex (Figure 8(e) ) has shown lesser and M PRO -rutin complex (Figure 8(c) ) has shown a similar pattern of correlated motions as compared to M PRO -inhibitor N3 (crystal structure) (Figure 8(b) ). 

The Gibbs free energy landscape study provides substantial information about the receptor-ligand complex. The conformational state changes during ligand binding can be most effectively predicted by calculating the Gibbs free energy analysis for PC1 and PC2. The Gibbs free energy landscape ( Figure  9 ) shows that the energy value ranging from 0 to 14.8, 0 

The exclusive study is done on the 3 D structural evolution of the protease-ligand complex over 10 ns interval of time scale by comparing 5 snapshots of the protein-ligand complexes ( Figure 10 ). The snapshots are observed to show very minimal changes in the protease structure during the 100 ns simulation. The secondary structure of the protease is conserved throughout the simulation when they are bound with the natural products, which shows the stability of the protease-ligand complex. Interestingly, the secondary structures are conserved during the MD simulations in both states. Additionally, we attribute the large RMSF peak in the residues in the ligand-binding process that signifies to the dynamical change of secondary structure in the bound form.

The overall secondary structural analysis of ligand-bound M pro (COVID-19 main protease) complexes and the apo-M pro has been performed using gmx_do_dssp option in GROMACS ( Figure 11 , Table 4 ). The secondary structural rigidity depends on the percentage of A-Helices/B-Sheets and flexibility depends on the percentage of coils/turns present in the structure. It has been observed that there are slight structural changes in M pro when it binds with We have found that apo-M pro and M pro -inhibitor N3 have a slightly higher percentage of coils as compared to other docked complexes (Table 4 ). Such data indicate that M pro structure is becoming more rigid and stable after docked with compounds like calendoflaside, isorhamnetin-3-O-b-D, and rutin.

Drug-likeness is a crucial consideration while choosing compounds during the initial stages of drug discovery. The pharmacological properties of the final three phytochemicals have been calculated to test their drug-likeness (Table 5 ).

The chemical structures of all these final phytochemicals are depicted in Figure 12 . Although there are few violations but all the phytochemicals can be used as a drug for proper application. Rutin on analysis has passed the Modern drug data report (MDDR Like Rule). MDDR comprises around 100000 drugs that are under development or are already launched in the market. These compounds are referenced in the various sources of literature (Kadam & Roy, 2007) . Rutin (DRUGBANK No.-DB01698) , in general, is also used in some commercially available drugs (Temerk et al., 2006) . Like Rutin, Calendoflaside also passed the MDDR Like Rule. Thus Calendoflaside can also be considered as an important candidate for drug development. On the other hand, Isorhamnetin-3-O-b-D passed the Ghose_Filter, Weighted QED test for drug-likeness. The Ghose rule is a valid rule for selecting drugs based on qualitative and quantitative characterization of known drug databases (Ghose et al., 1999) . A quantitative estimate of drug-likeness or QED depicts the fundamental distribution of molecular properties. Weighted QED is used majorly for molecular target druggability screening by arranging a large set of available bioactive compounds that are already published (Bickerton et al., 2012) . Thus we can consider Isorhamnetin-3-O-b-D as a lucrative candidate for drug design. It is also reported earlier that Isorhamnetin-3-O-b-D is used as a drug for the prevention and/or treatment of diabetes and its complications in the rat model (Lee et al., 2005) . Apart from this, we must consider PAINS Alert while drug design. Pan-assay interference compounds (PAINS) are chemical compounds that frequently give false positive consequences in high-throughput screening (Dahlin et al., 2015) . Rutin although shows 1 alert but the other two compound that is Isorhamnetin-3-O-b-D and Calendoflaside have 0 alert, indicating low or negligible chances for false-positive results.

The Coronaviruses are positive-strand RNA viruses. The process of viral replication and transcription is controlled by two overlapping polyproteins, the pp1a (replicase 1a, 450 kD) and pp1ab (replicase 1ab, 750 kD). Ribosomal frame shifting is a necessary process for the expression of the C-proximal portion of pp1ab polyprotein. Extensive proteolytic processing finally releases functional polypeptides from both the polyproteins which are in term controlled by the main protease (M pro ). Along with the papain-like protease(s) the main protease (Mpro) is essential for proteolytic processing of the polyproteins that are translated from the viral RNA (Xue et al., 2008) . The main protease slices the polyprotein at around 11 conserved sites which involve the Leu-Gln2 (#)(Ser, Ala, Gly) sequences. This slicing process is initiated by the enzyme's autolytic cleavage from pp1a and pp1ab. Thus we can come to the point that inhibiting the activity of the main protease would lead to blockage viral replication (Anand et al., 2003) . Herein, our computational study predicts that flavonoid based inhibitors of Calendula officinals (Rutin, Isorhamnetin-3-O-b-D, Calendoflaside) may block the active site of this protease and hence blocking its action. All the three phytochemicals bind with the amino acid residues Phe140, His164, Asn142, Gly143, Cys145, His163, Met165, Glu166, Gln189, Arg188, His41, Thr26, Thr25, Ser144, Leu141, Met49, Ser146, Asp187, Leu27 and thus blocking the major amino acid residue to initiate further proteolytic activity ( Figure 13 ). The major aim is to stop the viral replication which could be achieved by blocking the action of this main protease. As of now no human proteases possess similar cleavage specificity. Beside that inhibitors are likely to be non-toxic towards human body as they are from a wellestablished natural source and is used as homeopathy medicine for ages.

Our in silico study gives a clear idea that flavonoid based phytochemicals of calendula (rutin, isorhamnetin-3-O-b-D, calendoflaside) may be highly effective for inhibiting M pro which is the main protease for SARS-CoV-2 causing the deadly disease COVID-19. The virtual molecular docking study reveals that seven bioactive compounds from Calendula officinals have a higher binding affinity toward COVID-19 main protease (M Pro ) as compare to co-crystal ligand inhibitor N3. Molecular dynamics simulation (100 ns) studies of the top three docked compounds (rutin, isorhamnetin-3-O-b-D, and calendoflaside according to the docking score) and co-crystal ligand inhibitor N3 suggest that the top three docked compounds have a good amount of stability, flexibility and binding affinity toward COVID-19 main protease (M pro ) as compare to co-crystal ligand inhibitor N3. We have also observed the dynamics properties of the apo-form of COVID-19 main protease (apo-M pro ) and found that there is no significant amount of difference in terms of dynamic properties between apo-M pro and its docked complexes M pro -rutin, M pro -isorhamnetin-3-O-b-D, and M pro -calendoflaside. We have also evaluated the predicted ligands through ADMET descriptors. From the results, it can be concluded that the predicted ligands can be used for medicinal purposes and are also non-toxic.

Overall, this in silico study predicts that these three compounds from Calendula officinals, especially calendoflaside and rutin compounds due to their ability to form a stable complex with M pro , may have the potency to be evolved as an anti-COVID-19 main protease drug to fight against the novel coronavirus but before that, it must go through under the proper preclinical and clinical trials for further experimental and/or clinical scientific validation.

Opportunities for India

Using integrated computational approaches to identify safe and rapid treatment for SARS-CoV-2

Coronavirus main proteinase (3CLpro) structure: Basis for design of anti-SARS drugs

Therapeutic potential of Calendula officinalis

GROMACS: A message-passing parallel molecular dynamics implementation

The wealth of India raw materials (A Dictionary of Indian Raw Materials and Industrial Products)

Quantifying the chemical beauty of drugs

Novel 2019 coronavirus structure, mechanism of action, antiviral drug promises and rule out against its treatment

Phytochemical screening of Calendula officinalis Linn leaf extract by TLC

Prediction of the SARS-CoV-2 (2019-nCoV) 3C-like protease (3CLpro) structure: Virtual screening reveals velpatasvir, ledipasvir, and other drug repurposing candidates

Origin and evolution of pathogenic coronaviruses

PAINS in the assay: Chemical mechanisms of assay interference and promiscuous enzymatic inhibition observed during a sulfhydryl-scavenging HTS

SARS-CoV-2 RNA dependent RNA polymerase (RdRp) targeting: An in silico perspective

A knowledge-based approach in designing combinatorial or medicinal chemistry libraries for drug discovery. 1. A qualitative and quantitative characterization of known drug databases

In-silico approaches to detect inhibitors of the human severe acute respiratory syndrome coronavirus envelope protein ion channel

Structure of M pro from COVID-19 virus and discovery of its inhibitors

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

Recent trends in drug-likeness prediction: A comprehensive review of in silico methods

History and recent advances in coronavirus discovery

Marine natural compounds as potents inhibitors against the main protease of SARS-CoV-2. A molecular dynamic study

g_mmpbsa-a GROMACS tool for high-throughput MM-PBSA calculations

Inhibitory effects of isorhamnetin-3-O-$b$-D-glucoside from Salicornia herbacea on rat lens aldose reductase and sorbitol accumulation in streptozotocin-induced diabetic rat tissues

Genomic characterisation and epidemiology of 2019 novel coronavirus: Implications for virus origins and receptor binding

In vitro and in silico study of Aloe vera leaf extract against human breast cancer

Identification of potential molecules against COVID-19 main protease through structure-guided virtual screening approach

Phytochemical constituents and pharmacological activities of Calendula officinalis Linn (Asteraceae): A review

The expression level of angiotensin-converting enzyme 2 determine the severity of COVID-19: Lung and heart tissue as targets

In-silico homology assisted identification of inhibitor of RNA binding against 2019-nCoV N-protein (N terminal domain)

Targeting zoonotic viruses: Structure-based inhibition of the 3C-like protease from bat coronavirus HKU4 -The likely reservoir

Cathodic adsorptive stripping voltammetric determination of the antitumor drug rutin in pharmaceuticals, human urine, and blood serum

Cultivation of viruses from a high proportion of patients with colds

PRODRG, a program for generating molecular topologies and unique molecular descriptors from coordinates of small molecules

Hydrophobic interactions are a key to MDM2 inhibition by polyphenols as revealed by molecular dynamics simulations and MM/ PBSA free energy calculations

A new coronavirus associated with human respiratory disease in China

Structures of two coronavirus main proteases: Implications for substrate binding and antiviral drug design

We would like to acknowledge the "Centre of Excellence in Phase Transformation and Product Characterization, TEQIP-III", Jadavpur University, Kolkata, MHRD, Indian Institute of Technology, Kharagpur, India, Indian Council of Medical Research, Council of Scientific and Industrial Research and J C Bose Fellowship (DST, India) for providing fellowship and financial support during this work. We would also like to acknowledge every lab member of Biomaterial, Cell Culture, and Microbiology lab of School of Bioscience and Engineering, Jadavpur University. All the authors would like to thank all the medical workers and emergency service providers around the globe for their incomparable contribution towards the society and mankind while combating COVID-19.

No potential conflict of interest was reported by the author(s).

http://orcid.org/0000-0002-8462-4851 Ranabir Majumder http://orcid.org/0000-0003-3585-3276 Piyali Basak http://orcid.org/0000-0001-9245-0014