Pharmacophore-Based Virtual Screening of Potential SARS-CoV-2 Main Protease Inhibitors from Library of Natural Products

Abstract

Background: The SARS-CoV-2 main protease (M^pro) is an attractive target for drug discovery. Methods: A pharmacophore model was built using the three-dimensional (3D) pharmacophore generation algorithm HypoGen in Discovery Studio 2019. The best pharmacophore model was selected for validation using a test set of 24 compounds and was used as a 3D query for further screening of an in-house database of natural compounds. Lipinski's rule of five was used to assess the drug-like properties of the hit compounds. The filtered compounds were then subjected to bioactivity evaluations. The active compounds were docked into the active site of the SARS-CoV-2 M^pro crystal structure (PDB ID: 7D1M). Results: A suitable 3D pharmacophore model, Hypo1, was found to be the best model, consisting of four features (one hydrophobic feature, one hydrogen bond donor, and two hydrogen bond acceptors). Pharmacophore-based virtual screening with Hypo1 as the query to search an in-house database of 34 439 natural compounds resulted in 1502 hits. Among these, 255 compounds satisfied Lipinski's rule of five. The highest ranking 10 compounds were selected for further experimental testing, and one hit (W-7) illustrated inhibitory activity against SARS-CoV-2 M^pro with an IC₅₀ value of 75 μM. Docking studies revealed that this hit compound retained the necessary interactions within the active site of SARS-CoV-2 M^pro. Conclusion The identified lead natural compound could provide a scaffold for the further development of SARS-CoV-2 M^pro inhibitors.

Keywords

SARS-CoV-2 main protease natural inhibitor pharmacophore virtual screening bioactivity

Introduction

An outbreak of coronavirus disease 2019 (COVID-19), caused by a severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), has become an extraordinary event and a public health threat owing to its global transmission.¹ The SARS-CoV-2 main protease (M^pro; also known as 3-chymotrypsin-like protease) is involved in the processing of viral polyproteins to produce mature active proteins.² The inhibition of SARS-CoV-2 M^pro might be a therapeutic intervention for COVID-19.³ Hence, this protease presents an attractive target for drug-discovery efforts.

A large amount of activity data are available for existing inhibitors of SARS-CoV-2 M^pro that belong to diverse structural classes.^4–14 Several FDA-approved protease inhibitors such as ritonavir and danoprevir were found to inhibit SARS-CoV-2 M^pro at low micromolar concentrations.¹⁵ The pharmacophore models of SARS-CoV-2 M^pro inhibitors based on the structure of ligand-receptor complexes have been extracted.^16–19 However, to the best of our knowledge, to date, there is little information available regarding a highly predictive pharmacophore for SARS-CoV-2 M^pro inhibitors. Therefore, we aimed to generate chemical feature-based pharmacophores for SARS-CoV-2 M^pro inhibitors to guide the discovery of structurally novel inhibitors with improved efficacy.

In this study, we developed chemical feature-based pharmacophores for SARS-CoV-2 M^pro inhibitors. The best-validated pharmacophore model (Hypo1) was used as a query for virtual screening against natural molecules from an in-house database. The hit compounds were then filtered using Lipinski's rule of five as the filtration criteria. The study identified a potential hit molecule through an extensive activity assay in vitro and molecular docking studies.

Results and Discussion

The three-dimensional (3D)-QSAR pharmacophore protocol was used for pharmacophore model generation. The top 10 pharmacophore models were obtained, which contained two hydrogen bond acceptor features (A1 and A2), one hydrogen bond donor feature (D), and one general hydrophobic feature (H). The characteristics of the 10 generated pharmacophores, such as total cost, cost differences, root mean square (RMS) deviation, and correlation coefficient, are listed in Table 1. Cost analysis was used to determine the significance of the pharmacophore models.²⁰ For a significant pharmacophore model, the cost difference must be sufficiently large. Other parameters, such as the correlation coefficient and RMS deviation value, were used to evaluate the ability of a pharmacophore model to predict the activity of training set molecules.²¹ The first model (Hypo1) among the 10 generated pharmacophore models had the highest cost difference (735.01), lowest RMS deviation (1.5433), and best correlation coefficient (0.982). Therefore, Hypo1 was selected as the signified pharmacophore model. Figure 1 shows the 3D spatial arrangement of Hypo1 and its alignment with the most active compound (IC₅₀ = 0.0076 μM) in the training set.

