Comparative Performance of High-Throughput Methods for Protein pKa Predictions.

The medically relevant field of protein-based therapeutics has triggered a demand for protein engineering in different pH environments of biological relevance. In silico engineering workflows typically employ high-throughput screening campaigns that require evaluating large sets of protein residues and point mutations by fast yet accurate computational algorithms. While several high-throughput pKa prediction methods exist, their accuracies are unclear due to the lack of a current comprehensive benchmarking. Here, seven fast, efficient, and accessible approaches including PROPKA3, DeepKa, PKAI, PKAI+, DelPhiPKa, MCCE2, and H++ were systematically tested on a nonredundant subset of 408 measured protein residue pKa shifts from the pKa database (PKAD). While no method outperformed the null hypotheses with confidence, as illustrated by statistical bootstrapping, DeepKa, PKAI+, PROPKA3, and H++ had utility. More specifically, DeepKa consistently performed well in tests across multiple and individual amino acid residue types, as reflected by lower errors, higher correlations, and improved classifications. Arithmetic averaging of the best empirical predictors into simple consensuses improved overall transferability and accuracy up to a root-mean-square error of 0.76 pKa units and a correlation coefficient (R2) of 0.45 to experimental pKa shifts. This analysis should provide a basis for further methodological developments and guide future applications, which require embedding of computationally inexpensive pKa prediction methods, such as the optimization of antibodies for pH-dependent antigen binding.

Solvated interaction energy: from small-molecule to antibody drug design.

Scoring functions are ubiquitous in structure-based drug design as an aid to predicting binding modes and estimating binding affinities. Ideally, a scoring function should be broadly applicable, obviating the need to recalibrate and refit its parameters for every new target and class of ligands. Traditionally, drugs have been small molecules, but in recent years biologics, particularly antibodies, have become an increasingly important if not dominant class of therapeutics. This makes the goal of having a transferable scoring function, i.e., one that spans the range of small-molecule to protein ligands, even more challenging. One such broadly applicable scoring function is the Solvated Interaction Energy (SIE), which has been developed and applied in our lab for the last 15 years, leading to several important applications. This physics-based method arose from efforts to understand the physics governing binding events, with particular care given to the role played by solvation. SIE has been used by us and many independent labs worldwide for virtual screening and discovery of novel small-molecule binders or optimization of known drugs. Moreover, without any retraining, it is found to be transferrable to predictions of antibody-antigen relative binding affinities and as accurate as functions trained on protein-protein binding affinities. SIE has been incorporated in conjunction with other scoring functions into ADAPT (Assisted Design of Antibody and Protein Therapeutics), our platform for affinity modulation of antibodies. Application of ADAPT resulted in the optimization of several antibodies with 10-to-100-fold improvements in binding affinity. Further applications included broadening the specificity of a single-domain antibody to be cross-reactive with virus variants of both SARS-CoV-1 and SARS-CoV-2, and the design of safer antibodies by engineering of a pH switch to make them more selective towards acidic tumors while sparing normal tissues at physiological pH.

Antibody mutations favoring pH-dependent binding in solid tumor microenvironments: Insights from large-scale structure-based calculations.

Antibody-based therapeutics for treatment of various tumors have grown rapidly in recent years. Unfortunately, safety issues, attributed to off-tumor effects and cytotoxicity, are still a significant concern with the standard of care. Improvements to ensure targeted delivery of antitumor pharmaceuticals are desperately needed. We previously demonstrated that incorporating histidyl pH-switches in an anti-HER2 antibody induced selective antigen binding under acidic pH conditions (MAbs 2020;12:1682866). This led to an improved safety profile due to preferential targeting of the oncoprotein in the acidic solid tumor microenvironment. Following this success, we expanded this approach to a set of over 400 antibody structures complexed with over 100 different human oncoproteins, associated with solid tumors. Calculations suggested that mutations to His of certain residue types, namely Trp, Arg, and Tyr, could be significantly more successful for inducing pH-dependent binding under acidic conditions. Furthermore, 10 positions within the complementarity-determining region were also predicted to exhibit greater successes. Combined, these two accessible metrics could serve as the basis for a sequence-based engineering of pH-selective binding. This approach could be applied to most anticancer antibodies, which lack detailed structural characterization.

Mutations in Dynamic Structural Elements Alter the Kinetics and Fidelity of the Multifunctional Class II Lanthipeptide Synthetase, HalM2.

