PaCS-Toolkit: Optimized Software Utilities for Parallel Cascade Selection Molecular Dynamics (PaCS-MD) Simulations and Subsequent Analyses.

Parallel cascade selection molecular dynamics (PaCS-MD) is an enhanced conformational sampling method conducted as a "repetition of time leaps in parallel worlds", comprising cycles of multiple molecular dynamics (MD) simulations performed in parallel and selection of the initial structures of MDs for the next cycle. We developed PaCS-Toolkit, an optimized software utility enabling the use of different MD software and trajectory analysis tools to facilitate the execution of the PaCS-MD simulation and analyze the obtained trajectories, including the preparation for the subsequent construction of the Markov state model. PaCS-Toolkit is coded with Python, is compatible with various computing environments, and allows for easy customization by editing the configuration file and specifying the MD software and analysis tools to be used. We present the software design of PaCS-Toolkit and demonstrate applications of PaCS-MD variations: original targeted PaCS-MD to peptide folding; rmsdPaCS-MD to protein domain motion; and dissociation PaCS-MD to ligand dissociation from adenosine A2A receptor.

Engineering of an in-cell protein crystal for fastening a metastable conformation of a target miniprotein.

Protein crystals can be utilized as porous scaffolds to capture exogenous molecules. Immobilization of target proteins using protein crystals is expected to facilitate X-ray structure analysis of proteins that are difficult to be crystallized. One of the advantages of scaffold-assisted structure determination is the analysis of metastable structures that are not observed in solution. However, efforts to fix target proteins within the pores of scaffold protein crystals have been limited due to the lack of strategies to control protein-protein interactions formed in the crystals. In this study, we analyze the metastable structure of the miniprotein, CLN025, which forms a ß-hairpin structure in solution, using a polyhedra crystal (PhC), an in-cell protein crystal. CLN025 is successfully fixed within the PhC scaffold by replacing the original loop region. X-ray crystal structure analysis and molecular dynamics (MD) simulation reveal that CLN025 is fixed as a helical structure in a metastable state by non-covalent interactions in the scaffold crystal. These results indicate that modulation of intermolecular interactions can trap various protein conformations in the engineered PhC and provides a new strategy for scaffold-assisted structure determination.

Structure of MotA, a flagellar stator protein, from hyperthermophile.

Many motile bacteria swim and swarm toward favorable environments using the flagellum, which is rotated by a motor embedded in the inner membrane. The motor is composed of the rotor and the stator, and the motor torque is generated by the change of the interaction between the rotor and the stator induced by the ion flow through the stator. A stator unit consists of two types of membrane proteins termed A and B. Recent cryo-EM studies on the stators from mesophiles revealed that the stator consists of five A and two B subunits, whereas the low-resolution EM analysis showed that purified hyperthermophilic MotA forms a tetramer. To clarify the assembly formation and factors enhancing thermostability of the hyperthermophilic stator, we determined the cryo-EM structure of MotA from Aquifex aeolicus (Aa-MotA), a hyperthermophilic bacterium, at 3.42 Å resolution. Aa-MotA forms a pentamer with pseudo C5 symmetry. A simulated model of the Aa-MotA5MotB2 stator complex resembles the structures of mesophilic stator complexes, suggesting that Aa-MotA can assemble into a pentamer equivalent to the stator complex without MotB. The distribution of hydrophobic residues of MotA pentamers suggests that the extremely hydrophobic nature in the subunit boundary and the transmembrane region is a key factor to stabilize hyperthermophilic Aa-MotA.

Inhibition of the hexamerization of SARS-CoV-2 endoribonuclease and modeling of RNA structures bound to the hexamer.

Non-structural protein 15 (Nsp15) of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) forms a homo hexamer and functions as an endoribonuclease. Here, we propose that Nsp15 activity may be inhibited by preventing its hexamerization through drug binding. We first explored the stable conformation of the Nsp15 monomer as the global free energy minimum conformation in the free energy landscape using a combination of parallel cascade selection molecular dynamics (PaCS-MD) and the Markov state model (MSM), and found that the Nsp15 monomer forms a more open conformation with larger druggable pockets on the surface. Targeting the pockets with high druggability scores, we conducted ligand docking and identified compounds that tightly bind to the Nsp15 monomer. The top poses with Nsp15 were subjected to binding free energy calculations by dissociation PaCS-MD and MSM (dPaCS-MD/MSM), indicating the stability of the complexes. One of the identified pockets, which is distinctively bound by inosine analogues, may be an alternative binding site to stabilize viral RNA binding and/or an alternative catalytic site. We constructed a stable RNA structure model bound to both UTP and alternative binding sites, providing a reasonable proposed model of the Nsp15/RNA complex.

