Targeted Protein Degradation: Advances, Challenges, and Prospects for Computational Methods.

The therapeutic approach of targeted protein degradation (TPD) is gaining momentum due to its potentially superior effects compared with protein inhibition. Recent advancements in the biotech and pharmaceutical sectors have led to the development of compounds that are currently in human trials, with some showing promising clinical results. However, the use of computational tools in TPD is still limited, as it has distinct characteristics compared with traditional computational drug design methods. TPD involves creating a ternary structure (protein-degrader-ligase) responsible for the biological function, such as ubiquitination and subsequent proteasomal degradation, which depends on the spatial orientation of the protein of interest (POI) relative to E2-loaded ubiquitin. Modeling this structure necessitates a unique blend of tools initially developed for small molecules (e.g., docking) and biologics (e.g., protein-protein interaction modeling). Additionally, degrader molecules, particularly heterobifunctional degraders, are generally larger than conventional small molecule drugs, leading to challenges in determining drug-like properties like solubility and permeability. Furthermore, the catalytic nature of TPD makes occupancy-based modeling insufficient. TPD consists of multiple interconnected yet distinct steps, such as POI binding, E3 ligase binding, ternary structure interactions, ubiquitination, and degradation, along with traditional small molecule properties. A comprehensive set of tools is needed to address the dynamic nature of the induced proximity ternary complex and its implications for ubiquitination. In this Perspective, we discuss the current state of computational tools for TPD. We start by describing the series of steps involved in the degradation process and the experimental methods used to characterize them. Then, we delve into a detailed analysis of the computational tools employed in TPD. We also present an integrative approach that has proven successful for degrader design and its impact on project decisions. Finally, we examine the future prospects of computational methods in TPD and the areas with the greatest potential for impact.

A Conserved Local Structural Motif Controls the Kinetics of PTP1B Catalysis.

Protein tyrosine phosphatase 1B (PTP1B) is a negative regulator of the insulin and leptin signaling pathways, making it a highly attractive target for the treatment of type II diabetes. For PTP1B to perform its enzymatic function, a loop referred to as the "WPD loop" must transition between open (catalytically incompetent) and closed (catalytically competent) conformations, which have both been resolved by X-ray crystallography. Although prior studies have established this transition as the rate-limiting step for catalysis, the transition mechanism for PTP1B and other PTPs has been unclear. Here we present an atomically detailed model of WPD loop transitions in PTP1B based on unbiased, long-timescale molecular dynamics simulations and weighted ensemble simulations. We found that a specific WPD loop region─the PDFG motif─acted as the key conformational switch, with structural changes to the motif being necessary and sufficient for transitions between long-lived open and closed states of the loop. Simulations starting from the closed state repeatedly visited open states of the loop that quickly closed again unless the infrequent conformational switching of the motif stabilized the open state. The functional importance of the PDFG motif is supported by the fact that it is well conserved across PTPs. Bioinformatic analysis shows that the PDFG motif is also conserved, and adopts two distinct conformations, in deiminases, and the related DFG motif is known to function as a conformational switch in many kinases, suggesting that PDFG-like motifs may control transitions between structurally distinct, long-lived conformational states in multiple protein families.

Molecular Latent Space Simulators for Distributed and Multimolecular Trajectories.

All atom molecular dynamics (MD) simulations offer a powerful tool for molecular modeling, but the short time steps required for numerical stability of the integrator place many interesting molecular events out of reach of unbiased simulations. The popular and powerful Markov state modeling (MSM) approach can extend these time scales by stitching together multiple short discontinuous trajectories into a single long-time kinetic model but necessitates a configurational coarse-graining of the phase space that entails a loss of spatial and temporal resolution and an exponential increase in complexity for multimolecular systems. Latent space simulators (LSS) present an alternative formalism that employs a dynamical, as opposed to configurational, coarse graining comprising three back-to-back learning problems to (i) identify the molecular system's slowest dynamical processes, (ii) propagate the microscopic system dynamics within this slow subspace, and (iii) generatively reconstruct the trajectory of the system within the molecular phase space. A trained LSS model can generate temporally and spatially continuous synthetic molecular trajectories at orders of magnitude lower cost than MD to improve sampling of rare transition events and metastable states to reduce statistical uncertainties in thermodynamic and kinetic observables. In this work, we extend the LSS formalism to short discontinuous training trajectories generated by distributed computing and to multimolecular systems without incurring exponential scaling in computational cost. First, we develop a distributed LSS model over thousands of short simulations of a 264-residue proteolysis-targeting chimera (PROTAC) complex to generate ultralong continuous trajectories that identify metastable states and collective variables to inform PROTAC therapeutic design and optimization. Second, we develop a multimolecular LSS architecture to generate physically realistic ultralong trajectories of DNA oligomers that can undergo both duplex hybridization and hairpin folding. These trajectories retain thermodynamic and kinetic characteristics of the training data while providing increased precision of folding populations and time scales across simulation temperature and ion concentration.

