The membrane-embedded γ-secretase complex processively cleaves within the transmembrane domain of amyloid precursor protein (APP) to produce 37-to-43-residue amyloid β-peptides (Aβ) of Alzheimer’s disease (AD). Despite its importance in pathogenesis, the mechanism of processive proteolysis by γ-secretase remains poorly understood. Here, mass spectrometry and Western blotting were used to quantify the efficiency of tripeptide trimming of wild-type (WT) and familial AD (FAD) mutant Aβ49. In comparison to WT Aβ49, the efficiency of tripeptide trimming was similar for the I45F, A42T, and V46F Aβ49 FAD mutants but substantially diminished for the I45T and T48P mutants. In parallel with biochemical experiments, all-atom simulations using a novel peptide Gaussian accelerated molecular dynamics (Pep-GaMD) method were applied to investigate the tripeptide trimming of Aβ49 by γ-secretase. The starting structure was the active γ-secretase bound to Aβ49 and APP intracellular domain (AICD), as generated from our previous study that captured the activation of γ-secretase for the initial endoproteolytic cleavage of APP (BhattaraiA. Bhattarai, A. ACS Cent. Sci.20206969983). Pep-GaMD simulations captured remarkable structural rearrangements of both the enzyme and substrate, in which hydrogen-bonded catalytic aspartates and water became poised for tripeptide trimming of Aβ49 to Aβ46. These structural changes required a positively charged N-terminus of endoproteolytic coproduct AICD, which could dissociate during conformational rearrangements of the protease and Aβ49. The simulation findings were highly consistent with biochemical experimental data. Taken together, our complementary biochemical experiments and Pep-GaMD simulations have enabled elucidation of the mechanism of tripeptide trimming of Aβ49 by γ-secretase.

We introduce a Gaussian-accelerated molecular dynamics (GaMD), deep learning (DL), and free energy profiling workflow (GLOW) to predict molecular determinants and map free energy landscapes of biomolecules. All-atom GaMD-enhanced sampling simulations are first performed on biomolecules of interest. Structural contact maps are then calculated from GaMD simulation frames and transformed into images for building DL models using a convolutional neural network. Important structural contacts are further determined from DL models of attention maps of the structural contact gradients, which allow us to identify the system reaction coordinates. Finally, free energy profiles are calculated for the selected reaction coordinates through energetic reweighting of the GaMD simulations. We have also successfully demonstrated GLOW for the characterization of activation and allosteric modulation of a G protein-coupled receptor, using the adenosine A 1 receptor (A 1 AR) as a model system. GLOW findings are highly consistent with previous experimental and computational studies of the A 1 AR, while also providing further mechanistic insights into the receptor function. In summary, GLOW provides a systematic approach to mapping free energy landscapes of biomolecules. The GLOW workflow and its user manual can be downloaded at http://miaolab.org/GLOW.

Ligand binding thermodynamics and kinetics are critical parameters for drug design. However, it has proven challenging to efficiently predict ligand binding thermodynamics and kinetics from molecular simulations due to limited simulation timescales. Protein dynamics, especially in the ligand binding pocket, often plays an important role in ligand binding. Based on our previously developed Ligand Gaussian accelerated molecular dynamics (LiGaMD), here we present LiGaMD2 in which a selective boost potential was applied to both the ligand and protein residues in the binding pocket to improve sampling of ligand binding and dissociation. To validate the performance of LiGaMD2, the T4 lysozyme (T4L) mutants with open and closed pockets bound by different ligands were chosen as model systems. LiGaMD2 could efficiently capture repetitive ligand dissociation and binding within microsecond simulations of all T4L systems. The obtained ligand binding kinetic rates and free energies agreed well with available experimental values and previous modeling results. Therefore, LiGaMD2 provides an improved approach to sample opening of closed protein pockets for ligand dissociation and binding, thereby allowing for efficient calculations of ligand binding thermodynamics and kinetics.

