{'fraction_non_alphanumeric': 0.03519829738468465, 'fraction_numerical': 0.03413416281258953, 'mean_word_length': 4.870735223450264, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Fisetin 8-C-glucoside as entry inhibitor in SARS CoV-2 infection: molecular modelling study

Abstract Coronaviruses are RNA viruses that infect varied species including humans. TMPRSS2 is gateway for SARS CoV-2 entry into the host cell. It causes proteolytic activation of spike protein and discharge of the peptide into host cell. The TMPRSS2 inhibition could be one of the approaches to stop the viral entry, therefore, interaction pattern and binding energies for Fisetin and TMPRSS2 have been explored in the present study. TMPRSS2 peptide was used for homology modelling and then for further study. Molecular docking score and MMGBSA Binding energy of Fisetin was better than Nafamostat, a known inhibitor of TMPRSS2. Post docking MM-GBSA free energy for Fisetin and Nafamostat was −42.78 and −21.11 kcal/mol, respectively. Fisetin forms H bond with Val 25, His 41, Lys 42, Lys 45, Glu 44, Ser186. Nafamostat formed H bonds with Lys 85, Asp 90, Asp 203. RMSD plots of TMPRSS2, TMPRSS2-Fisetin and TMPRSS2-Nafamostat complex showed stable profile with very small fluctuation during entire simulation of 150 ns. Significant decrease in TMPRSS2-Fisetin and TMPRSS2-Nafamostat complex fluctuation occurred around His 41, Glu 44, Gly 136, Ser 186 in RMSF study. During simulation Fisetin interaction was observed with residues Val 25, His 41, Glu 44, Lys 45, Lys 87, Gly 136, Gln 183, Ser 186 likewise interaction of Nafamostat with Lys 85, Asp 90, Asn 163, Asp 203 and Ser 205. Post simulation MM-GBSA free energy was found to be −51.87 ± 4.3 and −48.23 ± 4.39 kcal/mol for TMPRSS2 with Fisetin and Nafamostat, respectively. Communicated by Ramaswamy H. Sarma


Introduction
Recent outbreak of Corona virus SARS CoV-2 infection has created a serious health problem to the whole world. There is urgent need to understand the exact mechanism of action of virus on human system and to develop a specific drug as no specific line of treatment is currently available. SARS-CoV-2 is a single-stranded RNA-enveloped virus, its spike (S) protein binds to the host angiotensin-converting enzyme 2 (ACE2) and primed by trans membrane serine protease, TMPRSS2 (Transmembrane serine protease 2), that leads to entry of virus into the cell. TMPRSS2 plays a role in the proteolytic activation of viral S protein and helps in respiratory disease caused by influenza (B€ ottcher et al., 2006), SARS-CoV and MERS (Iwata-Yoshikawa et al., 2019) viruses. TMPRSS2 may also regulate SARS-CoV-2 assembly and release from the cell membrane (Alireza et al., 2020;Shen et al., 2017). After entry into the cells there is synthesis of polyproteins and RNA with the help of replicase-transcriptase and RNAdependent RNA polymerase. These proteins like host TMPRSS2, ACE2 and viral proteins like M pro Protease and RNA-dependent RNA polymerase can be a potential drug targets for SARS-CoV-2 (Glowacka et al., 2011;Hoffmann et al., 2020;Shen et al., 2017;Shulla et al., 2011;Singh & Mishra, 2020;Sinha et al., 2020). TMPRSS2 consists of 492 amino acids which is attached to the cell membrane. The mature proteases after autocatalytic breakdown between Arg 255 and Ile 256 are generally membrane-bound. TMPRSS2 catalytic domain contains a catalytic triad consisting of the amino acid residues His 296, Asp 345 and Ser 441 (Lin et al., 1999;Park, 2010). The exact physiological function of TMPRSS2 in vivo remains unknown. TMPRSS2 protein is highly expressed in bronchial epithelial cells than in type II pneumocytes (Bertram et al., 2012).
There are mainly three classes of antiviral drugs under clinical trial for the effectiveness and safety against COVID-19 infection. First is RNA dependent RNA Polymerase Inhibitors (Remdesivir, Favipiravir, Ribavirin) and second is Protease Inhibitors (Lopinavir, Ritonavir) and various miscellaneous drugs like Interferons, Anti-cytokines, hydroxychloroquine comes in third group. Bromhexine is found to be an inhibitor of the TMPRSS2 with an IC50 of 0.75 mM (Lucas et al., 2014). The Camostat, Nafamostat, Aprotinine and bromhexine although not developed specifically for targeting TMPRSS2 were reported to inhibit TMPRSS2 protease activity (Lucas et al., 2014). Secondary metabolites, phytochemicals have therapeutic potential which has been explored in the present study against TMPRSS2 [Transmembrane serine protease 2] to prevent SARS CoV-2 infection. Ligand we selected for present study is a flavanol Fisetin 8-C-glucoside (Fisetin), generally present in fruits and vegetables such as apple, onion, strawberry, grape, cucumber etc. and are known to possess anticancer, antioxidant activities (Khan et al., 2013). Present study was initiated to understand how viral entry into target cells can be blocked by use of phytochemicals. Nafamostat, protease inhibitor, was taken for comparative study.

Method
Homology modelling method used for building a threedimensional model of the TMPRSS2 by comparing its amino acid sequence with an experimental model of related homologous template protein. The primary sequence of the target human TMPRSS2 was obtained from UniProtKB database with a sequence id: O15393, entry name TMPS2_HUMAN sequence length 492aa (UniProt). The amino acid sequence of TMPRSS2 subjected to BLASTp program with default parameters in Uniprot to find out template for sequence alignment (Altschul et al., 1990;Jakub ık et al., 2013). Peptide chain which consists of catalytic domain, from 256 to 492, was used for modelling. So, in our model Amino acid 1 corresponds to amino acid 256. The model was generated using Prime in Schr€ odinger Suite and loop refinement was done. To compare the model iTASSER model prediction tool was also used for model preparation. The validation of Target model obtained from Prime was performed by analysing the psi (w) and phi (u) torsion angles using Ramachandran plot. The 3D structure of target model obtained from homology modelling was optimized by hydrogen bond assignments (Ramachandran et al., 2016). PROCHECK was used to assess the quality of the model structure (Laskowski et al., 1993). The modelled structure was further minimized using 10 ns MD simulation using OPLS3 force field as per default protocol.

Active site prediction
Schr€ odinger suite was used to predict receptor sites of TMPRSS2 structure by using Sitemap module. It had generated potential site in the receptors with quantitative sitescore values for each site.

Molecular docking
Protein preparation wizard of Schr€ odinger suite was used to optimize TMPRSS2. Protonation of the TMPRSS2 at the physiological pH performed by Epik module of Schr€ odinger suite. Water molecules greater than 3 Å away from the Fisetin (Fisetin 8-C-glucoside:Pub chem i.d. 5281615) were deleted before subjecting to minimization and hydrogen bond optimization. Impact refinement module utilizing the OPLS3e force field for restrained molecular minimization to make relaxed structure after deleting all water molecules, heteroatom and adding hydrogen atoms. OPLS3e force field used in minimization studies to relieve steric clashes. LigPrep module was used for ligands preparation. Energy minimisation (in gas phase) was done using Macromodel with the OPLS3e force field. Glide software was used to carry out molecular docking and grid box was generated around catalytic triad His 41, Asp 90 and Serine 186 residues as centroid (Halgren, 2007). The docking was carried out using extra precision (XP) mode with default parameters to analyze the binding modes of compounds under study (Glide, 2015). The molecular docking generated H-bond between ligand and target, and the consequent glide energy. We have used 250 phytochemicals compiled from published literature with antimicrobial properties. Compounds were ranked based on GlideScore and the best pose of the ligand and chosen for further study (Friesner et al., , 2006Halgren et al., 2004;Jacobson et al., 2004;Kumar et al., 2016;Sastry et al., 2013).

Molecular dynamics simulation
Desmond module of Schrodinger Suite used for Molecular dynamics (MD) of ligand protein complex. The Fisetin-TMPRSS2 complex was solvated using SPC aqueous solvation of orthorhombic solvent box. Solvent buffer extended 10 Å beyond the TMPRSS2 in all directions and OPLS3 force field was applied. System was neutralized by adding the required number of counter ions. Initial two stages with steepest descent minimization comprised of 2000 steps with and without a restraint of 50 kcal/mol/Å 2 on the solute atoms carried out for relaxation protocol before production runs (Harder et al., 2016). The MD simulations done for short four times with 12 ps in NVT ensemble at 10 K with solute heavy atoms restrained with a force constant of 50 kcal/mol/Å 2 . NPT ensemble at10 K with the same restraint used for 12 ps simulation further NPT ensemble at 300 K with the same restraint for12 ps and NPT ensemble at 300 K without any restraints used finally for 24 ps simulation. Simulation with NPT ensemble at 300 K temperature and 1.01325 bar pressure run under unrestrained condition (Nos e-Hoover chain thermostat relaxation time of 1 ps and the isotropic Martyna-Tobias-Klein barostat with a relaxation time of 2 ps). Short range interactions studied with 9 Å cut off and long-range coulombic interactions analysed with smooth particle mesh Ewald method (PME). r-RESPA integrator used for non-bonded forces evaluation. After every three steps long range forces and at every step short-range force were recorded. MD simulations run for 150 ns. At every 150 ps trajectories were recorded and energy recording interval was maintained at 1.2 ps. OPLS3 (Optimized Kanhesia for Liquid Simulations) force fields used for analysing RMSD and RMSF of the complex structure. At every 150.0 ps intervals trajectories were recorded for evaluation. Molecular dynamic simulation study was done with Nafamostat and TMPRSS2 complex by the same protocol (MacroModel, 2015;Singh & Mishra, 2020).
Binding free energy was studied with Prime-MM/GBSA (Molecular Mechanics/Generalized Born Surface Area) (Mobley & Dill, 2009;Prime, 2015). Surface GB model with Gaussian surface in place of van der Waals surface used for better depiction of the solvent accessible surface area (Jakub ık et al., 2013). Equations used for binding energy (DG bind ) calculations were as follows: where E complex , E protein and E ligand representing minimized energies of the complex, TMPRSS2 and Fisetin, respectively (Das et al. 2009).
where solvation free energies were G solv-complex , G solv-protein and G solv-ligand for complex, TMPRSS2, and Fisetin, respectively.
The surface area energies were calculated for the complex, TMPRSS2, and Fisetin, by G SA-complex , G SA-protein and G SAligand , respectively; The scoring and interaction parameters of XP docking was used as selection parameter for selection of compound for MD simulation study (Kevin et al., 2006;Shivakumar et al., 2010;Wang et al., 2015). Prime program of Schr€ odinger software suite employed for all docking poses which further rescored with the MM/GBSA approach to allow more flexibility to protein and assess ligand affinity thereby providing more reliable results. The variable dielectric solvent model VSGB 2.0 was used for the study which comprises empirical corrections for modelling directionality of hydrogen-bond and p-stacking interactions. This method gives significant binding free energies for a wide range of protein-ligand complexes (Li et al., 2011;Mulakala & Viswanadhan, 2013). MM minimisation of the complex was carried out by allowing residues within 5.0 Å of the ligand to relax and keeping the rest of the structure fixed (Srivastava & Sastry, 2012).

Results and discussion
SARS CoV-2 entry into host cells involves two important steps like binding to host cell receptor and entry of viral genetic material into the host cell. Viral entry requires a main protease (TMPRSS2) which causes proteolytic processing of viral Spike (S) protein. The TMPRSS2 could be an important target to prevent viral entry into host cells. Hoffman et al. (2020) in an in vitro study showed that Camostat mesylate, an inhibitor of TMPRSS2, partially inhibit the entry of SARS-CoV-2 into primary pulmonary cells and cell lines. In another study it was reported that Camostat blocks the entry of SARS-CoV in human Calu-3 airway epithelial cells (Kawase et al., 2012). Serine protease inhibitor such as Nafamostat, inhibits TMPRSS2 and resists S-protein facilitated membrane fusion of MERS-CoV (Yamamoto et al., 2016). Bromhexine also reported to be an inhibitor of the TMPRSS2 with an IC-50 of 0.75 mM (Lucas et al., 2014). Inhibition of TMPRSS2 by Nafamostat was reported at 10-fold lower concentration than Camostat (Sielaff et al., 2011;Hoffman et al., 2020). There are reports on several other small molecule which inhibits TMPRSS2 with low nanomolar affinity (Yang & Honig, 2000). We tried to explore interaction pattern of some phytochemicals against TMPRSS2 and results were compared with Nafamostat.

Homology modelling (HM)
The three-dimensional structure of TMPRSS2 was not available in protein data bank (PDB) so we tried HM for Model preparation. With a similarity of 36.9% the sequence of Oryzias latipes enteropeptidase light chain (PDB ID: 3w94) has shown as best homology in the BLAST search for Human TMPRSS2 amino acid sequence. If the sequence similarity between two proteins is higher than 25% then their structures are considered to be similar (Rost, 1999;Yang & Honig, 2000). Therefore, the 3D structure of Human TMPRSS 2 built through homology modelling using the 3D structure of Oryzias latipes enteropeptidase light chain as the template for modelling TMPRSS2 protein. Model accuracy is related to the degree of sequence identity and similarity between template and TMPRSS2 (target). The crystal structure of Oryzias latipes enteropeptidase light chain was retrieved from PDB (PDB ID: 3w94) (Xu et al., 2014). PROCHECK analysis in Figure  1(a), showed that 81.6% of the residues are in the most favoured regions, 16.3% in additional allowed regions, 2% residues in generously allowed regions and 0.0% of the residue were obtained in disallowed regions (Laskowski et al., 1993). PROCHECK analysis of the model prepared by i-Tasser showed that 76.4% of the residues are in the most favoured regions, 16.6% in additional allowed regions, 3.5% residues in generously allowed regions and 3.5% of the residue were obtained in disallowed regions (Figure 1b).

Active site prediction
Active site residues of TMPRSS2 were predicted by Site Map. It identified the potential sites with residues with a site score of 1.06 considering parameters which include volume, size, contact, amino acid exposure, hydrophilicity, hydrophobicity and donor/acceptor ratio. The sites with site score of 1 and above are consider as suitable site for the ligand binding (Xie .& Huang 2012). The site included His 41 and Asp 90 also. The protease catalytic domain of TMPRSS2 contains a catalytic triad consisting of the amino acid residues His 296, Asp 345 and Ser 441 (Lin et al., 1999;Park, 2010

Molecular docking
An attempt was made in the present study to find inhibitor of TMPRSS2 through molecular docking and simulations. The receptor grid was formed around His 41, Asp 90 and Serine 186 in homology model of TMPRSS2 to conduct molecular docking. 250 phytochemicals with antimicrobial properties were chosen for docking against TMPRSS2. The phytochemicals were analysed by glide score, interaction with particular residue, His 41, Asp 90 and Ser 186 which forms catalytic triad of TMPRSS2 and Post docking MM-GBSA binding energies. We used a known inhibitor Nafamostat for comparative study. Out of 250 compound we identified 10 best phytochemicals from Extra precision (XP) docking and MMGBSA study (Table 1). Fisetin showed a good interaction at catalytic domain by forming hydrogen bonds and other non-bonded interactions with TMPRSS2 and lowest MM-GBSA binding energy. The Hydrogen bond interaction occurred at His 41, Ser 186 and with several other amino acids. His 41, Ser 186 are part of catalytic triad of TMPRSS2. Hydrophobic interactions were also seen with Cys 42, Tyr 82, Trp 206 (Figure 2a). The mode of binding and interaction of residues were similar to Nafamostat (Figure 2b). The Fisetin had produced much better MM-GBSA binding energies than Nafamostat. The post docking MM-GBSA binding free energy with TMPRSS2 was found to be for Fisetin and Nafamostat as À42.78 and À21.11 (kcal/mol), respectively (Table 1). Considering all factors, Fisetin was our choice phytochemical for molecular dynamic simulation study.

Molecular dynamics simulation
MD simulations of the TMPRSS2 complexed with Fisetin were performed using Desmond. It evaluates the dynamic interactions and analyse structural stability and binding site adaptations to the docked Fisetin during 150 ns simulation. MD simulations were run for the TMPRSS2 and Fisetin for 150 ns. The stability and fluctuations of the TMPRSS2 and Fisetin complex was investigated using carbon alpha (Ca) atoms and ligand during 150 ns simulation. RMSD is a measure of protein conformational stability. It can be observed that how fast and how far the protein deviates from the starting structure for a more quantitative analysis of the protein fluctuations. All protein frames are first aligned on the backbone of the reference frame, and then the RMSD was calculated on the basis of the selection of the atom. RMSD plot of the complex showed some fluctuations in the beginning then attained equilibrium at 25 ns and remained stable for next 125 ns, i.e. fluctuation remains within 1 Å. Changes of the order of 1-3 Å are considered satisfactory for small, globular proteins ( Figure 3).
Structural stability of protein gives an indication of the shape of the molecule at each time during entire simulation in radius of gyration plot (Rg). It was calculated to analyse structural changes of TMPRSS2, when the Fisetin was bound. The plot of radius of gyration in simulation time in different condition of the TMPRSS2 has given in Figure 4. TMPRSS2-Fisetin complex during 150 ns molecular dynamics simulation did not bring significant conformational changes in its structure. Figure 5(a) and 5(b) showed H-bond formation between TMPRSS2-Fisetin and TMPRSS2-Nafamostat, respectively, during 150 ns MD simulation study. The average H-bond between TMPRSS2 and Fisetin was 3.04 while 2.67 found with TMPRSS2 and Nafamostat. Upon general hydrogen bond analysis throughout 150 ns simulation, Fisetin was found to form more than 76% hydrogen bonds with Glu44, Lys45 and Lys87. The average distance was also less than 3 Å signifies moderate type of binding ( Table 2).
The residues involved in the interactions with a Fisetin during simulation were depicted in Root mean square fluctuations (RMSF) (Figure 6). Root mean square fluctuations  (RMSF) was calculated to identify the residues involved in the key interactions with a ligand. Significant movement was observed around Amino acid His 41, Glu44, Lys45, Pro46, Leu47, Lys87, Gly 136 and Ser 186 which are critical residues in bonds formation. These amino acids might be playing significant role in complex formation ( Figure 6). The secondary structures like alpha helices and beta strands displayed less fluctuations as compared to loop regions. The essential hydrophobic residues in the protein binding site for the ligand showed rigid behaviour with very less amounts of fluctuations, which reflects the potential of ligand for forming stable interactions with protein during 150 ns simulation ( Figure 6). Figure 7 shows H bond, hydrophobic, ionic and water bridges formation between TMPRSS2 and Fisetin during entire period of simulation. Top panel indicates interaction that TMPRSS2 makes with the Fisetin in entirety. Interaction between specific residue and ligand in each   trajectory frame was given in bottom panel. Dark orange shade indicates formation of more than one specific contacts (Figure 7). Figures 8a, 8b, (Figure 8a and 8b). It was also evident that Nafamostat was making    similar pattern of hydrogen bonds, hydrophobic interactions and water bridges with TMPRSS2 (Figure 9a and 9b). Free energy determines most of the molecular processes, like molecular interaction, chemical reaction, protein folding and others. The enthalpy component influences strength of the interactions in terms of free energy. These interactions include hydrogen, van der Waals bonding, ionic, electrostatic and polarisation of the interacting groups. The entropy term of the binding free energy greatly influenced by solvation effects (due to release of bound water molecules on Fisetin binding to TMPRSS2). Molecular dynamics simulation (MD) helps to calculate free energies of molecular interactions by using theoretical calculations (Christ et al. 2009). In binding free energy estimation, docking and scoring are generally efficient but not very accurate. These methods can differentiate between binders and non-binders to certain extent  Schematic representation of TMPRSS2 and Fisetin interaction during 150 ns simulation (green-hydrophobic, blue-polar, red charged-salt bridge; amino acid 1 will correspond to amino acid 256). (Gohlke & Klebe 2002). Alchemical perturbation (AP) method is another alternative which is accurate but requires large data sets and computationally intensive. Therefore, not often used in drug design (Homeyer et al. 2014;Shirts et al. 2007). There are third group of methods with intermediate performance which is again statistical method but only of the end states. These end point methods are MM/PBSA (molecular mechanics [MM] with Poisson-Boltzmann [PB] and surface area solvation) and GBSA, which is inexpensive than AP, more precise than the scoring functions (Kollman et al., 2000). MM/PB(GB)SA have been extensively studied for free energy calculations (Wang et al., 2001, Gohlke et al. 2003, Hou et al., 2008. where DG is binding free energy, DH enthalpy, DS entropy, and T represent temperature. Surface tension of the solvent related with c a coefficient and b is fitting parameter. DG sol, DE MM and -TDS, are the variations of the solvation free energy, gas phase energy and the conformational entropy, respectively. DE MM calculate DE vdw (van der Waals) energies and DE internal (bond, angle, and dihedral energies) and DE electrostatic (electrostatic). DG solv is the sum of the nonelectrostatic solvation component (nonpolar contribution), DG SA electrostatic solvation energy (polar contribution). DG PB/GB , GB or PB model were used for polar contributions. Solvent accessible surface area (SASA) is used for nonpolar energy contribution. Conformational snapshots of MD simulations were used for calculation of -TDS. Post simulation MM-GBSA binding energy estimation for 1000 snapshots was performed for the Complex. The binding energy of Fisetin to TMPRSS2 was found to be À51.87 ± 4.3 kcal/mol and range as À67.22 to À34.25 kcal/mol. The binding energy of Nafamostat to TMPRSS2 was observed as À48.23 ± 4.39 kcal/ mol and range as À61.46 to À30.69 kcal/mol (Table 3). In an another study on TMPRSS2 with Camostat mesylate interaction the average MMGBSA free binding energy was reported as À65.20 ± 1.64 kcal/mol in 100 snapshots (Kumar et al., 2020).

Conclusion
SARS CoV-2 entry into host cells involves two important steps like binding to host cell receptor and release of virus genome into the host cell's cytoplasm. Viral entry requires a main protease (TMPRSS2) which causes proteolytic processing of viral Spike (S) protein. In the present study we demonstrated that Fisetin 8-C-glucoside stably interact with TMPRESS2. MMGBSA binding energy after docking for Fisetin 8-C-glucoside and Nafamostat was À42.78 and À21.11 kcal/ mol, respectively. Post simulation MMGBSA for Fisetin and Nafamostat for 1000 snapshots was À51.87 ± 4.3 and À48.23 ± 4.39 kcal/mol. Molecular docking and simulation studies revealed that Fisetin 8-C-glucoside TMPRSS2 protein complex showed better hydrogen bonds, hydrophobic interactions than known inhibitor Nafamostat. The findings suggest a probable role of Fisetin as potential lead compound in drug development against host protease TMPRSS2. Thus, the Fisetin 8-C-glucoside can be helpful in the management of COVID-19.

Disclosure statement
No potential conflict of interest was reported by the authors.

Funding
Banaras Hindu University.