Deep learning has driven major breakthroughs in protein structure prediction; however, one of the next critical steps forward is accurately predicting how proteins interact with small-molecule ligands, to enable real-world applications such as drug discovery. Recent cofolding methods aim to address this challenge, but evaluating their performance has been inconclusive because of the lack of relevant benchmarking datasets. Here we present a comprehensive evaluation of four leading all-atom cofolding methods using our newly introduced benchmark dataset, Runs N' Poses. Runs N' Poses comprises 2,600 high-resolution protein-ligand systems released after the training cutoff used by these methods. We demonstrate that current cofolding approaches largely memorize ligand poses from their training data, hindering their use for de novo drug design. With this assessment and benchmark dataset, we aim to accelerate progress in the field by allowing for a more realistic assessment of the current state-of-the-art deep learning methods for predicting protein-ligand interactions.

Deep learning has driven major breakthroughs in protein structure prediction, however the next critical advance is accurately predicting how proteins interact with other molecules, especially small molecule ligands, to enable real-world applications such as drug discovery and design. Recent deep learning all-atom methods have been built to address this challenge, but evaluating their performance on the prediction of protein-ligand complexes has been inconclusive due to the lack of relevant benchmarking datasets. Here we present a comprehensive evaluation of four leading all-atom cofolding deep learning methods using our newly introduced benchmark dataset Runs N' Poses, which comprises 2,600 high-resolution protein-ligand systems released after the training cutoff used by these methods. We demonstrate that current co-folding approaches largely memorise ligand poses from their training data, hindering their use for de novo drug design. This limitation is especially pronounced for ligands that have only been seen binding in one pocket, whereas more promiscuous ligands such as cofactors show moderately improved performance. With this work and benchmark dataset, we aim to accelerate progress in the field by allowing for a more realistic assessment of the current state-of-the-art deep learning methods for predicting protein-ligand interactions. ### Competing Interest Statement The authors have declared no competing interest.

Molecule parametrization is an essential requirement to guarantee the accuracy of docking calculations. Parametrization includes a proper perception of chemical properties such as bonds, formal charges and protonation states. This includes large biological macromolecules, such as proteins and nucleic acids, and small molecules, such as ligands and cofactors. The structures of proteins and nucleic acids are challenging due to omission of several atoms from the structural model, and from the lack of connectivity and bond order information in the PDB and mmCIF file formats. For small molecules, the very large chemical diversity poses challenges for both validating correctness and providing accurate parameters. These challenges affect various modeling approaches like molecular docking and molecular dynamics. Moreover, several specialized methods (particularly in molecular docking) leverage specific chemical properties to add custom potentials, pseudoatoms, or manipulate atomic connectivity. To address these challenges, we developed Meeko, a molecular parametrization Python package that leverages the widely used RDKit cheminformatics library for a chemically accurate description of the molecular representation. Small molecules are modeled as single RDKit molecules, and biological macromolecules as multiple RDKit molecules, one for each residue. Meeko is highly customizable and designed to be easily scriptable for high-throughput processing, replacing MGLTools for receptor and ligand preparation.

Cysteine residues play key roles in protein structure and function and can serve as targets for chemical probes and even drugs. Chemoproteomic studies have revealed that heightened cysteine reactivity toward electrophilic probes, such as iodoacetamide alkyne (IAA), is indicative of likely residue functionality. However, while the cysteine coverage of chemoproteomic studies has increased substantially, these methods still provide only a partial assessment of proteome-wide cysteine reactivity, with cysteines from low-abundance proteins and tough-to-detect peptides still largely refractory to chemoproteomic analysis. Here, we integrate cysteine chemoproteomic reactivity data sets with structure-guided computational analysis to delineate key structural features of proteins that favor elevated cysteine reactivity toward IAA. We first generated and aggregated multiple descriptors of cysteine microenvironment, including amino acid content, solvent accessibility, residue proximity, secondary structure, and predicted pKa. We find that no single feature is sufficient to accurately predict the reactivity. Therefore, we developed the CIAA (Cysteine reactivity toward IodoAcetamide Alkyne) method, which utilizes a Random Forest model to assess cysteine reactivity by incorporating descriptors that characterize the three-dimensional (3D) structural properties of thiol microenvironments. We trained the CIAA model on existing and newly generated cysteine chemoproteomic reactivity data paired with high-resolution crystal structures from the Protein Data Bank (PDB), with cross-validation against an external data set. CIAA analysis reveals key features driving cysteine reactivity, such as backbone hydrogen bond donor atoms, and reveals still underserved needs in the area of computational predictions of cysteine reactivity, including challenges surrounding protein structure selection data set curation. Thus, our work provides a strong foundation for deploying artificial intelligence (AI) on cysteine chemoproteomic data sets.