Using molecular dynamics simulations to prioritize and understand AI-generated cell penetrating peptides.

Cell-penetrating peptides have important therapeutic applications in drug delivery, but the variety of known cell-penetrating peptides is still limited. With a promise to accelerate peptide development, artificial intelligence (AI) techniques including deep generative models are currently in spotlight. Scientists, however, are often overwhelmed by an excessive number of unannotated sequences generated by AI and find it difficult to obtain insights to prioritize them for experimental validation. To avoid this pitfall, we leverage molecular dynamics (MD) simulations to obtain mechanistic information to prioritize and understand AI-generated peptides. A mechanistic score of permeability is computed from five steered MD simulations starting from different initial structures predicted by homology modelling. To compensate for variability of predicted structures, the score is computed with sample variance penalization so that a peptide with consistent behaviour is highly evaluated. Our computational pipeline involving deep learning, homology modelling, MD simulations and synthesizability assessment generated 24 novel peptide sequences. The top-scoring peptide showed a consistent pattern of conformational change in all simulations regardless of initial structures. As a result of wet-lab-experiments, our peptide showed better permeability and weaker toxicity in comparison to a clinically used peptide, TAT. Our result demonstrates how MD simulations can support de novo peptide design by providing mechanistic information supplementing statistical inference.

Delineating the conformational landscape of the adenosine A2A receptor during G protein coupling.

G-protein-coupled receptors (GPCRs) represent a ubiquitous membrane protein family and are important drug targets. Their diverse signaling pathways are driven by complex pharmacology arising from a conformational ensemble rarely captured by structural methods. Here, fluorine nuclear magnetic resonance spectroscopy (19F NMR) is used to delineate key functional states of the adenosine A2A receptor (A2AR) complexed with heterotrimeric G protein (Gαsß1γ2) in a phospholipid membrane milieu. Analysis of A2AR spectra as a function of ligand, G protein, and nucleotide identifies an ensemble represented by inactive states, a G-protein-bound activation intermediate, and distinct nucleotide-free states associated with either partial- or full-agonist-driven activation. The Gßγ subunit is found to be critical in facilitating ligand-dependent allosteric transmission, as shown by 19F NMR, biochemical, and computational studies. The results provide a mechanistic basis for understanding basal signaling, efficacy, precoupling, and allostery in GPCRs.

Generating Ampicillin-Level Antimicrobial Peptides with Activity-Aware Generative Adversarial Networks.

Antimicrobial peptides are a potential solution to the threat of multidrug-resistant bacterial pathogens. Recently, deep generative models including generative adversarial networks (GANs) have been shown to be capable of designing new antimicrobial peptides. Intuitively, a GAN controls the probability distribution of generated sequences to cover active peptides as much as possible. This paper presents a peptide-specialized model called PepGAN that takes the balance between covering active peptides and dodging nonactive peptides. As a result, PepGAN has superior statistical fidelity with respect to physicochemical descriptors including charge, hydrophobicity, and weight. Top six peptides were synthesized, and one of them was confirmed to be highly antimicrobial. The minimum inhibitory concentration was 3.1 µg/mL, indicating that the peptide is twice as strong as ampicillin.

Edge expansion parallel cascade selection molecular dynamics simulation for investigating large-amplitude collective motions of proteins.