Figure 1.

Best pharmacophore model Hypo 1. (A) Hypo 1 illustrating hydrophobic (H), hydrogen bond donor (D), and hydrogen bond acceptor (A). Distance-1, 10.386 Å; Distance-2, 5.491 Å; Distance-3, 7.343 Å; Distance-4, 5.552 Å; Distance-5, 6.232 Å; and Distance-6 6.538 Å. (B) Hypo 1 mapping with training set compound 1 (IC₅₀ = 0.0076 μM).

Table 1.

Summary of Hypotheses Run.^a

Hypothesis no.	Total cost	Δ Cost	RMS deviation	Correlation (r)
1	108.460	735.01	1.54330	0.981689
2	135.887	707.59	2.17407	0.963293
3	160.114	683.36	2.58981	0.947487
4	169.282	674.19	2.74897	0.940611
5	173.294	670.18	2.8232	0.937242
6	175.300	668.17	2.81124	0.937822
7	175.728	667.74	2.86416	0.935342
8	176.379	667.09	2.73165	0.941515
9	177.810	665.66	2.85273	0.935905
10	181.365	662.11	2.89909	0.933734

Null cost of 10 top-scored hypotheses is 843.473. Fixed cost value is 77.1876. Configuration cost is 19.2586.

The significance of the best pharmacophore model, Hypo1, was evaluated based on its ability to predict the biological activity of the compound. First, Hypo1 showed a reliable ability to predict the activities of the compounds from the training set (Table 2). Of the 24 compounds in the training set, 23 had errors of less than 5%. Second, the predictive ability of Hypo1 was validated using a test set of 24 compounds (Figure 2). Most of the IC₅₀ values were predicted accurately, indicating that Hypo1 has excellent predictive validity. The experimental and estimated activities of the test set are listed in Table 3. Third, the crystal structure of GC-376 (an effective inhibitor of SARS-CoV-2 M^pro) was mapped to the pharmacophore Hypo1. The ligand GC-376 was mapped satisfactorily to the four essential features of Hypo1 with a crystallized conformation (PDB ID: 7D1M) (Figure 3A). The predicted orientation of GC-376 generated by Hypo1 was similar to the experimental binding orientation characterized by X-ray crystallography, with an RMS deviation of 0.8708 (Figure 3B). The ligand GC-376 was used to assess whether Hypo1 could accurately predict the experimental activity of the compound. Our results indicated that the experimental and predicted IC₅₀ values of GC-376 were 0.15 μM and 0.19 μM, respectively, with an error factor of 1.25%. These results further indicate that our pharmacophore model is reliable.

Figure 2.

Molecular structures and activity data (IC₅₀ values, μM) of the 24 test set compounds.

Figure 3.

Validation of the best pharmacophore model Hypo 1 with the crystal structure of compound GC-376. (A) Hypo 1 mapping with the crystal structure of compound GC-376. Pharmacophore model illustrating hydrophobic (H), hydrogen bond donor (D), and hydrogen bond acceptor (A). (B) The conformation of GC-376 proposed by Hypo 1 (C1) mapping with that taken from the crystal structure (C2).

Table 2.

Experimental Activity and Estimated Activity for the Training Set Based on the Pharmacophore Model Hypo 1.