The protein-ligand component of the 16th Critical Assessment of Structure Prediction (CASP16) challenged participants to predict both binding poses and affinities of small molecules to protein targets, with a focus on drug-like compounds from pharmaceutical discovery projects. Thirty research groups submitted predictions for 229 protein-ligand pose targets and 140 affinity targets across five protein systems. Among the submitted predictions, template-based pose-prediction methods did particularly well, with the best groups achieving mean LDDT-PLI values of 0.69 (scale of 0-1 with 1 best). For comparison, we also ran a set of automated baseline pose-prediction methods, including ones using deep neural networks. Of these, AlphaFold 3 did particularly well, with a mean LDDT-PLI of 0.8, thus outscoring the best CASP16 predictor. The CASP affinity predictions showed modest correlation with experimental data (maximum Kendall's τ = 0.42), well below the theoretical maximum possible given experimental uncertainty (~0.73). As seen in prior challenges, providing experimental structures did not improve affinity predictions in the second stage of the challenge, suggesting that the scoring functions used here are a key limiting factor. Overall, the accuracy achieved by CASP participants is similar to that observed in the prior Drug Design Data Resource (D3R) blinded prediction challenges. The present results highlight the progress and persistent challenges in computational protein-ligand modeling and provide valuable benchmarks for the field of computer-aided drug design.

Cosolvent molecular dynamics (MD) are an increasingly popular form of simulations where small molecule cosolvents are added to water-solvated protein systems. These simulations can perform diverse target characterization tasks, including cryptic and allosteric pocket identification and pharmacophore profiling, and supplement suites of enhanced sampling methods to explore protein conformational landscapes. The behavior of these systems is tied to the cosolvents used, so the ability to define diverse and complex mixtures is critical in dictating the outcome of the simulations. However, existing methods for preparing cosolvent simulations only support a limited number of predefined cosolvents and concentrations. Here we present CosolvKit, a tool for the preparation and analysis of systems composed of user-defined cosolvents and concentrations. This tool is modular and agnostic of the MD engine and force field used, offering access to a variety of generalizable small molecule force fields. To the best of our knowledge, CosolvKit represents the first generalized approach for the construction of these simulations.

Water desolvation is one of the key components of the free energy of binding of small molecules to their receptors. Thus, understanding the energetic balance of solvation and desolvation resulting from individual water molecules can be crucial when estimating ligand binding, especially when evaluating different molecules and poses as done in High-Throughput Virtual Screening (HTVS). Over the most recent decades, several methods were developed to tackle this problem, ranging from fast approximate methods (usually empirical functions using either discrete atom-atom pairwise interactions or continuum solvent models) to more computationally expensive and accurate ones, mostly based on Molecular Dynamics (MD) simulations, such as Grid Inhomogeneous Solvation Theory (GIST) or Double Decoupling. On one hand, MD-based methods are prohibitive to use in HTVS to estimate the role of waters on the fly for each ligand. On the other hand, fast and approximate methods show an unsatisfactory level of accuracy, with low agreement with results obtained with the more expensive methods. Here we introduce WaterKit, a new grid-based sampling method with explicit water molecules to calculate thermodynamic properties using the GIST method. Our results show that the discrete placement of water molecules is successful in reproducing the position of crystallographic waters with very high accuracy, as well as providing thermodynamic estimates with accuracy comparable to more expensive MD simulations. Unlike these methods, WaterKit can be used to analyze specific regions on the protein surface, (such as the binding site of a receptor), without having to hydrate and simulate the whole receptor structure. The results show the feasibility of a general and fast method to compute thermodynamic properties of water molecules, making it well-suited to be integrated in high-throughput pipelines such as molecular docking.

Prediction categories in the Critical Assessment of Structure Prediction (CASP) experiments change with the need to address specific problems in structure modeling. In CASP15, four new prediction categories were introduced: RNA structure, ligand-protein complexes, accuracy of oligomeric structures and their interfaces, and ensembles of alternative conformations. This paper lists technical specifications for these categories and describes their integration in the CASP data management system.

The prediction of protein-ligand complexes (PLC), using both experimental and predicted structures, is an active and important area of research, underscored by the inclusion of the Protein-Ligand Interaction category in the latest round of the Critical Assessment of Protein Structure Prediction experiment CASP15. The prediction task in CASP15 consisted of predicting both the three-dimensional structure of the receptor protein as well as the position and conformation of the ligand. This paper addresses the challenges and proposed solutions for devising automated benchmarking techniques for PLC prediction. The reliability of experimentally solved PLC as ground truth reference structures is assessed using various validation criteria. Similarity of PLC to previously released complexes are employed to judge PLC diversity and the difficulty of a PLC as a prediction target. We show that the commonly used PDBBind time-split test-set is inappropriate for comprehensive PLC evaluation, with state-of-the-art tools showing conflicting results on a more representative and high quality dataset constructed for benchmarking purposes. We also show that redocking on crystal structures is a much simpler task than docking into predicted protein models, demonstrated by the two PLC-prediction-specific scoring metrics created. Finally, we introduce a fully automated pipeline that predicts PLC and evaluates the accuracy of the protein structure, ligand pose, and protein-ligand interactions.

CASP15 introduced a new category, ligand prediction, where participants were provided with a protein or nucleic acid sequence, SMILES line notation, and stoichiometry for ligands and tasked with generating computational models for the three-dimensional structure of the corresponding protein-ligand complex. These models were subsequently compared with experimental structures determined by x-ray crystallography or cryoEM. To assess these predictions, two novel scores were developed. The Binding-Site Superposed, Symmetry-Corrected Pose Root Mean Square Deviation (BiSyRMSD) evaluated the absolute deviations of the models from the experimental structures. At the same time, the Local Distance Difference Test for Protein-Ligand Interactions (lDDT-PLI) assessed the ability of models to reproduce the protein-ligand interactions in the experimental structures. The ligands evaluated in this challenge range from single-atom ions to large flexible organic molecules. More than 1800 submissions were evaluated for their ability to predict 23 different protein-ligand complexes. Overall, the best models could faithfully reproduce the geometries of more than half of the prediction targets. The ligands' size and flexibility were the primary factors influencing the predictions' quality. Small ions and organic molecules with limited flexibility were predicted with high fidelity, while reproducing the binding poses of larger, flexible ligands proved more challenging.

This study introduces a novel Bayesian Optimization (BO) method to support the design and optimization of bioactive peptide sequences in the context of a fully automated closed-loop Design-Make-Test (DMT) pipeline. Using the major histocompatibility complex class I receptor system as test case, we showed that BO is capable to efficiently navigate vast sequence spaces. Starting from a single peptide-lead sequence in the $\mu$M IC50 range, the method is able to optimize a peptide sequence to its optimal binding affinity in less than 5 DMT cycles, with 96 peptide sequences per batch. We extensively evaluated its performance, in various conditions and with different parameters, providing valuable insights for peptide optimization tasks in future closed-loop DMT environments. Different sequence- and structure-based initialization strategies were also tested, to generate the initial batch of peptide sequences, as well as different molecular fingerprints and protein language models. Additionally, the method developed here can natively handle various peptide sequence lengths and scaffolds (e.g. macrocycles) and support any arbitrary non-standard amino acids or residue modifications. The source code of our method, Mobius, is publicly available under the Apache license at https://git.scicore.unibas.ch/schwede/mobius.

AutoDock Vina is arguably one of the fastest and most widely used open-source docking engines. However, compared to other docking engines in the AutoDock Suite, it lacks features that support modeling of specific systems such as macrocycles or modeling water explicitly. Here, we describe the implementation of these functionality in AutoDock Vina 1.2.0. Additionally, AutoDock Vina 1.2.0 supports the AutoDock4.2 scoring function, simultaneous docking of multiple ligands, and a batch mode for docking a large number of ligands. Furthermore, we implemented Python bindings to facilitate scripting and the development of docking workflows. This work is an effort toward the unification of the features of the AutoDock4 and AutoDock Vina docking engines. The source code is available at https://github.com/ccsb-scripps/AutoDock-Vina

Activity-based protein profiling (ABPP) has been used extensively to discover and optimize selective inhibitors of enzymes. Here, we show that ABPP can also be implemented to identify the converse-small-molecule enzyme activators. Using a kinetically controlled, fluorescence polarization-ABPP assay, we identify compounds that stimulate the activity of LYPLAL1-a poorly characterized serine hydrolase with complex genetic links to human metabolic traits. We apply ABPP-guided medicinal chemistry to advance a lead into a selective LYPLAL1 activator suitable for use in vivo. Structural simulations coupled to mutational, biochemical and biophysical analyses indicate that this compound increases LYPLAL1's catalytic activity likely by enhancing the efficiency of the catalytic triad charge-relay system. Treatment with this LYPLAL1 activator confers beneficial effects in a mouse model of diet-induced obesity. These findings reveal a new mode of pharmacological regulation for this large enzyme family and suggest that ABPP may aid discovery of activators for additional enzyme classes.

AUTODOCK is a molecular docking software widely used in computational drug design. Its time-consuming executions have motivated the development of AUTODOCK-GPU, an OpenCL-accelerated version that can run on GPUs and CPUs. This work discusses the development of AUTODOCK-GPU from a programming perspective, detailing how our design addresses the irregularity of AUTODOCK while pushing towards higher performance. Details on required data transformations, re-structuring of complex functionality, as well as the performance impact of different configurations are also discussed. While AUTODOCK-GPU reaches speedup factors of 341x on a Titan V GPU and 51x on a 48-core Xeon Platinum 8175M CPU, experiments show that performance gains are highly dependent on the molecular complexity under analysis. Finally, we summarize our preliminary experiences when porting AUTODOCK onto FPGAs.

Caspases are a critical class of proteases involved in regulating programmed cell death and other biological processes. Selective inhibitors of individual caspases, however, are lacking, due in large part to the high structural similarity found in the active sites of these enzymes. We recently discovered a small-molecule inhibitor, 63-R, that covalently binds the zymogen, or inactive precursor (pro-form), of caspase-8, but not other caspases, pointing to an untapped potential of procaspases as targets for chemical probes. Realizing this goal would benefit from a structural understanding of how small molecules bind to and inhibit caspase zymogens. There have, however, been very few reported procaspase structures. Here, we employ X-ray crystallography to elucidate a procaspase-8 crystal structure in complex with 63-R, which reveals large conformational changes in active-site loops that accommodate the intramolecular cleavage events required for protease activation. Combining these structural insights with molecular modeling and mutagenesis-based biochemical assays, we elucidate key interactions required for 63-R inhibition of procaspase-8. Our findings inform the mechanism of caspase activation and its disruption by small molecules and, more generally, have implications for the development of small molecule inhibitors and/or activators that target alternative (e.g., inactive precursor) protein states to ultimately expand the druggable proteome.

Molecular docking has been successfully used in computer-aided molecular design projects for the identification of ligand poses within protein binding sites. However, relying on docking scores to rank different ligands with respect to their experimental affinities might not be sufficient. It is believed that the binding scores calculated using molecular mechanics combined with the Poisson–Boltzman surface area (MM-PBSA) or generalized Born surface area (MM-GBSA) can predict binding affinities more accurately. In this perspective, we decided to take part in Stage 2 of the Drug Design Data Resource (D3R) Grand Challenge 4 (GC4) to compare the performance of a quick scoring function, AutoDock4, to that of MM-GBSA in predicting the binding affinities of a set of $$\beta$$-Amyloid Cleaving Enzyme 1 (BACE-1) ligands. Our results show that re-scoring docking poses using MM-GBSA did not improve the correlation with experimental affinities. We further did a retrospective analysis of the results and found that our MM-GBSA protocol is sensitive to details in the protein-ligand system: (i) neutral ligands are more adapted to MM-GBSA calculations than charged ligands, (ii) predicted binding affinities depend on the initial conformation of the BACE-1 receptor, (iii) protonating the aspartyl dyad of BACE-1 correctly results in more accurate binding affinity predictions.

The thiazolidinedione (TZD) pioglitazone (Pio) is an FDA-approved drug for type 2 diabetes mellitus that binds and activates the nuclear receptor peroxisome proliferator-activated receptor gamma (PPARγ). Although TZDs have potent antidiabetic effects, they also display harmful side effects that have necessitated a better understanding of their mechanisms of action. In particular, little is known about the effect of in vivo TZD metabolites on the structure and function of PPARγ. Here, we present a structure-function comparison of Pio and a major in vivo metabolite, 1-hydroxypioglitazone (PioOH). PioOH displayed a lower binding affinity and reduced potency in coregulator recruitment assays compared to Pio. To determine the structural basis of these findings, we solved an X-ray crystal structure of PioOH bound to PPARγ ligand-binding domain (LBD) and compared it to a published Pio-bound crystal structure. PioOH exhibited an altered hydrogen bonding network that could underlie its reduced affinity and potency compared to Pio. Solution-state structural analysis using NMR spectroscopy and hydrogen/deuterium exchange mass spectrometry (HDX-MS) analysis revealed that PioOH stabilizes the PPARγ activation function-2 (AF-2) coactivator binding surface better than Pio. In support of AF-2 stabilization, PioOH displayed stabilized coactivator binding in biochemical assays and better transcriptional efficacy (maximal transactivation response) in a cell-based assay that reports on the activity of the PPARγ LBD. These results, which indicate that Pio hydroxylation affects both its potency and efficacy as a PPARγ agonist, contribute to our understanding of PPARγ-binding drug metabolite interactions and may assist in future PPARγ drug design efforts.

Allosteric regulation plays an important role in many biological processes, such as signal transduction, transcriptional regulation, andmetabolism. Allostery is rooted in the fundamental physical properties of macromolecular systems, but its underlying mechanisms are still poorly understood. A collection of contributions to a recent interdisciplinary CECAM(Center Europeen de Calcul Atomique et Moleculaire) workshop is used here to provide an overview of the progress and remaining limitations in the understanding of the mechanistic foundations of allostery gained from computational and experimental analyses of real protein systems and model systems. The main conceptual frameworks instrumental in driving the field are discussed. We illustrate the role of these frameworks in illuminating molecular mechanisms and explaining cellular processes, and describe some of their promising practical applications in engineering molecular sensors and informing drug design efforts.

Molecular docking has been successfully used in computer-aided molecular design projects for the identification of ligand poses within protein binding sites. However, relying on docking scores to rank different ligands with respect to their experimental affinities might not be sufficient. It is believed that the binding scores calculated using molecular mechanics combined with the Poisson-Boltzman surface area (MM-PBSA) or generalized Born surface area (MM-GBSA) can more accurately predict binding affinities. In this perspective, we decided to take part in Stage 2 in the Drug Design Data Resource (D3R) Grand Challenge 4 (GC4) to compare the performance of a quick scoring function, Autodock4, to that of MM-GBSA in predicting the binding affinities of a set of Beta-Amyloid Cleaving Enzyme 1 (BACE-1) ligands. Our results show that re-scoring docking poses using MM-GBSA did not improve the correlation with experimental affinities. We further did a retrospective analysis of the results and found that our MM-GBSA protocol is sensitive to details in the protein-ligand system: i) neutral ligands are more adapted to MM-GBSA calculations than charged ligands, ii) predicted binding affinities depend on the initial conformation of the BACE-1 receptor, iii) protonating the aspartyl dyad of BACE-1 correctly results in more accurate binding pose and affinity predictions.