Predicting the structural basis of targeted protein degradation by integrating molecular dynamics simulations with structural mass spectrometry.

Targeted protein degradation (TPD) is a promising approach in drug discovery for degrading proteins implicated in diseases. A key step in this process is the formation of a ternary complex where a heterobifunctional molecule induces proximity of an E3 ligase to a protein of interest (POI), thus facilitating ubiquitin transfer to the POI. In this work, we characterize 3 steps in the TPD process. (1) We simulate the ternary complex formation of SMARCA2 bromodomain and VHL E3 ligase by combining hydrogen-deuterium exchange mass spectrometry with weighted ensemble molecular dynamics (MD). (2) We characterize the conformational heterogeneity of the ternary complex using Hamiltonian replica exchange simulations and small-angle X-ray scattering. (3) We assess the ubiquitination of the POI in the context of the full Cullin-RING Ligase, confirming experimental ubiquitinomics results. Differences in degradation efficiency can be explained by the proximity of lysine residues on the POI relative to ubiquitin.

Accelerated molecular dynamics simulations of the octopamine receptor using GPUs: discovery of an alternate agonist-binding position.

Octopamine receptors (OARs) perform key biological functions in invertebrates, making this class of G-protein coupled receptors (GPCRs) worth considering for insecticide development. However, no crystal structures and very little research exists for OARs. Furthermore, GPCRs are large proteins, are suspended in a lipid bilayer, and are activated on the millisecond timescale, all of which make conventional molecular dynamics (MD) simulations infeasible, even if run on large supercomputers. However, accelerated Molecular Dynamics (aMD) simulations can reduce this timescale to even hundreds of nanoseconds, while running the simulations on graphics processing units (GPUs) would enable even small clusters of GPUs to have processing power equivalent to hundreds of CPUs. Our results show that aMD simulations run on GPUs can successfully obtain the active and inactive state conformations of a GPCR on this reduced timescale. Furthermore, we discovered a potential alternate active-state agonist-binding position in the octopamine receptor which has yet to be observed and may be a novel GPCR agonist-binding position. These results demonstrate that a complex biological system with an activation process on the millisecond timescale can be successfully simulated on the nanosecond timescale using a simple computing system consisting of a small number of GPUs. Proteins 2016; 84:1480-1489. © 2016 Wiley Periodicals, Inc.

A comparison of weighted ensemble and Markov state model methodologies.

Computation of reaction rates and elucidation of reaction mechanisms are two of the main goals of molecular dynamics (MD) and related simulation methods. Since it is time consuming to study reaction mechanisms over long time scales using brute force MD simulations, two ensemble methods, Markov State Models (MSMs) and Weighted Ensemble (WE), have been proposed to accelerate the procedure. Both approaches require clustering of microscopic configurations into networks of "macro-states" for different purposes. MSMs model a discretization of the original dynamics on the macro-states. Accuracy of the model significantly relies on the boundaries of macro-states. On the other hand, WE uses macro-states to formulate a resampling procedure that kills and splits MD simulations for achieving better efficiency of sampling. Comparing to MSMs, accuracy of WE rate predictions is less sensitive to the definition of macro-states. Rigorous numerical experiments using alanine dipeptide and penta-alanine support our analyses. It is shown that MSMs introduce significant biases in the computation of reaction rates, which depend on the boundaries of macro-states, and Accelerated Weighted Ensemble (AWE), a formulation of weighted ensemble that uses the notion of colors to compute fluxes, has reliable flux estimation on varying definitions of macro-states. Our results suggest that whereas MSMs provide a good idea of the metastable sets and visualization of overall dynamics, AWE provides reliable rate estimations requiring less efforts on defining macro-states on the high dimensional conformational space.