We propose edge expansion parallel cascade selection molecular dynamics (eePaCS-MD) as an efficient adaptive conformational sampling method to investigate the large-amplitude motions of proteins without prior knowledge of the conformational transitions. In this method, multiple independent MD simulations are iteratively conducted from initial structures randomly selected from the vertices of a multi-dimensional principal component subspace. This subspace is defined by an ensemble of protein conformations sampled during previous cycles of eePaCS-MD. The edges and vertices of the conformational subspace are determined by solving the "convex hull problem." The sampling efficiency of eePaCS-MD is achieved by intensively repeating MD simulations from the vertex structures, which increases the probability of rare event occurrence to explore new large-amplitude collective motions. The conformational sampling efficiency of eePaCS-MD was assessed by investigating the open-close transitions of glutamine binding protein, maltose/maltodextrin binding protein, and adenylate kinase and comparing the results to those obtained using related methods. In all cases, the open-close transitions were simulated in ∼10 ns of simulation time or less, offering 1-3 orders of magnitude shorter simulation time compared to conventional MD. Furthermore, we show that the combination of eePaCS-MD and accelerated MD can further enhance conformational sampling efficiency, which reduced the total computational cost of observing the open-close transitions by at most 36%.

Kinetic Selection and Relaxation of the Intrinsically Disordered Region of a Protein upon Binding.

Here, we investigate the association and dissociation mechanisms of a typical intrinsically disordered region (IDR), transcriptional activation subdomain of tumor suppressor protein p53 (TAD-p53), with murine double-minute clone 2 protein (MDM2). Using a combination of cycles of association and dissociation parallel cascade molecular dynamics, multiple standard molecular dynamics (MD), and the Markov state model, we were successful in obtaining the lowest free energy structure of the MDM2/TAD-p53 complex as the structure closest to the crystal structure without prior knowledge of the crystal structure. This method also reproduced the experimentally measured standard binding free energy, and the association and dissociation rate constants, requiring only an accumulated MD simulation cost of 11.675 µs even though that actual dissociation occurs on the order of seconds. We identified few complex intermediates with similar free energies; yet TAD-p53 first binds MDM2 as the second lowest free energy intermediate kinetically with >90% of the flux, adopting a conformation similar to that of one of these few intermediates in its monomeric state. Even though the mechanism of the first step has a conformational-selection-type aspect, the second step shows induced-fit-like features and occurs as concomitant dehydration of the interface, side-chain π-π stacking, and main-chain hydrogen-bond formation to complete binding as an α-helix. In addition, dehydration is a key process for the final relaxation process around the complex interface. These results demonstrate that TAD-p53 kinetically selects its initial binding form and then relaxes to complete the binding.

Enhancing Biomolecular Sampling with Reinforcement Learning: A Tree Search Molecular Dynamics Simulation Method.

This paper proposes a novel molecular simulation method, called tree search molecular dynamics (TS-MD), to accelerate the sampling of conformational transition pathways, which require considerable computation. In TS-MD, a tree search algorithm, called upper confidence bounds for trees, which is a type of reinforcement learning algorithm, is applied to sample the transition pathway. By learning from the results of the previous simulations, TS-MD efficiently searches conformational space and avoids being trapped in local stable structures. TS-MD exhibits better performance than parallel cascade selection molecular dynamics, which is one of the state-of-the-art methods, for the folding of miniproteins, Chignolin and Trp-cage, in explicit water.

Dissociation Process of a MDM2/p53 Complex Investigated by Parallel Cascade Selection Molecular Dynamics and the Markov State Model.

Recently, we efficiently generated dissociation pathways of a protein-ligand complex without applying force bias with parallel cascade selection molecular dynamics (PaCS-MD) and showed that PaCS-MD in combination with the Markov state model (MSM) yielded a binding free energy comparable to experimental values. In this work, we applied the same procedure to a complex of MDM2 protein and the transactivation domain of p53 protein (TAD-p53), the latter of which is known to be very flexible in the unbound state. Using 30 independent MD simulations in PaCS-MD, we successfully generated 25 dissociation pathways of the complex, which showed complete or partial unfolding of the helical region of TAD-p53 during the dissociation process within an average simulation time of 154.8 ± 46.4 ns. The standard binding free energy obtained in combination with one-dimensional-, three-dimensional (3D)- or Cα-MSM was in good agreement with those determined experimentally. Using 3D-MSM based on the center of mass position of TAD-p53 relative to MDM2, the dissociation rate constant was calculated, which was comparable to those measured experimentally. Cα-MSM based on all Cα coordinates of TAD-p53 reproduced the experimentally measured standard binding free energy, and dissociation and association rate constants. We conclude that the combination of PaCS-MD and MSM offers an efficient computational procedure to calculate binding free energies and kinetic rates.