Compound	Measured IC₅₀ (μM)	Estimated IC₅₀ (μM)	Error factor	Ref.
1	0.0076	0.0058	−1.31	⁹
2	0.012	0.0186	1.55	⁹
3	0.036	0.032	−1.11	⁹
4	0.04	0.064	1.60	⁷
5	0.15	0.13	−1.15	⁴
6	0.15	0.19	1.25	¹⁴
7	0.18	0.349	1.90	⁵
8	0.40	0.80	1.99	¹⁰
9	0.55	0.67	1.22	¹²
10	0.67	0.88	1.31	⁵
11	1.55	4.21	2.72	⁶
12	2.39	2.60	1.09	⁵
13	7.29	5.73	−1.27	¹²
14	9.45	4.77	−1.98	¹²
15	10.69	6.61	−1.62	¹³
16	15.75	2.46	−6.40	⁶
17	20.31	7.50	−2.71	¹²
18	20.96	20.33	−1.03	¹²
19	23.50	31.70	1.35	⁸
20	100.70	61.87	−1.63	⁸
21	124.93	223.97	1.79	⁴
22	349	368.88	1.05	¹¹
23	370	922.75	2.49	¹¹
24	635	378.54	−1.68	¹¹

Table 3.

Experimental Activity and Estimated Activity for the Test Set Based on the Pharmacophore Model Hypo 1.

Compound	Measured IC₅₀ (μM)	Estimated IC₅₀ (μM)	Error factor	Ref.
1	0.0076	0.0077	1.02	⁶
2	0.017	0.016	−1.09	⁶
3	0.052	0.052	1.00	⁶
4	0.069	0.07	1.02	⁶
5	0.093	0.09	−1.04	⁶
6	0.29	0.26	−1.13	⁶
7	0.45	0.58	1.29	⁵
8	0.75	0.52	−1.44	⁶
9	0.97	0.29	−3.30	⁵
10	1.82	1.83	1.00	³
11	2.78	0.83	−3.35	⁹
12	3.9	6.98	1.79	¹⁰
13	6.4	2.54	−2.51	⁵
14	7.16	368.91	51.52	⁹
15	42	368.74	8.78	⁸
16	43.13	44.24	1.02	⁹
17	68	7.84	−8.67	⁸
18	76	75.53	−1.01	⁸
19	85	99.95	1.17	⁸
20	111	111.45	1.00	⁸
21	234	73.93	−3.16	⁸
22	438	110.04	−3.98	⁸
23	497	499.17	1.00	⁸
24	3047	378.54	−1.68	⁸

We used Hypo1 as a 3D query to screen an in-house library of natural compounds, and 1502 hits were captured. To assess drug-like and pharmacokinetic properties, these hit compounds were further subjected to Lipinski's rule of five. The results showed that 255 compounds met the standards. Subsequently, the 10 compounds with the highest ranking were selected for determination of their inhibitory activity against SARS-CoV-2 M^pro in vitro using a fluorescence resonance energy transfer protease assay. The results showed that one compound (W-7) demonstrated a significant IC₅₀ value of 75 μM.

Docking was performed using the CDOCKER algorithm to predict the binding mode between the hit compound W-7 and SARS-CoV-2 M^pro. The co-crystallized compound GC-376 was first re-docked into the active site of SARS-CoV-2 M^pro. The results showed that CDOCKER was able to correctly dock compound GC-376 to the active center of SARS-CoV-2 M^pro with an RMS deviation of 0.8708. Subsequently, the hit compound W-7 was docked into the active site of SARS-CoV-2 M^pro. As shown in Figure 4A, compound W-7 was found to dock well and was located at the center of the binding pocket. The docking results (Figure 4B) also indicated that two of the hydroxyl groups of the hit compound and Glu 166 of SARS-CoV-2 M^pro formed optimal hydrogen bonds. In addition, a hydrogen bond was formed between one of the hydroxyl oxygens of the hit compound and Gln 189. Furthermore, compound W-7 satisfied the expected hydrophobic interaction, as defined by our query Hypo1, with the surrounding residues as His41, Met49, Phe140, Ser144, and Met165. This binding motif is in agreement with the results of a previous study.²² The alignment of Hypo1 with compound W-7 was performed. The results indicated that W-7 could map Hypo1 well (Figure 4C). The RMS deviation between the orientation of compound W-7 mapped with Hypo1 and its orientation after docking was 0.807 (Figure 4D). These results clearly confirmed that the binding mode between SARS-CoV-2 M^pro and W-7 is consistent with that proposed by Hypo1.

