Inhibitory effect of thymoquinone from Nigella sativa against SARS-CoV-2 main protease. An in-silico study SARS-CoV-2.

Nigella sativa is known for the safety profile, containing a wealth of useful antiviral compounds. The main protease (Mpro, 3CLpro) of severe acute respiratory syndrome 2 (SARS-CoV-2) is being considered as one of the most attractive viral target, processing the polyproteins during viral pathogenesis and replication. In the current investigation we analyzed the potency of active component, thymoquinone (TQ) of Nigella sativa against SARS-CoV-2 Mpro. The structures of TQ and Mpro was retrieved from PubChem (CID10281) and Protein Data Bank (PDB ID 6MO3) respectively. The Mpro and TQ were docked and the complex was subjected to molecular dynamic (MD) simulations for a period 50ns. Protein folding effect was analyzed using radius of gyration (Rg) while stability and flexibility was measured, using root means square deviations (RMSD) and root means square fluctuation (RMSF) respectively. The simulation results shows that TQ is exhibiting good binding activity against SARS-CoV-2 Mpro, interacting many residues, present in the active site (His41, Cys145) and also the Glu166, facilitating the pocket shape. Further, experimental approaches are needed to validate the role of TQ against virus infection. The TQ is interfering with pocket maintaining residues as well as active site of virus Mpro which may be used as a potential inhibitor against SARS-CoV-2 for better management of COVID-19.