Class II lanthipeptide synthetases (LanM enzymes) catalyze the multistep post-translational modification of genetically encoded precursor peptides into macrocyclic (often antimicrobial) lanthipeptides. The reaction sequence involves dehydration of serine/threonine residues, followed by intramolecular addition of cysteine thiols onto the nascent dehydration sites to construct thioether bridges. LanMs utilize two separate active sites in an iterative yet highly coordinated manner to maintain a remarkable level of regio- and stereochemical control over the multistep maturation. The mechanisms underlying this biosynthetic fidelity remain enigmatic. We recently demonstrated that proper function of the haloduracin ß synthetase (HalM2) requires dynamic structural elements scattered across the surface of the enzyme. Here, we perform kinetic simulations, structural analysis of reaction intermediates, hydrogen-deuterium exchange mass spectrometry studies, and molecular dynamics simulations to investigate the contributions of these dynamic HalM2 structural elements to biosynthetic efficiency and fidelity. Our studies demonstrate that a large, conserved loop (HalM2 residues P349-P405) plays essential roles in defining the precursor peptide binding site, facilitating efficient peptide dehydration, and guiding the order of thioether ring formation. Moreover, mutations near the interface of the HalM2 dehydratase and cyclase domains perturb cyclization fidelity and result in aberrant thioether topologies that cannot be corrected by the wild type enzyme, suggesting an element of kinetic control in the normal cyclization sequence. Overall, this work provides the most comprehensive correlation of the structural and functional properties of a LanM enzyme reported to date and should inform mechanistic studies of the biosynthesis of other ribosomally synthesized and post-translationally modified peptide natural products.

Use of Extended-Hückel Descriptors for Rapid and Accurate Predictions of Conjugated Torsional Energy Barriers.

Over the past few decades, virtual high-throughput screening (vHTS) and molecular dynamics simulations have become effective and widely used tools in the initial stages of drug discovery efforts. These methods allow a great number of druglike molecules to be screened quickly and inexpensively. Unfortunately, however, the accuracies of both these methods rely on the quality of the underlying molecular mechanics force fields (FFs), which are often poor. This major weakness originates from the reliance of FFs on a finite list of specific parameters, called atom types, which have low transferability between molecules. In particular, the torsional energy barriers of druglike molecules are notoriously difficult to predict. Continuing our endeavor to understand factors affecting the torsional energy barriers of small molecules and quantify them, we showed that descriptors calculated using the extended-Hückel method could be used to rapidly assign accurate torsion parameters for conjugated molecules. This method, called H-TEQ 4.5, was developed using a set of 684 conjugated molecules. It was subsequently validated on a test set of 200 diverse molecules and produced an average root-mean-square error (rmse) of 1.01 kcal·mol-1, with respect to the reference quantum mechanic torsional profiles. For comparison, GAFF2, MMFF94, and MAB produced average rmse's of 3.49, 1.50, and 1.77 kcal·mol-1, respectively. H-TEQ 4.5 is also computationally inexpensive, running just under 0.25 ms for a biphenyl molecule on a home computer, allowing it to be used for vHTS of large libraries of compounds. Overall, H-TEQ 4.5 solved the problems associated with the transferability of torsion parameters for conjugated molecules. This method was incorporated into the Molecular Operating Environment and will be available for a wide variety of applications.

Atom Type Independent Modeling of the Conformational Energy of Benzylic, Allylic, and Other Bonds Adjacent to Conjugated Systems.

Applications of computational methods to predict binding affinities for protein/drug complexes are routinely used in structure-based drug discovery. Applications of these methods often rely on empirical force fields (FFs) and their associated parameter sets and atom types. However, it is widely accepted that FFs cannot accurately cover the entire chemical space of drug-like molecules, due to the restrictive cost of parametrization and the poor transferability of existing parameters. To address these limitations, initiatives have been carried out to develop more transferable methods, in order to allow for more rigorous descriptions of any drug-like molecule. We have previously reported H-TEQ, a method which does not rely on atom types and incorporates well established chemical principles to assign parameters to organic molecules. The previous implementation of H-TEQ (a torsional barrier prediction method) only covered saturated and lone pair containing molecules; here, we report our efforts to incorporate conjugated systems into our model. The next step was the evaluation of the introduction of unsaturations. The developed model (H-TEQ3.0) has been validated on a wide variety of molecules containing heteroaromatic groups, alkyls, and fused ring systems. Our method performs on par with one of the most commonly used FFs (GAFF2), without relying on atom types or any prior parametrization.

Torsional Energy Barriers of Biaryls Could Be Predicted by Electron Richness/Deficiency of Aromatic Rings; Advancement of Molecular Mechanics toward Atom-Type Independence.