Figure 4.

The hit compound W-7. (A) A 3D structure of W-7 in SARS-CoV-2 M^pro binding groove (PDB ID: 7D1M). SARS-CoV-2 M^pro is represented as a molecular surface colored by electrostatic potential. The negative, neutral, and positive regions are represented by opaque red, white, and blue, respectively. (B) Model of W-7 bound to SARS-CoV-2 M^pro (PDB ID: 7D1M). Hydrogen bonds shown in dotted lines. (C) Hypo 1 mapping with W-7. Pharmacophore model illustrating hydrophobic (H), hydrogen bond donor (D), and hydrogen bond acceptor (A). (D) Conformation of W-7 proposed by Hypo 1 (C1) mapping with that created by the docking (C2).

Conclusions

In summary, pharmacophore models for SARS-CoV-2 M^pro inhibitors were constructed based on 24 diverse compounds with varying levels of activity. The best pharmacophore (Hypo1) was selected based on the highest cost difference and correlation coefficient, and the lowest total cost and RMS value. The validated Hypo1 was then employed as a query for virtual screening of molecules from an in-house natural product database, and 1502 hits were identified. The hit compounds were further analyzed using drug-like filters and SARS-CoV-2 M^pro-inhibitory activity analysis. Finally, one hit (W-7) showed effective biological activity, with an IC₅₀ value of 75 μM. Molecular docking was performed to understand the interaction between this hit compound and SARS-CoV-2 M^pro. Taken together, these results show that the hit W-7 may act as a good lead against SARS-CoV-2 M^pro.

Materials and Methods

Generation of Pharmacophore Models

A set of 24 structurally diverse compounds with IC₅₀ values ranging from 0.0076 μM to 635 μM was used as the training set compounds (Figure 5). The structures of the selected molecular datasets were built using Discovery Studio 2019 (DS). Every compound was minimized to a local energy minimum using the CHARMM force field.²³ A maximum of 250 conformers within a threshold of 20 kcal/mol above the global energy minimum was generated for each compound using the Poling algorithm.²⁴ Three-dimensional QSAR pharmacophore methodologies in DS were used to create hypotheses. The chemical features used to construct the hypotheses were the hydrogen bond acceptor (A), hydrogen bond donor (D), and general hydrophobic features (H).¹⁸

Figure 5.

Molecular structures and activity data (IC₅₀ values, μM) of the 24 training set compounds.

Database Search

The validated pharmacophore Hypo1 was utilized as a 3D query to search an in-house library containing 34 439 natural compounds to obtain hit molecules that matched the chemical features of the queried pharmacophore. The Search 3D database protocol implemented within DS was employed with a best/flexible search approach. All hit compounds were further evaluated using Lipinski's rule of five.

Biological Testing

A fluorescence resonance energy transfer protease assay (Beyotime Biotechnology) was used to determine the inhibitory activity of the compounds against SARS-CoV-2 M^pro. The recombinant SARS-CoV-2 M^pro was mixed with an assay buffer, serial dilutions of each compound, or the positive control inhibitor ebselen with a final volume of 96 μL, and incubated. The reaction was initiated by the addition of 4 μL fluorogenic substrate. Subsequently, the fluorescence signal at 325 nm (excitation) and 393 nm (emission) was measured every 30 s for 10 min using a SpectraMax i3x multifunctional enzyme marker (Molecular Devices). The relative fluorescence units of compounds at various concentrations compared to those added with DMSO were used to generate IC₅₀ curves.