Evaluation of conformational changes in diabetes-associated mutation in insulin a chain: a molecular dynamics study.

Insulin plays a central role in the regulation of metabolism in humans. Mutations in the insulin gene can impair the folding of its precursor protein, proinsulin, and cause permanent neonatal-onset diabetes mellitus known as Mutant INS-gene induced Diabetes of Youth (MIDY) with insulin deficiency. To gain insights into the molecular basis of this diabetes-associated mutation, we perform molecular dynamics simulations in wild-type and mutant (Cys(A7) to Tyr or C(A7)Y) insulin A chain in aqueous solutions. The C(A7)Y mutation is one of the identified mutations that impairs the protein folding by substituting the cysteine residue which is required for the disulfide bond formation. A comparative analysis reveals structural differences between the wild-type and the mutant conformations. The analyzed mutant insulin A chain forms a metastable state with major effects on its N-terminal region. This suggests that MIDY mutant involves formation of a partially folded intermediate with conformational change in N-terminal region in A chain that generates flexible N-terminal domain. This may lead to the abnormal interactions with other proinsulins in the aggregation process.

Characterization of the Anopheles gambiae octopamine receptor and discovery of potential agonists and antagonists using a combined computational-experimental approach.

BACKGROUND: Octopamine receptors (OARs) perform key functions in the biological pathways of primarily invertebrates, making this class of G-protein coupled receptors (GPCRs) a potentially good target for insecticides. However, the lack of structural and experimental data for this insect-essential GPCR family has promoted the development of homology models that are good representations of their biological equivalents for in silico screening of small molecules. METHODS: Two Anopheles gambiae OARs were cloned, analysed and functionally characterized using a heterologous cell reporter system. Four antagonist- and four agonist-binding homology models were generated and virtually screened by docking against compounds obtained from the ZINC database. Resulting compounds from the virtual screen were tested experimentally using an in vitro reporter assay and in a mosquito larvicide bioassay. RESULTS: Six An. gambiae OAR/tyramine receptor genes were identified. Phylogenetic analysis revealed that the OAR (AGAP000045) that encodes two open reading frames is an α-adrenergic-like receptor. Both splice variants signal through cAMP and calcium. Mutagenesis analysis revealed that D100 in the TM3 region and S206 and S210 in the TM5 region are important to the activation of the GPCR. Some 2,150 compounds from the virtual screen were structurally analysed and 70 compounds were experimentally tested against AgOAR45B expressed in the GloResponse™CRE-luc2P HEK293 reporter cell line, revealing 21 antagonists, 17 weak antagonists, 2 agonists, and 5 weak agonists. CONCLUSION: Reported here is the functional characterization of two An. gambiae OARs and the discovery of new OAR agonists and antagonists based on virtual screening and molecular dynamics simulations. Four compounds were identified that had activity in a mosquito larva bioassay, three of which are imidazole derivatives. This combined computational and experimental approach is appropriate for the discovery of new and effective insecticides.

AWE-WQ: fast-forwarding molecular dynamics using the accelerated weighted ensemble.

A limitation of traditional molecular dynamics (MD) is that reaction rates are difficult to compute. This is due to the rarity of observing transitions between metastable states since high energy barriers trap the system in these states. Recently the weighted ensemble (WE) family of methods have emerged which can flexibly and efficiently sample conformational space without being trapped and allow calculation of unbiased rates. However, while WE can sample correctly and efficiently, a scalable implementation applicable to interesting biomolecular systems is not available. We provide here a GPLv2 implementation called AWE-WQ of a WE algorithm using the master/worker distributed computing WorkQueue (WQ) framework. AWE-WQ is scalable to thousands of nodes and supports dynamic allocation of computer resources, heterogeneous resource usage (such as central processing units (CPU) and graphical processing units (GPUs) concurrently), seamless heterogeneous cluster usage (i.e., campus grids and cloud providers), and support for arbitrary MD codes such as GROMACS, while ensuring that all statistics are unbiased. We applied AWE-WQ to a 34 residue protein which simulated 1.5 ms over 8 months with peak aggregate performance of 1000 ns/h. Comparison was done with a 200 µs simulation collected on a GPU over a similar timespan. The folding and unfolded rates were of comparable accuracy.