Introduction
The main protease (Mpro, 3CLpro) is an important drug targets of sever acute respiratory syndrome 2 (SARS-CoV-2) proteome, processing the polyproteins. A number of studies has well characterized the Mpro (Anand et al., level mechanism underlined. The MD simulations are time efficient and may accurately predict how interact (Liu et al., 2018;Hollingsworth and Dror, 2018). These MD studies are useful to capture a variety of biomolecular interactions, including ligand binding, protein folding, and changes in proteins behavior over time.
Knowing the strength of MD simulations and the importance of pharmacological characteristics of N. sativa, we performed the current study on thymoquinone (TQ) against Mpro to analyze the behavior of target protein at molecular level and their molecular effect on the SARS-CoV-2 target proteins.

Protein preparation
The worldwide biomolecular structural information is being archived at Brookhaven National Laboratories, called Protein Data Bank (PDB) (Bernstein et al., 1977;Berman et al., 2000). The researchers around the world can easily retrieve the crystal structures of biomolecules. The crystal structure of COVID-19 Mpro (PDB ID: 6M03) was retrieved from PDB. Prior to further analysis, the structure was subjected to Preparation, using MOE (molecular operating environment) (Vilar et al., 2008). The partial charges and missing hydrogen were assigned. The thymoquinone (CID:10281) was also prepared. Protein and ligand were docked and the complex was subjected to dynamics study.

Molecular Dynamics (MD) simulation
The molecular dynamics (MD) simulations in drug discovery capture the behavior of proteins in full molecular and atomic detail to an extent of fine temporal resolution (Liu et al., 2018;Hollingsworth and Dror, 2018). The Mpro and TQ docking complex were subjected to molecular dynamics (MD) simulation as described in our previous study using Amber package (Berendsen et al., 1995;Khan et al., 2020). Briefly, MD simulation was performed on docking complex with the ff14SB force field through Amber14 package (Salomon-Ferrer et al., 2013;Sun et al., 2014b, a). To solvate each system the TIP3P water model was applied while system was neutralized with counterions (Jorgensen et al., 1983). The system was energy minimized and conjugate gradient followed by heating upto 300K and 1atm pressure to equilibrate the system. Temperature regulations was achieved with the Langevin thermostat while Particle Mesh Ewald algorithm was applied for long-range electrostatic interactions (Essmann et al., 1995;Darden et al., 1993). The MD simulation production step was carried with pmemd code 30 (Götz et al., 2012).

Results and Discussion
The current study shows that TQ may be effective against Mpro of SARS-CoV-2. The calculation of drug-likeness may help to understand the pharmacokinetic of a novel compound as well the pharmaceutical properties before its clinical application. The TQ drug likeness properties, 2005;Liang, 2006;Gan et al., 2006;Xue et al., 2008;Pillaiyar et al., 2016) and its potential role. This viral protein acts on 11 positions along polyproteins. Human proteases do not share the cleavage specificity with NSP5 of SARS-COV-2 (Lee et al., 2020, p. 2).
The substrate-binding sites residues 10-99 (Domains I) and 100-182 (Domain II) in picornavirus, are six-stranded antiparallel β-barrels while domain III (198-303) forming five helices, regulating the dimerization of Mpro (Shi and Song, 2006). Residues Cys145 and His41 form the catalytic site. The catalytic activity depends on the dimerization of the enzyme, as the N-finger interacts with Glu166 to facilitate the S1 pocket shape of the substrate-binding site (Anand et al., 2002). The residue T285 and I286in CoV-2 Mpro is substituted by A285 and L286 respectively (Zhang et al., 2020) leading to a threefold upsurge Mpro (Lim et al., 2014). Inhibitors may be useful to reduce the catalytic degree against these locations.
In the drug discovery and their mechanism of action are important for better understanding the insight of the molecules. The molecular biologist desire to know that how a protein and small molecules works. An atomic level information is typically generating significant insight information of biomolecular interactions. The intermolecular interactions could be explored through the dynamic's studies. Unfortunately, such kind of information is difficult to obtained through experimental approaches. An alternative to such approaches is computational molecular dynamics simulation (MD) of proteins and natural compounds to understand the atomic calculated through Swiss ADME, is also in accordance the drug likeness rule (Walters and Murcko, 2002;Daina et al., 2017). The drug likeness properties and absorption has been given in Table 1. The pharmacokinetics of TQ seems in accordance with desire drug like compound. Similarly, drug likeness properties also favor its clinical application (Brüstle et al., 2002;Vistoli et al., 2008;Ursu et al., 2011). The TQ docked against SARS-CoV-2 Mpro seems potent, altering the stability of protein. The RMSD of TQ and Mpro in Figure 1 shows stability. The Mpro is not stable in the whole simulation period which might be useful for effective inhibition of viral activity.
The RMSD graph during the 50ns simulation period shows that the Mpro exhibited an unstable fluctuation. The RMSD at the start (1.06Å) is rising to 2.9Å at 18ns simulations. A downfall fall in fluctuation was again observed to 1.4Å at the end of 50ns MD simulation. The RMSD of Mpro seems highly unstable due to TQ which may alter the stability of target, assisting in the inhibition of viral proteins. Flexibility is also one of the important thermodynamic properties, maintaining the optimal functions of proteins (Nagasundaram et al., 2015). A large change in this property may alter the biomolecules optimal function. The TQ may cause an increase in the flexibility of Mpro (Figure 2) which might be useful for better management of SARS-CoV-2 infections. The Mpro exhibited the RMSF among 5Å and 25Å at residues position 48 and 310 respectively. Residues at location 145 to 160 also attained a high RMSF (20.4Å) which contains the active site residue Cys145. MD simulations may explore the insight mechanisms of changes at the molecular level (Liu and Yao, 2010;Liu et al., 2018;He et al., 2018) which might be difficult through experimental work. Several studies reported that any change in protein function might be due the change in RMSF (Berhanu and Masunov, 2011;Chong et al., 2011;Bavi et al., 2016).
The degree of folding stability could be measured through Rg. Fluctuations in Rg with respect to time shows unstable folding while a straight value reveals stability in folding (Lobanov et al., 2008;Smilgies and Folta-Stogniew, 2015;Khan et al., 2019Khan et al., , 2021. A protein with misfolding shows variations in Rg over time (Figure 3). The plot shows large variations between 22Å and 22.8Å. Majority of the variations have been found from 811ps to 2026ps.    The lowest Rg was detected at 811ps (22Å) while the highest at 1216ps (22.8Å). This shows that the TQ may affect the folding of Mpro which might be important to inhibit the protein activity. The Rg plot of Mpro is not stable during the simulation peried which shows the potential activity of TQ. Change in protein stability may be due to the alteration of thermodynamic property (Chen and Shen, 2009). This includes protein RMSD, fluctuations and also the protein folding. Destabilization in folding and thermodynamic stability may affect biomolecules function.
The TQ is fitting in the pocket, interacting with active site (C145, H41) (Figure 4), altering the catalytic activity of viral protein. The residues located in the binding pocket and its surrounding (T24, L27, H41, F140, C145, H163, M165, P168, and His172) are imported for a natural compound to interact with. The phytocompound TQ form a hydrogen bond with Glu166, facilitating the pocket shape of the substrate-binding site (Anand et al., 2002) and many hydrophobic interactions with active site (His41, Cys145) and its surrounding residues ( Figure 5).

Conclusion
TQ shows good binding affinity with SARS-CoV-2 NSP5, interacting with active site residues and also with Glue166, maintaining the pocket shape for viral enzymatic activity. This phytomedicine alters the overall thermodynamics properties of SARS-CoV-2 Mpro which may useful for better management of COVID-19 in future. Further experimental validation is required to observe the TQ effect in vivo. The TQ may be used as therapeutic compound against SARS-CoV-2 after experimental confirmation.