Biaryl molecules are ubiquitous pharmacophores found in natural products and pharmaceuticals. In spite of this, existing molecular mechanics force fields are unable to accurately reproduce their torsional energy profiles, except for a few well-parametrized cases. This effectively limits the ability of structure-based drug design methods to correctly identify hits involving biaryls with confidence (e.g., during virtual screening, employing docking and/or molecular dynamics simulations). Continuing in our endeavor to quantify organic chemistry principles, we showed that the torsional energy profile of biaryl compounds could be computed on-the-fly based on the electron richness/deficiency of the aromatic rings. This method, called H-TEQ 4.0, was developed using a set of 131 biaryls. It was subsequently validated on a separate set of 100 diverse biaryls, including multisubstituted, bicyclic and tricyclic druglike molecules, and produced an average root-mean-square error (RMSE) of 0.95 kcal·mol-1. For comparison, GAFF2 produced an RMSE of 3.88 kcal·mol-1, owing to problems associated with the transferability of torsion parameters. The success of H-TEQ 4.0 provided further evidence that force fields could transition to become atom-type independent, providing that the correct chemical principles are used. Overall, this method solved the problem of transferability of biaryl torsion parameters, while simultaneously improving the overall accuracy of the force field.

Predicting Positions of Bridging Water Molecules in Nucleic Acid-Ligand Complexes.

Over the past two decades, interests in DNA and RNA as drug targets have been growing rapidly. Following the trends observed with protein drug targets, computational approaches for drug design have been developed for this new class of molecules. Our efforts toward the development of a universal docking program, Fitted, led us to focus on nucleic acids. Throughout the development of this docking program, efforts were directed toward displaceable water molecules which must be accurately located for optimal docking-based drug discovery. However, although there is a plethora of methods to place water molecules in and around protein structures, there is, to the best of our knowledge, no such fully automated method for nucleic acids, which are significantly more polar and solvated than proteins. We report herein a new method, Splash'Em (Solvation Potential Laid around Statistical Hydration on Entire Macromolecules) developed to place water molecules within the binding cavity of nucleic acids. This fast method was shown to have high agreement with water positions in crystal structures and will therefore provide essential information to medicinal chemists.

Challenges and current status of computational methods for docking small molecules to nucleic acids.

Since the development of the first docking program in 1982, the use of docking-based in silico screening for potentially bioactive molecule discovery has become a common strategy in academia and pharmaceutical industry. Up until recently, application of docking programs has largely focused on drugs binding to proteins. However, with the discovery of promising drug targets in nucleic acids, including RNA riboswitches, DNA G-quadruplexes, and extended repeats in RNA, there has been greater interests in developing drugs for nucleic acids. However, due to major biochemical and physical differences in charges, binding pockets, and solvation, existing docking programs, developed for proteins, face difficulties when adopted directly for nucleic acids. In this review, we cover the current field of in silico docking to nucleic acids, available programs, as well as challenges faced in the field.

Glyceraldehyde 3-phosphate dehydrogenase (GAPDH) is inactivated by S-sulfuration in vitro.

Hydrogen sulfide (H2S) is produced enzymatically by cystathionine ß-synthase (CBS) and cystathionine γ-lyase (CSE), as well as other enzymes in mammalian tissues. These discoveries have led to the crowning of H2S as yet another toxic gas that serves as a gasotransmitter like NO and CO. H2S is thought to exert its biological effects through its reaction with cysteine thiols in proteins, yielding sulfurated thiol (-SSH) derivatives. One of the first proteins shown to be modified by H2S was glyceraldehyde 3-phosphate dehydrogenase (GAPDH) [1] where the S-sulfuration of the active site cysteine (Cys 152) resulted in ~7-fold increase in the activity of the enzyme. In the present study we have attempted to reproduce this result with no success. GAPDH in its reduced, or hydrogen peroxide, or glutathione disulfide, or nitrosonium oxidized forms was reacted with sulfide or polysulfides. Sulfide had no effect on reduced GAPDH activity, while polysulfides inhibited GAPDH to ~42% of control. S-sulfuration of GAPDH occurred at Cys 247 after sulfide treatment, Cys 156 and Cys 247 after polysulfide treatment. No evidence of S-sulfuration at active site Cys 152 was discovered. Both sulfide and polysulfide was able to restore the activity of glutathione disulfide oxidized GAPDH, but not to control untreated levels. Treatment of glutathione disulfide oxidized GAPDH with polysulfide also produced S-sulfuration of Cys 156. Treatment of a C156S mutant of GAPDH with sulfide and polysulfide resulted in S-sulfuration of Cys 152, which also caused a decrease and not an increase in enzymatic activity. Computational chemistry shows S-sulfuration of Cys 156 may affect the position of catalytic Cys 152, raising its pKa by 0.5, which may affect the nucleophilicity of Cys 152. The current study raises significant questions about the reported ability of H2S to activate GAPDH by the sulfuration of its active site thiol, and indicates that polysulfide is a stronger protein S-sulfurating agent than sulfide.