Identification of novel arthropod vector G protein-coupled receptors.

BACKGROUND: The control of vector-borne diseases, such as malaria, dengue fever, and typhus fever is often achieved with the use of insecticides. Unfortunately, insecticide resistance is becoming common among different vector species. There are currently no chemical alternatives to these insecticides because new human-safe classes of molecules have yet to be brought to the vector-control market. The identification of novel targets offer opportunities for rational design of new chemistries to control vector populations. One target family, G protein-coupled receptors (GPCRs), has remained relatively under explored in terms of insecticide development. METHODS: A novel classifier, Ensemble*, for vector GPCRs was developed. Ensemble* was validated and compared to existing classifiers using a set of all known GPCRs from Aedes aegypti, Anopheles gambiae, Apis Mellifera, Drosophila melanogaster, Homo sapiens, and Pediculus humanus. Predictions for unidentified sequences from Ae. aegypti, An. gambiae, and Pe. humanus were validated. Quantitative RT-PCR expression analysis was performed on previously-known and newly discovered Ae. aegypti GPCR genes. RESULTS: We present a new analysis of GPCRs in the genomes of Ae, aegypti, a vector of dengue fever, An. gambiae, a primary vector of Plasmodium falciparum that causes malaria, and Pe. humanus, a vector of epidemic typhus fever, using a novel GPCR classifier, Ensemble*, designed for insect vector species. We identified 30 additional putative GPCRs, 19 of which we validated. Expression of the newly discovered Ae. aegypti GPCR genes was confirmed via quantitative RT-PCR. CONCLUSION: A novel GPCR classifier for insect vectors, Ensemble*, was developed and GPCR predictions were validated. Ensemble* and the validation pipeline were applied to the genomes of three insect vectors (Ae. aegypti, An. gambiae, and Pe. humanus), resulting in the identification of 52 GPCRs not previously identified, of which 11 are predicted GPCRs, and 19 are predicted and confirmed GPCRs.

Long Timestep Molecular Dynamics on the Graphical Processing Unit.

Molecular dynamics (MD) simulations now play a key role in many areas of theoretical chemistry, biology, physics, and materials science. In many cases, such calculations are significantly limited by the massive amount of computer time needed to perform calculations of interest. Herein, we present Long Timestep Molecular Dynamics (LTMD), a method to significantly speed MD simulations. In particular, we discuss new methods to calculate the needed terms in LTMD as well as issues germane to a GPU implementation. The resulting code, implemented in the OpenMM MD library, can achieve a significant 6-fold speed increase, leading to MD simulations on the order of 5 µs/day using implicit solvent models.

Folding Proteins at 500 ns/hour with Work Queue.

Molecular modeling is a field that traditionally has large computational costs. Until recently, most simulation techniques relied on long trajectories, which inherently have poor scalability. A new class of methods is proposed that requires only a large number of short calculations, and for which minimal communication between computer nodes is required. We considered one of the more accurate variants called Accelerated Weighted Ensemble Dynamics (AWE) and for which distributed computing can be made efficient. We implemented AWE using the Work Queue framework for task management and applied it to an all atom protein model (Fip35 WW domain). We can run with excellent scalability by simultaneously utilizing heterogeneous resources from multiple computing platforms such as clouds (Amazon EC2, Microsoft Azure), dedicated clusters, grids, on multiple architectures (CPU/GPU, 32/64bit), and in a dynamic environment in which processes are regularly added or removed from the pool. This has allowed us to achieve an aggregate sampling rate of over 500 ns/hour. As a comparison, a single process typically achieves 0.1 ns/hour.

Modeling conformational ensembles of slow functional motions in Pin1-WW.