Molecular Docking

The molecules were subjected to docking in the SARS-CoV-2 M^pro active site using CDOCKER within DS.²⁵ The 3D structure of SARS-CoV-2 M^pro (PDB ID: 7D1M) was obtained from the Protein Data Bank (http://www.rcsb.org/pdb/). The active site was defined as a sphere of 12 Å around the geometric centroid of the co-crystal ligand GC-376. Each ligand was allowed to generate 10 different orientations, and the most desirable orientation was selected based on the scoring function (-CDOCKER interaction energy) and binding interaction.

Footnotes

Acknowledgments

We would like to thank Editage () for English language editing.

Authors’ Notes

Jing Wang, Yu Jiang, Zhanli Wang, and Yuheng Ma contributed equally to this work.

Authors’ Contributions

JW and YJ conducted the study, performed data analysis, and prepared the manuscript. YW and HY performed experiments. ZW and YM designed and supported the study, and edited the manuscript. All authors approved the final version of the manuscript.

Declaration of Conflicting Interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article. This study was supported by the National Natural Science Foundation of China (grant number 82170296), and Scientific Research Fund Project of Baotou Medical College (grant number BYJJ-XG 202001).

ORCID iD

MiZhanli Wang

References

Banerjee

Perera

Tillekeratne

LMV

. Potential SARS-CoV-2 main protease inhibitors. Drug Discov Today. 2021;26(3):804-816. doi: 10.1016/j.drudis.2020.12.005

Zhu

Shyr

Zheng

. Viral proteases as targets for coronavirus disease 2019 drug development. J Pharmacol Exp Ther. 2021;378(2):166-172. doi: 10.1124/jpet.121.000688

Ullrich

Nitsche

. The SARS-CoV-2 main protease as drug target. Bioorg Med Chem Lett. 2020;30(17):127377. doi: 10.1016/j.bmcl.2020.127377

Faheem, Kumar

Sekhar

KVGC

, et al. Druggable targets of SARS-CoV-2 and treatment opportunities for COVID-19. Bioorg Chem. 2020;104:104269. doi: 10.1016/j.bioorg.2020.104269

Zhang

Lin

Sun

, et al. Crystal structure of SARS-CoV-2 main protease provides a basis for design of improved α-ketoamide inhibitors. Science. 2020;368(6489):409-412. doi: 10.1126/science.abb3405

Jin

, et al. Structure of M^pro from SARS-CoV-2 and discovery of its inhibitors. Nature. 2020;582(7811):289-293. doi: 10.1038/s41586-020-2223-y

Dai

Zhang

Jiang

, et al. Structure-based design of antiviral drug candidates targeting the SARS-CoV-2 main protease. Science. 2020;368(6497):1331-1335. doi: 10.1126/science.abb4489

Kang

. Progress in developing inhibitors of SARS-CoV-2 3C-like protease. Microorganisms. 2020;8(8):1250. doi: 10.3390/microorganisms8081250

Gimeno

Mestres-Truyol

Ojeda-Montes

, et al. Prediction of novel inhibitors of the main protease (M-pro) of SARS-CoV-2 through consensus docking and drug reposition. Int J Mol Sci. 2020;21(11):3793. doi: 10.3390/ijms21113793

10.

Sacco

Hurst

, et al. Boceprevir, GC-376, and calpain inhibitors II, XII inhibit SARS-CoV-2 viral replication by targeting the viral main protease. Cell Res. 2020;30(8):678-692. doi: 10.1038/s41422-020-0356-z

11.

Amin

Banerjee

Ghosh

Gayen

Jha

. Protease targeted COVID-19 drug discovery and its challenges: insight into viral main protease (M^pro) and papain-like protease (PLpro) inhibitors. Bioorg Med Chem. 2021;29:115860. doi: 10.1016/j.bmc.2020.115860

12.

Qiao

Zeng

, et al. SARS-CoV-2 M^pro inhibitors with antiviral activity in a transgenic mouse model. Science. 2021;371(6536):1374-1378. doi: 10.1126/science.abf1611

13.

Vatansever

Yang

Drelich

, et al. Bepridil is potent against SARS-CoV-2 in vitro. Proc Natl Acad Sci U S A. 2021;118(10):e2012201118. doi: 10.1073/pnas.2012201118

14.

Bafna

White

Harish

, et al. Hepatitis C virus drugs that inhibit SARS-CoV-2 papain-like protease synergize with remdesivir to suppress viral replication in cell culture. Cell Rep. 2021;35(7):109133. doi: 10.1016/j.celrep.2021.109133

15.

Citarella

Scala

Piperno

Micale

. SARS-CoV-2 M^pro: a potential target for peptidomimetics and small-molecule inhibitors. Biomolecules. 2021;11(4):607. doi: 10.3390/biom11040607

16.

Mosayebnia

Hajiagha Bozorgi

Rezaeianpour

Kobarfard

. In silico prediction of SARS-CoV-2 main protease and polymerase inhibitors: 3D-pharmacophore modelling. J Biomol Struct Dyn. 2022;40(14):6569-6586. doi: 10.1080/07391102.2021.1886991

17.

Pathak

Chen

Hsu

, et al. Uncovering flexible active site conformations of SARS-CoV-2 3CL proteases through protease pharmacophore clusters and COVID-19 drug repurposing. ACS Nano. 2021;15(1):857-872. doi: 10.1021/acsnano.0c07383

18.

Halimi

Bararpour

. Natural inhibitors of SARS-CoV-2 main protease: structure based pharmacophore modeling, molecular docking and molecular dynamic simulation studies. J Mol Model. 2022;28(9):279. doi: 10.1007/s00894-022-05286-6

19.

Vázquez-Mendoza

Mendoza-Figueroa

García-Vázquez

Correa-Basurto

García-Machorro

. In silico drug repositioning to target the SARS-CoV-2 main protease as covalent inhibitors employing a combined structure-based virtual screening strategy of pharmacophore models and covalent docking. Int J Mol Sci. 2022;23(7):3987. doi: 10.3390/ijms23073987

20.

John

Thangapandian

Arooj

Hong

Kim

Lee

. Development, evaluation and application of 3D QSAR pharmacophore model in the discovery of potential human renin inhibitors. BMC Bioinformatics. 2011;12(Suppl 14):S4. doi: 10.1186/1471-2105-12-S14-S4

21.

Huang

Wang

Qiang

, et al. The discovery of novel BCR-ABL tyrosine kinase inhibitors using a pharmacophore modeling and virtual screening approach. Front Cell Dev Biol. 2021;9:649434. doi: 10.3389/fcell.2021.649434

22.

Feng

, et al. Both boceprevir and GC376 efficaciously inhibit SARS-CoV-2 by targeting its main protease. Nat Commun. 2020;11(1):4417. doi: 10.1038/s41467-020-18233-x

23.

Vanommeslaeghe

Hatcher

Acharya

, et al. CHARMM General force field: a force field for drug-like molecules compatible with the CHARMM all-atom additive biological force fields. J Comput Chem. 2010;31(4):671-690. doi: 10.1002/jcc.21367

24.

Wang

Zhang

Huang

. Pharmacophore modeling and in silico screening for new KDR kinase inhibitors. Bioorg Med Chem Lett. 2007;17(8):2126-2133. doi: 10.1016/j.bmcl.2007.01.089

25.

Lipinski

. Lead- and drug-like compounds: the rule-of-five revolution. Drug Discov Today Technol. 2004;1(4):337-341. doi: 10.1016/j.ddtec.2004