key: cord-0709365-tym691he authors: Zev, Shani; Raz, Keren; Schwartz, Renana; Tarabeh, Reem; Gupta, Prashant Kumar; Major, Dan T. title: Benchmarking the Ability of Common Docking Programs to Correctly Reproduce and Score Binding Modes in SARS-CoV-2 Protease Mpro date: 2021-05-28 journal: J Chem Inf Model DOI: 10.1021/acs.jcim.1c00263 sha: 16c12490f2c1f5c0372b4bb80fb8dcdcc8336278 doc_id: 709365 cord_uid: tym691he [Image: see text] The coronavirus SARS-CoV-2 main protease, M(pro), is conserved among coronaviruses with no human homolog and has therefore attracted significant attention as an enzyme drug target for COVID-19. The number of studies targeting M(pro) for in silico screening has grown rapidly, and it would be of great interest to know in advance how well docking methods can reproduce the correct ligand binding modes and rank these correctly. Clearly, current attempts at designing drugs targeting M(pro) with the aid of computational docking would benefit from a priori knowledge of the ability of docking programs to predict correct binding modes and score these correctly. In the current work, we tested the ability of several leading docking programs, namely, Glide, DOCK, AutoDock, AutoDock Vina, FRED, and EnzyDock, to correctly identify and score the binding mode of M(pro) ligands in 193 crystal structures. None of the codes were able to correctly identify the crystal structure binding mode (lowest energy pose with root-mean-square deviation < 2 Å) in more than 26% of the cases for noncovalently bound ligands (Glide: top performer), whereas for covalently bound ligands the top score was 45% (EnzyDock). These results suggest that one should perform in silico campaigns of M(pro) with care and that more comprehensive strategies including ligand free energy perturbation might be necessary in conjunction with virtual screening and docking. Coronaviruses are positive-stranded RNA viruses that infect humans and animals and cause common and severe respiratory diseases, including severe acute respiratory syndrome (SARS) and Middle East respiratory syndrome (MERS). 1,2 These viruses rely heavily on functional polypeptides that are generated by proteolytic cleavage of polyproteins that are translated from viral RNA. The principal coronavirus proteases responsible for this polypeptide formation are mainly protease and papain protease. The coronavirus SARS-CoV-2 main protease, M pro (henceforth denoted simply as M pro ), has garnered significant attention in the past year as an enzyme drug target due to the COVID-19 pandemic outbreak. 3 M pro is a druggable target 4,5 as it is conserved among coronaviruses and has no human homolog. The first M pro structures were published in early 2020. 3, 6, 7 The first crystal structures revealed that the active form of M pro is a homodimer containing two protomers, each composed of three domains ( Figure 1A ). The active site in M pro is located in a cleft between domains I and II ( Figure 1B) and features a noncanonical catalytic Cys−His dyad. The active site is composed of four regions: S1′, S1, S2, and S4, with the catalytic Cys145 located in S1′ and His41 located in S1′ and S2 ( Figure 1C ). To date, well over 200 three-dimensional holostructures have been resolved at a resolution of 3 Å or better. 3,6−9 In these structures, ligands bind to a variety of binding sites, including covalent and noncovalent binding in the main catalytic site, noncovalent binding in pockets in between the two proteomers, and weakly bound at the protein surface (i.e., in between crystallographic homodimer units). This wealth of structural information makes M pro an attractive benchmarking system for testing the ability of docking programs to correctly identify and rank the correct poses. This becomes particularly important in light of the severity of the COVID-19 pandemic 2 and the significant number of mutant forms of the virus that are rapidly appearing and might render current and future vaccines ineffective. To date, M pro has been studied extensively using ligand docking and screening tools 10−19 and computational enzymology tools, such as hybrid quantum mechanics−molecular mechanics (QM/MM). 20−29 The number of studies targeting M pro for in silico screening has grown, 30 and it would be of great interest to know in advance how well docking methods can reproduce the correct ligand binding modes and rank these correctly. Clearly, current attempts at designing drugs targeting M pro with the aid of computational docking are problematic if programs struggle to predict correct binding modes and score these correctly. This is true even if common docking programs have undergone extensive testing since each protein target comes with its own challenges due to the complexity of binding pocket crevices, nature of interactions, and solvent exposure. Therefore, critical evaluation of how well common docking programs perform for M pro is important. Since the development of the first automated docking program DOCK, 31 57 and OPLS 58,59 (i.e., molecular mechanics, MM), and knowledge-based scoring functions, such as DRUG-SCORE, 60 IT-SCORE, 61 DSX, 62 CHEMSCORE, 63 and SMoG. 64 Specialized docking programs addressing enzymes, 65, 66 such as EnzyDock, have also been developed. 67 Current challenges for docking methods include protein flexibility, 68 ligand solvation, and binding-site hydration. 47, 69 Thus far, several docking approaches have been employed to screen M pro for potential drugs in virtual screening and drug repurposing campaigns, including Glide, 7,10,13,17,18,70−72 Autodock, 11, 13, 73, 74 Autodock Vina, 11, 13, 14, 19, 71, 73, 75, 76 Surflex, 77 PLANT, 78 DockThor, 76 fast pulling of ligands, 14 deep docking, 70 algebraic topology and deep learning, 79 and virtual reality-based docking. 16 However, to the best of our knowledge, no rigorous benchmark study addressing the ability of such docking tools to reproduce and correctly rank known ligand binding modes has been published, in spite of the known inherent challenges in docking. 80−83 In the current work, we tested the ability of several leading docking programs to correctly identify and score the binding mode of M pro ligands in 193 crystal structures ( Figure S1 ). We tested the following docking programs: Glide, 40 39 FRED, 84−87 and EnzyDock. 67 The current results suggest that care should be taken in applying docking programs to a challenging protein target such as M pro . Preparing Mpro Structure Database. The available crystal structures of ligand-containing SARS-CoV-2 M pro were downloaded from the RCSB PDB website (March-December 2020). 88 In total, we collected 193 different structures, including covalent and noncovalent ligands (Table S1 and Figure S1 ). All structures were aligned relative to one reference structure (PDB ID: 5R84) for easy comparison. To perform docking, we separated each protein and ligand into separate files, removing crystal waters, ions, and cosolvents. Missing residues were added using the Modeller homology program. 89 Hydrogens were added using the CHARMM simulation platform (using HBUILD) for the protein structure and using Openbabel for the ligands. 90, 91 Visual inspection was also performed. For systems including only one monomer, the complementary unit was generated using the crystallographic information included in the PDB file using CHARMM. Protonation states of His residues were determined based on hydrogen bonding patterns and knowledge of the chemistry catalyzed by M pro (Table S2) , and they match the protonation states of key His residues recently published. 23 All docking simulations described below commenced with the CHARMM prepared systems. Clustering of Ligands and Water Molecules. Chemical descriptors were calculated for all ligands from the 193 PDB files using RDKit libraries in Python. Features thought to be important for ligand binding were chosen. Specifically, we computed the number of rotatable torsions, molecular weight, number of H-bond donors and acceptors, number of aromatic rings, the fraction of carbon atoms in sp 3 hybridization (relative to all carbon atoms in the molecule), and log P values. We applied principal component analysis (PCA) using these ligands descriptors, followed by k-means clustering to cluster the ligands into groups. We selected the number of clusters by silhouette analysis of the k-means clustering results. Ligand clustering was performed using Python 3.7. Density-Based Spatial Clustering of Applications with Noise (DBSCAN) clustering was performed to analyze the water molecules from all crystal structures. These water molecules were not included in the docking studies. Subsite-Binding Pocket Binding. To classify the binding patterns of all the 193 protein−ligand complexes, we categorized the ligands as bound on the surface, at the dimer interface, or in the active site. The latter was characterized according to the subsites S1, S1′, S2, and S4 ( Figure 1C , Table S3 ). 7,77 A ligand is considered to occupy a binding pocket if any ligand atom is within 4.0 Å of any pocket atom and also within 3.0 Å of the geometric center of the pocket (defined as geometric center of all pocket atoms). Moreover, a ligand is considered to be in the proximity of a binding pocket if any ligand atom is within 4.0 Å of any pocket atom and also within 5.0 Å of the geometric center of the pocket. The criteria were designed to account for both small and bulky ligands and to distinguish between binding poses where ligand groups are well docked inside a subpocket and poses where ligand groups are located at the periphery of a subpocket. Cutoff values are suitable for various nonbonded interactions (e.g., π−π stacking, hydrogen bond, ionic, and hydrophobic interactions). The final values were obtained by trial and error and validated by means of visual inspection. The subsite-binding pocket occupancy analysis was implemented as a CHARMM 90, 91 script. Docking Protocols. To compare the performance of selected docking programs for use with M pro (search algorithm and scoring function), we performed noncovalent docking using DOCK, 31−35 Autodock Vina, 39 and FRED 84−87 and noncovalent and covalent docking using Glide, 40−42 Autodock, 36−38 and EnzyDock 67 to the systems described above. In all docking simulations described below, the ligand was fully flexible, while the protein was fixed (except for the covalently connected complexes, where the appropriate Cys residue was flexible). Ligand Docking with Glide. Proteins and ligands were prepared using Schrodinger's Maestro (version 11.4, 2017-4 release) Prep Wiz and LigPrep modules, respectively, with default settings for docking with Glide. All covalent docking simulations were performed using the CovDock module available in Glide. 92 For the noncovalent simulations, the grid was generated using XGlide, which enables creation of different grids in parallel. The grids were centered around the ligand's centroid. The dimensions of the enclosing box and the bounding box were set to 12 × 12 × 12 Å 3 and 26 × 26 × 26 Å 3 , respectively, for all cases. The ligand stereochemistry was kept during all docking simulations. The number of poses written per ligand was set to 10,000. The scaling factors of the vdW radii and the partial atomic charge cutoff were set to the default values 0.80 and 0.15, respectively. Standard precision (SP) mode was chosen for all ligand docking runs. The selection of the best-docked ligand structure among the proposed poses is made based on several model energies implemented with Glide (docking score, Prime energy and E-model energy, and cdock affinity). Solvent effects were incorporated using MMGBSA. All reported energies herein used the docking score function for noncovalent docking and the cdock affinity scoring function for covalent docking as these performed best [i.e., produced the highest number of top-ranking structures with root-mean-square deviation (rmsd) < 2 Å]. In all Glide docking simulations (ligand preparation, protein preparation, grid generation, covalent, and noncovalent docking), the OPLS3 force field 59,93 was used. Ligand Docking with DOCK (Version 6.9). Proteins and ligands were prepared for docking simulations using the DockPrep option of Chimera v.14. 94 The grid was generated according to the center of mass of the crystal structure ligand with a grid spacing of 0.4 Å. The maximum and minimum radius of the sphere generated was set to 4 and 1.4 Å. All the spheres within 10 Å of the ligand were selected for docking. The box length surrounding the ligand was set to 6 Å, that is, the edge of the box from any atom of the ligand was at least 6 Å away, which easily accommodates all the selected spheres. Ligand Docking with Autodock (Version 4.2). Proteins and ligand pdbqt files were prepared using standard AutoDock tools (prepare_flexreceptor4.py and prepare_ligand4.py). These files include Cartesian coordinates and Gasteiger atomic charges 95 for each atom. AutoDock employs a united atom method, and thus, no nonpolar hydrogens are present. The center of mass of the crystallographic ligand was used to determine the center of the grid. AutoDock uses one grid box to perform the docking calculations, and the dimensions of this box were set to 37.5 × 37.5 × 37.5 Å 3 and the spacing was set to 0.375 Å for all systems. We performed flexible ligand docking into a rigid protein environment using GA, with default settings. For covalent docking, 96 each ligand was prepared with the active Cys residue already present in the input file using AutoDock tools (prepare_receptor4.py and prepare_flexreceptor4.py). For covalent docking, the ligand flexible torsional angles were presampled using MC simulations with CHARMM prior to docking. Ligand Docking with Autodock Vina. Protein and ligand input pdbqt files were prepared in the same way as for Autodock4.2 (see above). The size of the grid was set to 30.0 Ligand Docking with FRED. 85−87 FRED is one of the docking programs available within the OpenEye scientific library. For the docking process, proteins and ligands were prepared using the graphical user interface "Make Receptor" provided with OpenEye. FRED creates a potential field around the binding site by producing a negative image, which complements the shape of the protein site. This potential field is represented on a contour, which completely surrounds the ligand. OMEGA, an internal program within OpenEye is used to generate an ensemble of conformers for each ligand. A total of 200 different conformers were generated for each ligand for docking, and the 50 lowest energy docked structures were used to select the best pose in terms of lowest rmsd or Chemgauss energy scoring relative to the crystal ligand structure. The proteins and ligands were held fixed during the docking process. Ligand Docking with EnzyDock. Protein and ligand files were prepared as described at the beginning of Methods. CHARMM topology (RTF) and parameters (PRM) files for the ligands were generated using the CHARMM General Force Field (CGenFF) program. 53 To better understand the nature of the 193 M pro complexes prior to discussing the docking results, we present analyses of the ligands and their binding modes. The ligands were clustered into seven main groups based on their chemical features by PCA and k-means clustering. The features of each cluster were normalized, and the average value for each cluster was calculated ( Figure 2 ). The relative amount of the 193 ligands composing each cluster can be seen in the inserted pie chart ( Figure 2 ). For instance, cluster 6 is characterized by 17 ligands of low molecular weight and high fraction of carbon atoms with sp 3 hybridization, while cluster 3 is composed of 78 low-molecularweight ligands that are rather rigid and slightly hydrophobic. Cluster 4 has 15 ligands of high molecular weight, flexible chains with sp 3 carbons, and several hydrogen bond donors and acceptors, while cluster 7 has medium-molecular-weight ligands that are highly hydrophobic with aromatic rings, yet has some hydrogen bond donors and acceptors. Next, we analyzed the binding modes of the ligand clusters in M pro (Tables 1 and S1). In Table 1 we present the fraction of each cluster that is bound in a specific subpocket of the active site, at the surface, or at the interface between the two monomers. Note that ligands can bind in more than one pocket, and hence, the fractions do not add up to unity for each cluster. Inspection of the binding data clearly shows that ligands of low molecular weight (clusters 1, 3, 5, and 6) tend to bind at the surface of the protein (i.e., clusters 3, 5, and 6 are caught in between crystal units) or at the interface between the homodimer units (cluster 1). Still many low-molecular-weight ligands occupy pockets S1 and S1′ as these are covalently attached to C145. Ligands rich in aromatic rings and correspondingly high log P values (cluster 7) tend to occupy all pockets more than average, specifically sites S1 and S2. This is due to favorable π−π interactions with F140, H163, and H172 located in S1 or H41 and Y54 in S2. Another important observation is that ligands more likely to bind to S4 (which is rarely occupied) belong to clusters 2 and 4, which tend to include large, flexible molecules that are rich in H-bond donors and acceptors. We also clustered the water molecules in all crystal structures using DBSCAN clustering. Following clustering, we removed all waters that overlap any bound ligand ( Figure S2 ). These waters form a set of active site features that can be included in docking studies (these waters were not included in the current docking study). Docking of Ligands in M pro . We docked all ligands from 193 crystal M pro structures into their respective crystal protein structure (Table S4 ). These crystal structures include ligands bound noncovalently to the main binding pocket, surface and dimer interface, as well as covalently attached ligands. In all results below, we present the success rate of different docking programs in reproducing the crystal bound poses. For DOCK, AutoDock Vina, and FRED, the results reflect the noncovalent complexes only. In Figure 3A we show the overall success of all programs. Glide and EnzyDock reproduce the correct crystal structure pose (rmsd < 2 Å) for over 50% of the structures, with success rates of 64 and 70%, respectively, while for AutoDock, this rate falls to 40%. However, in many cases, even if the correct pose is identified, it is not scored as lowest in energy, and the success rate reduces to 33% (Glide), AutoDock (30%), and EnzyDock (35%). The overall success rates of Glide, AutoDock, and EnzyDock are in part due to the covalent complexes, whose poses are easier to reproduce than the noncovalent ones. If we analyze the success rates of identifying only the covalently bound complexes, Glide, AutoDock, and EnzyDock identify the correct poses 70, 42, and 71% of the cases, while the correct pose is also scored as the best one in 38, 36, and 45% of the cases ( Figure 3B ). For the noncovalent complexes, Glide, DOCK, AutoDock, AutoDock Vina, FRED, and EnzyDock identify the correct poses in 55, 61, 37, 29, 46 , and 68% of the cases, respectively, while these are ranked as the lowest energy poses in 26, 20, 24, 14, 14 , and 22% of the cases, respectively ( Figure 3C ). Finally, if we remove the complexes with surface-bound ligands (i.e., ligands bound in between crystal units), all programs perform significantly better ( Figure 3D ). The correct poses are identified (and scored as best) as follows (%): Glide 74 (39), DOCK 77 (29), (18), and EnzyDock 80 (35) , respectively. Next, we analyze how the different programs perform as a function of binding site locations on M pro . In Figures S3 and 4, we present box plots of the best rmsd values and the rmsd values for the lowest scoring pose for noncovalently bound ligands, respectively. All methods struggle with ligands bound at the interface between monomers and on the protein surface, Journal of Chemical Information and Modeling pubs.acs.org/jcim Article with Glide and FRED producing the best results. Additionally, most methods perform better for ligands bound to more than a single pocket (i.e., S1 + S2), and this trend is particularly clear for Glide and EnzyDock. Similarly, we analyze how the different covalent docking programs perform as a function of binding site locations on M pro . In Figures S4 and 5 , we present box plots of the best rmsd values and the rmsd values for the lowest scoring pose for covalently bound ligands, respectively. Also here, we observe a general trend, where the docking methods perform better for ligands occupying more pockets. In wake of the growing threat emerging from the SARS-CoV-2 pandemic, the modeling community has rushed to study a variety of potential pharmaceutical targets. One of these targets, M pro , is particularly attractive as this enzyme has no human analogue and is conserved among coronaviruses. A large number of studies have addressed ligand docking and virtual screening of ligand libraries against M pro in search of promising leads. A prerequisite for such studies is the ability of the docking programs to correctly identify and score ligand poses. Due to the intense efforts by the scientific community, there is already a wealth of structural biology information on M pro , hence enabling comparative studies of different docking approaches against this target. To date, the available crystal structures of M pro include ligands bound covalently and noncovalently to the main catalytic site, surface, and in between the two monomers. Here, we studied several leading docking codes, namely, Glide, DOCK, AutoDock, AutoDock Vina, FRED, and EnzyDock, and evaluate their ability to correctly reproduce and score the crystal structure ligand configuration for 193 M pro crystal structures. None of the codes are able to correctly identify and score the crystal structure in more than 26% of the cases for noncovalently bound ligands (Glide top performer), whereas for covalently bound ligands, the top score was 45% (EnzyDock best performer). Additionally, a general trend, where several of the docking methods (e.g., Glide and EnzyDock) perform better for larger, bulkier ligands occupying more than a single pocket, is observed. All docking methods struggle with prediction of small ligands. In the original crystal structures, many of the smaller ligands are surrounded by numerous explicit water molecules, dimethyl sulfoxide molecules, and Cl − ions that were removed prior to docking. Thus, these redocking trends may be ascribed to difficulty in accurately scoring docking poses, where a delicate balance between intra-and intermolecular terms and solvation terms must be stricken. In conclusion, the current results suggest that one should perform docking studies and virtual screening campaigns of M pro with care and that more comprehensive strategies might be necessary. Such strategies might include initial virtual screening (e.g., using FRED or AutoDock Vina) or docking (e.g., Glide or EnzyDock), followed by more rigorous ligand free energy binding calculations 98, 99 and in-depth QM/MM studies. 20, 24, 26, 28, 29 Inclusion of conserved water molecules, as identified in this study, may also be of help in guiding the docking process. Indeed, MD simulations have pointed to several water molecules, as important for M pro . 11 S.Z., K.R., and P.K.G. contributed equally. The docking simulations and system preparations were performed by all authors. The manuscript was written through contributions of all authors. All authors have given approval to the final version of the manuscript. The authors declare no competing financial interest. All preprepared 193 Mpro systems and accompanying Python and CHARMM scripts are available at https://github.com/ shanizev/Benchmarking-SARS-CoV-2. EnzyDock is freely available for noncommercial use on https://github.com/ majordt/EnzyDock. This work was supported by the Israel Ministry of Science, Technology and Space (grant 3-16310) and Israeli Science Foundation (grant # 1683/18). Coronavirus main proteinase (3CLpro) structure: Basis for design of anti-sars drugs Research and development on therapeutic agents and vaccines for COVID-19 and related human coronavirus diseases Crystal structure of SARS-CoV-2 main protease provides a basis for design of improved α-ketoamide inhibitors Potent SARS-CoV-2 direct-acting antivirals provide an important complement to COVID-19 vaccines Structure-based design of antiviral drug candidates targeting the SARS-CoV-2 main protease Structure of m(pro) from SARS-CoV-2 and discovery of its inhibitors Structural plasticity of SARS-CoV-2 3CL M pro active site cavity revealed by room temperature x-ray crystallography Fast identification of possible drug treatment of coronavirus disease-19 (COVID-19) through computational drug repurposing study Alpha-ketoamides as broad-spectrum inhibitors of coronavirus and enterovirus replication: Structure-based design, synthesis, and activity assessment Identification of 14 known drugs as inhibitors of the main protease of SARS-CoV-2 Rapid prediction of possible inhibitors for SARS-CoV-2 main protease using docking and fpl simulations Structural and evolutionary analysis indicate that the SARS-CoV-2 M pro is a challenging target for small-molecule inhibitor design Interactive molecular dynamics in virtual reality is an effective tool for flexible substrate and inhibitor docking to the SARS-CoV-2 main protease Structure-based virtual screening to discover potential lead molecules for the SARS-CoV-2 main protease Discovery of new hydroxyethylamine analogs against 3CL(pro) protein target of SARS-CoV-2: Molecular docking, molecular dynamics simulation, and structure-activity relationship studies Computational determination of potential inhibitors of SARS-CoV-2 main protease Revealing the molecular mechanisms of proteolysis of SARS-CoV-2 M pro by QM/MM computational methods Comprehensive insights into the catalytic mechanism of middle east respiratory syndrome 3c-like protease and severe acute respiratory syndrome 3c-like protease Dynamical properties of enzyme−substrate complexes disclose substrate specificity of the SARS-CoV-2 main protease as characterized by the electron density descriptors Proton-coupled conformational activation of sars coronavirus main proteases and opportunity for designing small-molecule broad-spectrum targeted covalent inhibitors Exploring the mechanism of covalent inhibition: Simulating the binding free energy of α-ketoamide inhibitors of the main protease of sars-cov-2 Tunón, I. A microscopic description of SARS-CoV-2 main protease inhibition with michael acceptors. Strategies for improving inhibitors design High-resolution mining of SARS-CoV-2 main protease conformational space: Supercomputer-driven unsupervised adaptive sampling Covalent and non-covalent binding free energy calculations for peptidomimetic inhibitors of SARS-CoV-2 main protease Mechanism of inhibition of SARS-CoV-2 M pro by n3 peptidyl michael acceptor explained by QM/ MM simulations and design of new derivatives with tunable chemical reactivity Crowdsourcing drug discovery for pandemics A geometric approach to macromolecule-ligand interactions Protein docking and complementarity Docking small-molecule ligands into active-sites Molecular recognition and docking algorithms Covalent docking of large libraries for the discovery of chemical probes Automated docking of substrates to proteins by simulated annealing Distributed automated docking of flexible ligands to proteins: Parallel applications of autodock 2.4 Autodock4 and autodocktools4: Automated docking with selective receptor flexibility Autodock vina: Improving the speed and accuracy of docking with a new scoring function, efficient optimization, and multithreading Glide: A new approach for rapid, accurate docking and scoring. 2. Enrichment factors in database screening Extra precision glide: Docking and scoring incorporating a model of hydrophobic enclosure for protein-ligand complexes Docking performance of the glide program as evaluated on the astex and dud datasets: A complete set of glide sp results and selected results for a new scoring function integrating watermap and glide Rosettaligand docking with full ligand and receptor flexibility Blind docking of pharmaceutically relevant compounds using rosettaligand Development and validation of a genetic algorithm for flexible docking Improved protein-ligand docking using GOLD Modeling water molecules in protein-ligand docking using GOLD III Assessing search strategies for flexible docking III Assessing energy functions for flexible docking Detailed analysis of grid-based molecular docking: A case study of CDOCKERa CHARMm-based md docking algorithm III Flexible CDOCKER: Development and application of a pseudo-explicit structure-based docking method within CHARMM All-atom empirical potential for molecular modeling and dynamics studies of proteins CHARMM general force field: A force field for drug-like molecules compatible with the CHARMM all-atom additive biological force fields CHARMM36 all-atom additive protein force field: Validation based on comparison to nmr data Charmm36m: An improved force field for folded and intrinsically disordered proteins Improved modeling of halogenated ligand−protein interactions using the drude polarizable and CHARMM additive empirical force fields A second generation force field for the simulation of proteins, nucleic acids, and organic molecules The OPLS force field for proteins. Energy minimizations for crystals of cyclic peptides and crambin Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids Drugscore(CSD)-knowledge-based scoring function derived from small molecule crystal data with superior recognition rate of near-native ligand poses and better affinity prediction Automated large-scale file preparation, docking, and scoring: Evaluation of itscore and stscore using the 2012 community structureactivity resource benchmark Towards the development of universal, fast and highly accurate docking/scoring methods: A long way to go An improved knowledge-based scoring function for protein-ligand interactions Mechanistically informed predictions of binding modes for carbocation intermediates of a sesquiterpene synthase reaction Predicting productive binding modes for substrates and carbocation intermediates in terpene synthases-bornyl diphosphate synthase as a representative case Protein-ligand docking of multiple reactive states along a reaction coordinate in enzymes Receptor flexibility in small-molecule docking calculations Automated docking to multiple target structures: Incorporation of protein mobility and structural water heterogeneity in autodock Rapid identification of potential inhibitors of SARS-CoV-2 main protease by deep docking of 1.3 billion compounds Targeting the SARS-CoV-2 main protease using FDAapproved isavuconazonium, a p2-p3 alpha-ketoamide derivative and pentagastrin: An in-silico drug discovery approach Putative inhibitors of SARS-CoV-2 main protease from a library of marine natural products: A virtual screening and molecular modeling study Molecular docking, validation, dynamics simulations, and pharmacokinetic prediction of natural compounds against the SARS-CoV-2 main-protease Rational approach toward COVID-19 main protease inhibitors via molecular docking, molecular dynamics simulation and free energy calculation Ligand and structure-based virtual screening applied to the SARS-CoV-2 main protease: An in silico repurposing study Optimization rules for SARS-CoV-2 Mpro antivirals: Ensemble docking and exploration of the coronavirus protease active site In silico drug repurposing for SARS-CoV-2 main proteinase and spike proteins Unveiling the molecular mechanism of SARS-CoV-2 main protease inhibition from 137 crystal structures using algebraic topology and deep learning Diverse, high-quality test set for the validation of protein-ligand native and cross-docking performance Can we trust docking results? Evaluation of seven commonly used programs on PDBbind database Beware of docking! Predictive power of different types of experimental restraints in small molecule docking: A review Comparison of shape-matching and docking as virtual screening tools FRED pose prediction and virtual screening accuracy FRED and hybrid docking performance on standardized datasets POSIT: Flexible shape-guided docking for pose prediction Comparative protein modelling by satisfaction of spatial restraints CHARMM -a program for macromolecular energy, minimization, and dynamics calculations The biomolecular simulation program Docking covalent inhibitors: A parameter free approach to pose prediction and scoring OPLS3: A force field providing broad coverage of drug-like small molecules and proteins UCSF Chimeraa visualization system for exploratory research and analysis Iterative partial equalization of orbital electronegativitya rapid access to atomic charges Covalent docking using autodock: Two-point attractor and flexible side chain methods Extension of the CHARMM general force field to sulfonylcontaining compounds and its utility in biomolecular simulations On achieving high accuracy and reliability in the calculation of relative protein−ligand binding affinities Potent noncovalent inhibitors of the main protease of SARS-CoV-2 from molecular sculpting of the drug perampanel guided by free energy perturbation calculations