Protein-protein interactions are often mediated by flexible loops that experience conformational dynamics on the microsecond to millisecond time scales. NMR relaxation studies can map these dynamics. However, defining the network of inter-converting conformers that underlie the relaxation data remains generally challenging. Here, we combine NMR relaxation experiments with simulation to visualize networks of inter-converting conformers. We demonstrate our approach with the apo Pin1-WW domain, for which NMR has revealed conformational dynamics of a flexible loop in the millisecond range. We sample and cluster the free energy landscape using Markov State Models (MSM) with major and minor exchange states with high correlation with the NMR relaxation data and low NOE violations. These MSM are hierarchical ensembles of slowly interconverting, metastable macrostates and rapidly interconverting microstates. We found a low population state that consists primarily of holo-like conformations and is a "hub" visited by most pathways between macrostates. These results suggest that conformational equilibria between holo-like and alternative conformers pre-exist in the intrinsic dynamics of apo Pin1-WW. Analysis using MutInf, a mutual information method for quantifying correlated motions, reveals that WW dynamics not only play a role in substrate recognition, but also may help couple the substrate binding site on the WW domain to the one on the catalytic domain. Our work represents an important step towards building networks of inter-converting conformational states and is generally applicable.

MDLab: a molecular dynamics simulation prototyping environment.

Molecular dynamics (MD) simulation involves solving Newton's equations of motion for a system of atoms, by calculating forces and updating atomic positions and velocities over a timestep Deltat. Despite the large amount of computing power currently available, the timescale of MD simulations is limited by both the small timestep required for propagation, and the expensive algorithm for computing pairwise forces. These issues are currently addressed through the development of efficient simulation methods, some of which make acceptable approximations and as a result can afford larger timesteps. We present MDLab, a development environment for MD simulations built with Python which facilitates prototyping, testing, and debugging of these methods. MDLab provides constructs which allow the development of propagators, force calculators, and high level sampling protocols that run several instances of molecular dynamics. For computationally demanding sampling protocols which require testing on large biomolecules, MDL includes an interface to the OpenMM libraries of Friedrichs et al. which execute on graphical processing units (GPUs) and achieve considerable speedup over execution on the CPU. As an example of an interesting high level method developed in MDLab, we present a parallel implementation of the On-The-Fly string method of Maragliano and Vanden-Eijnden. MDLab is available at http://mdlab.sourceforge.net.

A separable shadow Hamiltonian hybrid Monte Carlo method.

Hybrid Monte Carlo (HMC) is a rigorous sampling method that uses molecular dynamics (MD) as a global Monte Carlo move. The acceptance rate of HMC decays exponentially with system size. The shadow hybrid Monte Carlo (SHMC) was previously introduced to reduce this performance degradation by sampling instead from the shadow Hamiltonian defined for MD when using a symplectic integrator. SHMC's performance is limited by the need to generate momenta for the MD step from a nonseparable shadow Hamiltonian. We introduce the separable shadow Hamiltonian hybrid Monte Carlo (S2HMC) method based on a formulation of the leapfrog/Verlet integrator that corresponds to a separable shadow Hamiltonian, which allows efficient generation of momenta. S2HMC gives the acceptance rate of a fourth order integrator at the cost of a second-order integrator. Through numerical experiments we show that S2HMC consistently gives a speedup greater than two over HMC for systems with more than 4000 atoms for the same variance. By comparison, SHMC gave a maximum speedup of only 1.6 over HMC. S2HMC has the additional advantage of not requiring any user parameters beyond those of HMC. S2HMC is available in the program PROTOMOL 2.1. A Python version, adequate for didactic purposes, is also in MDL (http://mdlab.sourceforge.net/s2hmc).

Inferring protein-protein interactions from multiple protein domain combinations.

The ever accumulating wealth of knowledge about protein interactions and the domain architecture of involved proteins in different organisms offers ways to understand the intricate interplay between interactome and proteome. Ultimately, the combination of these sources of information will allow the prediction of interactions among proteins where only domain composition is known. Based on the currently available protein-protein interaction and domain data of Saccharomyces cerevisiae and Drosophila melanogaster we introduce a novel method, Maximum Specificity Set Cover (MSSC), to predict potential protein-protein interactions. Utilizing interactions and domain architectures of domains as training sets, this algorithm employs a set cover approach to partition domain pairs, which allows the explanation of the underlying protein interaction to the largest degree of specificity. While MSSC in its basic version only considers domain pairs as the driving force between interactions, we also modified the algorithm to account for combinations of more than two domains that govern a protein-protein interaction. This approach allows us to predict the previously unknown protein-protein interactions in S. cerevisiae and D. melanogaster, with a degree of sensitivity and specificity that clearly outscores other approaches. As a proof of concept we also observe high levels of co-expression and decreasing GO distances between interacting proteins. Although our results are very encouraging, we observe that the quality of predictions significantly depends on the quality of interactions, which were utilized as the training set of the algorithm. The algorithm is part of a Web portal available at http://ppi.cse.nd.edu .