In silico study of Bombyx mori fibroin enhancement by graphene in acidic environment.

Bombyx mori fibroin has been widely used since a long time ago and has become a popular material. Here, we carry out a molecular dynamics simulation-based docking simulation of a small fragment of graphene in order to seek the best binding position on the N-termini domain of Bombyx mori fibroin. We report the best binding position, of which binding free energy falls at -54.8 kJ mol-1, indicating the strong binding. The further analysis of the binding pathway shows that this position is selective for single layered graphene rather than multi-layered graphene within our limited simulation times. Via comparing the RAMAN spectra of the corresponding binding pose of atomic clusters, we report the change in the bands compared with free standing graphene fragments, implying the change in molecular orbitals.

Vibrational Energy Transfer from Heme through Atomic Contacts in Proteins.

A pathway of vibrational energy flow in myoglobin was studied by time-resolved anti-Stokes ultraviolet resonance Raman spectroscopy combined with site-directed mutagenesis. Our previous study suggested that atomic contacts in proteins provide the dominant pathway for energy transfer while covalent bonds do not. In the present study, we directly examined the contributions of covalent bonds and atomic contacts to the pathway of vibrational energy flow by comparing the anti-Stokes resonance Raman spectra of two myoglobin mutants: one lacked a covalent bond between heme and the polypeptide chain, and the other retained the intact bond. The two mutants showed no significant difference in temporal changes in the anti-Stokes Raman intensities of the tryptophan bands, implying that the dominant channel of vibrational energy transfer is not through the covalent bond but rather through van der Waals atomic contacts between heme and the protein moiety. The obtained insights contribute to our general understanding of energy transfer in the condensed phase.

Protein-Ligand Dissociation Simulated by Parallel Cascade Selection Molecular Dynamics.

We investigated the dissociation process of tri-N-acetyl-d-glucosamine from hen egg white lysozyme using parallel cascade selection molecular dynamics (PaCS-MD), which comprises cycles of multiple unbiased MD simulations using a selection of MD snapshots as the initial structures for the next cycle. Dissociation was significantly accelerated by PaCS-MD, in which the probability of rare event occurrence toward dissociation was enhanced by the selection and rerandomization of the initial velocities. Although this complex was stable during 1 µs of conventional MD, PaCS-MD easily induced dissociation within 100-101 ns. We found that velocity rerandomization enhances the dissociation of triNAG from the bound state, whereas diffusion plays a more important role in the unbound state. We calculated the dissociation free energy by analyzing all PaCS-MD trajectories using the Markov state model (MSM), compared the results to those obtained by combinations of PaCS-MD and umbrella sampling (US), steered MD (SMD) and US, and SMD and the Jarzynski equality, and experimentally determined binding free energy. PaCS-MD/MSM yielded results most comparable to the experimentally determined binding free energy, independent of simulation parameter variations, and also gave the lowest standard errors.

Ligand binding to anti-cancer target CD44 investigated by molecular simulations.

CD44 is a cell-surface glycoprotein and receptor for hyaluronan, one of the major components of the tumor extracellular matrix. There is evidence that the interaction between CD44 and hyaluronan promotes breast cancer metastasis. Recently, the molecule F-19848A was shown to inhibit hyaluronan binding to receptor CD44 in a cell-based assay. In this study, we investigated the mechanism and energetics of F-19848A binding to CD44 using molecular simulation. Using the molecular mechanics/Poisson Boltzmann surface area (MM-PBSA) method, we obtained the binding free energy and inhibition constant of the complex. The van der Waals (vdW) interaction and the extended portion of F-19848A play key roles in the binding affinity. We screened natural products from a traditional Chinese medicine database to search for CD44 inhibitors. From combining pharmaceutical requirements with docking and molecular dynamics simulations, we found ten compounds that are potentially better or equal to the F-19848A ligand at binding to CD44 receptor. Therefore, we have identified new candidates of CD44 inhibitors, based on molecular simulation, which may be effective small molecules for the therapy of breast cancer.