Protein−protein interactions (PPIs) play key roles in many fundamental biological processes such as cellular signaling and immune responses. However, it has proven challenging to simulate repetitive protein association and dissociation in order to calculate binding free energies and kinetics of PPIs due to long biological timescales and complex protein dynamics. To address this challenge, we have developed a new computational approach to all-atom simulations of PPIs based on a robust Gaussian accelerated molecular dynamics (GaMD) technique. The method, termed "PPI-GaMD", selectively boosts interaction potential energy between protein partners to facilitate their slow dissociation. Meanwhile, another boost potential is applied to the remaining potential energy of the entire system to effectively model the protein's flexibility and rebinding. PPI-GaMD has been demonstrated on a model system of the ribonuclease barnase interactions with its inhibitor barstar. Six independent 2 μs PPI-GaMD simulations have captured repetitive barstar dissociation and rebinding events, which enable calculations of the protein binding thermodynamics and kinetics simultaneously. The calculated binding free energies and kinetic rate constants agree well with the experimental data. Furthermore, PPI-GaMD simulations have provided mechanistic insights into barstar binding to barnase, which involves long-range electrostatic interactions and multiple binding pathways, being consistent with previous experimental and computational findings of this model system. In summary, PPI-GaMD provides a highly efficient and easy-to-use approach for binding free energy and kinetics calculations of PPIs.

Gaussian accelerated molecular dynamics (GaMD) is a computational technique that provides both unconstrained enhanced sampling and free energy calculations of biomolecules. Here, we present the implementation of GaMD in the OpenMM simulation package and validate it on model systems of alanine dipeptide and RNA folding. For alanine dipeptide, 30 ns GaMD production simulations reproduced free energy profiles of 1000 ns conventional molecular dynamics (cMD) simulations. In addition, GaMD simulations captured the folding pathways of three hyperstable RNA tetraloops (UUCG, GCAA, and CUUG) and binding of the rbt203 ligand to the HIV-1 Tar RNA, both of which involved critical electrostatic interactions such as hydrogen bonding and base stacking. Together with previous implementations, GaMD in OpenMM will allow for wider applications in simulations of proteins, RNA, and other biomolecules.

Biomolecular binding kinetics including the association (k on ) and dissociation (k of f ) rates are critical parameters for therapeutic design of smallmolecule drugs, peptides, and antibodies. Notably, the drug molecule residence time or dissociation rate has been shown to correlate with their efficacies better than binding affinities. A wide range of modeling approaches including quantitative structure-kinetic relationship models, Molecular Dynamics simulations, enhanced sampling, and Machine Learning has been developed to explore biomolecular binding and dissociation mechanisms and predict binding kinetic rates. Here, we review recent advances in computational modeling of biomolecular binding kinetics, with an outlook for future improvements.

Biomolecular recognition such as binding of small molecules, nucleic acids, peptides and proteins to their target receptors plays key roles in cellular function and has been targeted for therapeutic drug design. Molecular dynamics (MD) is a computational approach to analyze these binding processes at an atomistic level, which provides valuable understandings of the mechanisms of biomolecular recognition. However, the rather slow biomolecular binding events often present challenges for conventional MD (cMD), due to limited simulation timescales (typically over hundreds of nanoseconds to tens of microseconds). In this regard, enhanced sampling methods, particularly accelerated MD (aMD), have proven useful to bridge the gap and enable all-atom simulations of biomolecular binding events. Here, we will review the recent method developments of Gaussian aMD (GaMD), ligand GaMD (LiGaMD) and peptide GaMD (Pep-GaMD), which have greatly expanded our capabilities to simulate biomolecular binding processes. Spontaneous binding of various biomolecules to their receptors has been successfully simulated by GaMD. Microsecond LiGaMD and Pep-GaMD simulations have captured repetitive binding and dissociation of small-molecule ligands and highly flexible peptides, and thus enabled ligand/peptide binding thermodynamics and kinetics calculations. We will also present relevant application studies in simulations of important drug targets and future perspectives for rational computer-aided drug design.

The RNA-binding protein Hu antigen R (HuR) is a post-transcriptional regulator critical in several types of diseases, including cancer, making it a promising therapeutic target. We have identified small-molecule inhibitors of HuR through a screening approach used in combination with fragment analysis. A total of 36 new compounds originating from fragment linking or structural optimization were studied to establish structure–activity relationships in the set. Two top inhibitors, 1c and 7c, were further validated by binding assays and cellular functional assays. Both block HuR function by directly binding to the RNA-binding pocket, inhibit cancer cell growth dependence of HuR, and suppress cancer cell invasion. Intraperitoneal administration of inhibitor 1c inhibits tumor growth as a single agent and shows a synergistic effect in combination with chemotherapy docetaxel in breast cancer xenograft models. Mechanistically, 1c interferes with the HuR–TGFB/THBS1 axis.