Normal mode partitioning of Langevin dynamics for biomolecules.

We propose a novel normal mode multiple time stepping Langevin dynamics integrator called NML. The aim is to approximate the kinetics or thermodynamics of a biomolecule by a reduced model based on a normal mode decomposition of the dynamical space. Our basis set uses the eigenvectors of a mass reweighted Hessian matrix calculated with a biomolecular force field. This particular choice has the advantage of an ordering according to the eigenvalues, which have a physical meaning of being the square of the mode frequency. Low frequency eigenvalues correspond to more collective motions, whereas the highest frequency eigenvalues are the limiting factor for the stability of the integrator. In NML, the higher frequency modes are overdamped and relaxed near their energy minimum while respecting the subspace of low frequency dynamical modes. Our numerical results confirm that both sampling and rates are conserved for an implicitly solvated alanine dipeptide model, with only 30% of the modes propagated, when compared to the full model. For implicitly solvated systems, NML gives a twofold improvement in efficiency over plain Langevin dynamics for sampling a small 22 atom (alanine dipeptide) model and in excess of an order of magnitude for sampling an 882 atom (bovine pancreatic trypsin inhibitor) model, with good scaling with system size subject to the number of modes propagated. NML has been implemented in the open source software PROTOMOL.

A parallel implementation of the Cellular Potts Model for simulation of cell-based morphogenesis.

The Cellular Potts Model (CPM) has been used in a wide variety of biological simulations. However, most current CPM implementations use a sequential modified Metropolis algorithm which restricts the size of simulations. In this paper we present a parallel CPM algorithm for simulations of morphogenesis, which includes cell-cell adhesion, a cell volume constraint, and cell haptotaxis. The algorithm uses appropriate data structures and checkerboard subgrids for parallelization. Communication and updating algorithms synchronize properties of cells simulated on different processor nodes. Tests show that the parallel algorithm has good scalability, permitting large-scale simulations of cell morphogenesis (10(7) or more cells) and broadening the scope of CPM applications. The new algorithm satisfies the balance condition, which is sufficient for convergence of the underlying Markov chain.

Accelerating the replica exchange method through an efficient all-pairs exchange.

The authors accelerate the replica exchange method through an efficient all-pairs replica exchange. A proof of detailed balance is shown along with an analytical estimate of the enhanced exchange efficiency. The new method provides asymptotically four fold speedup of conformation traversal for replica counts of 8 and larger with typical exchange rates. Experimental tests using the blocked alanine dipeptide demonstrate the method's correctness and show an approximate sampling efficiency improvement of 100% according to potential energy cumulative averages and an ergodic measure. An explicitly solvated PIN1 WW domain system of 4958 atoms is sampled using our new method, yielding a cluster sampling rate almost twice that of the single exchange near neighbor implementation. Computational software and scripts along with input and output data sets are available at.

Predicting protein-protein interactions from protein domains using a set cover approach.

One goal of contemporary proteome research is the elucidation of cellular protein interactions. Based on currently available protein-protein interaction and domain data, we introduce a novel method, Maximum Specificity Set Cover (MSSC), for the prediction of protein-protein interactions. In our approach, we map the relationship between interactions of proteins and their corresponding domain architectures to a generalized weighted set cover problem. The application of a greedy algorithm provides sets of domain interactions which explain the presence of protein interactions to the largest degree of specificity. Utilizing domain and protein interaction data of S. cerevisiae, MSSC enables prediction of previously unknown protein interactions, links that are well supported by a high tendency of coexpression and functional homogeneity of the corresponding proteins. Focusing on concrete examples, we show that MSSC reliably predicts protein interactions in well-studied molecular systems, such as the 26S proteasome and RNA polymerase II of S. cerevisiae. We also show that the quality of the predictions is comparable to the Maximum Likelihood Estimation while MSSC is faster. This new algorithm and all data sets used are accessible through a Web portal at http://ppi.cse.nd.edu.